Thanks a lot dokken, your code works just fine. After that I made a few more modifications to my original code snippet.
The whole idea here is to pass a dict from python to c++ (especially if the mesh/functionspace remains same) and use c++ to calculate dof_component_map and get it back in python so that can be used later in a fluid structure-interaction code.
Note 1: I need this because:
- It is very expensive to call a recursive function each time
- I am interpolating multiple times in each time-step in the solver on a 1 million DOF function-space. However my fluid mesh remain constant with time.
Note 2: That is the reason I intend to pass it by reference. The whole idea of this interpolation script is copied from the fenicstools/interpolate module.
The modified code is given below. It works perfectly fine, except the fact that the dof_component_map is calculated correctly in c++ and does not get carried back to python. (One can check that by uncommenting the “printing map” lines in interpolate.cpp).
test.py
from dolfin import *
import cppimport
compiled_cpp_module = cppimport.imp('delta_interpolation')
mesh = UnitCubeMesh(1, 1, 1)
V = VectorFunctionSpace(mesh, 'CG', 2)
dof_component_map = {}
compiled_cpp_module.extract_dof_component_map_user(dof_component_map, V)
for k,v in dof_component_map.items():
print (k, ":", v)
delta_interpolation.cpp
/*
<%
from dolfin.jit.jit import dolfin_pc
setup_pybind11(cfg)
cfg['include_dirs'] = dolfin_pc['include_dirs']
cfg['library_dirs'] = dolfin_pc['library_dirs']
cfg['compiler_args'] = ['-std=c++11', '-DHAS_MPI']
%>
*/
#include <pybind11/pybind11.h>
#include <pybind11/stl_bind.h>
#include <pybind11/eigen.h>
#include <pybind11/numpy.h>
#include <pybind11/stl.h>
#include <dolfin/la/GenericVector.h>
#include <dolfin/function/Function.h>
#include <dolfin/function/FunctionSpace.h>
#include <dolfin/fem/GenericDofMap.h>
#include <dolfin/fem/FiniteElement.h>
#include <dolfin/common/RangedIndexSet.h>
#include <dolfin/geometry/BoundingBoxTree.h>
#include <dolfin/mesh/Edge.h>
#include <dolfin/mesh/Mesh.h>
#include <dolfin/mesh/Cell.h>
#include <dolfin/mesh/CellType.h>
#include <dolfin/mesh/MeshEntityIterator.h>
#include <dolfin/mesh/Vertex.h>
#include <dolfin/geometry/Point.h>
#include <dolfin/mesh/MeshEntity.h>
PYBIND11_MAKE_OPAQUE(std::unordered_map<std::size_t, std::size_t> &);
using namespace dolfin;
namespace py = pybind11;
void extract_dof_component_map(std::unordered_map<std::size_t,
std::size_t>& dof_component_map,
const FunctionSpace& V,
int* component)
{
// Extract sub dofmaps recursively and store dof to component map
if (V.element()->num_sub_elements() == 0)
{
std::unordered_map<std::size_t, std::size_t> collapsed_map;
std::shared_ptr<GenericDofMap> dummy
= V.dofmap()->collapse(collapsed_map, *V.mesh());
(*component)++;
for (const auto &map_it : collapsed_map)
dof_component_map[map_it.second] = (*component);
}
else
{
for (std::size_t i = 0; i < V.element()->num_sub_elements(); ++i)
{
const std::vector<std::size_t> comp = {i};
std::shared_ptr<FunctionSpace> Vs = V.extract_sub_space(comp);
extract_dof_component_map(dof_component_map, *Vs, component);
}
}}
void extract_dof_component_map_user(std::unordered_map<std::size_t,
std::size_t> &dof_component_map,
const FunctionSpace &V)
{
int component = -1;
extract_dof_component_map(dof_component_map, V, &component);
// Printing map
//for (auto &pair: dof_component_map)
//{
// std::cout << "{" << pair.first << ": " << pair.second << "}\n";
//}
}
PYBIND11_MODULE(delta_interpolation, m)
{
py::bind_map<std::unordered_map<std::size_t, std::size_t>>(m, "UnorderedMapSize_tSize_t", py::module_local(false));
m.def("extract_dof_component_map_user", (void (*)(std::unordered_map<std::size_t,
std::size_t> & , const FunctionSpace &)) &
extract_dof_component_map_user);
m.def("extract_dof_component_map_user", [](py::dict d, py::object U){
auto _d = d.cast<std::unordered_map<std::size_t, std::size_t>>();
auto _U = U.attr("_cpp_object").cast<const FunctionSpace &>();
extract_dof_component_map_user(_d, _U);
});
}
If we look here and also at (test_stl_binders.cpp), it looks like we should be able to pass dict by reference, maybe I am wrong, but I tried my best, for the life of me, I can’t seem to figure it out. They also mention that using “PYBIND11_MAKE_OPAQUE” macro helps us to do that. What makes it even more complicated is the fact that std:unordered_map uses “dynamic vectors - size_t” (because of the size of Function space).
Any help or even a workaround will be very helpful. Again I am using FEniCS 2019.1.0. Thanks in advance.
Also dokken, I dont have permission to delete the former post.