Skip to content
Snippets Groups Projects
Commit 0c4d5696 authored by Christian Engwer's avatar Christian Engwer
Browse files

[python] simplify setup of yaspGrid with custom communicator, using new dune-common infrastructure

This uses the new feature of dune-common!1230
parent 9476506f
No related branches found
No related tags found
1 merge request!696[python] support yaspGrid construction with user defined mpi4py.MPI.Comm
......@@ -2,17 +2,10 @@
# SPDX-License-Identifier: LicenseRef-GPL-2.0-only-with-DUNE-exception
from dune.common.checkconfiguration import assertCMakeHave, ConfigurationError
from dune.common import comm as defaultCommunicator
from dune.typeregistry import generateTypeName
from dune.packagemetadata import getCMakeFlags
if getCMakeFlags()["HAVE_MPI"]:
from mpi4py import MPI
withMPI = True
defaultCommunicator = MPI.COMM_WORLD
else:
withMPI = False
defaultCommunicator = None
class CartesianDomain(tuple):
@staticmethod
def bndDomain(lower, upper, division, tol = 1e-4 ):
......@@ -312,24 +305,16 @@ def yaspGrid(constructor, dimgrid=None, coordinates="equidistant", ctype=None,
raise ValueError("yaspGrid: unsupported constructor parameter " + str(constructor))
# compile YaspGrid for a given dimension & coordinate type
includes = ["dune/grid/yaspgrid.hh", "dune/grid/io/file/dgfparser/dgfyasp.hh"]
includes = ["dune/grid/yaspgrid.hh", "dune/grid/io/file/dgfparser/dgfyasp.hh", "dune/common/parallel/mpihelper.hh"]
gridTypeName, _ = generateTypeName("Dune::YaspGrid", str(dimgrid), coordinates_type)
setupPeriodic = [ 'std::bitset<'+str(dimgrid)+'> periodic_;',
'for (int i=0;i<'+str(dimgrid)+';++i) periodic_.set(i,periodic[i]);']
if withMPI:
includes.append("mpi4py/mpi4py.h")
generator = setupPeriodic + ['import_mpi4py();', # import C symbol of mpi4py (see https://enccs.se/news/2021/03/mpi-hybrid-c-python-code/)
'MPI_Comm* comm_p = PyMPIComm_Get(py_comm.ptr());',
'Dune::Communication<MPI_Comm> communicator(*comm_p);',
'return new DuneType(coordinates,periodic_,overlap,communicator);']
else:
generator = setupPeriodic + [ 'return new DuneType(coordinates,periodic_,overlap);']
generator = setupPeriodic + [ 'return new DuneType(coordinates,periodic_,overlap,communicator);' ]
ctor = Constructor(
[ "const " + coordinates_type + "& coordinates",
'std::array<bool, '+str(dimgrid)+'> periodic',
'int overlap',
'pybind11::object py_comm'],
'Dune::Communication< Dune::MPIHelper::MPICommunicator > communicator'],
generator,
[ '"coordinates"_a', '"periodic"_a', '"overlap"_a', '"communicator"_a' ])
gridModule = module(includes, gridTypeName, ctor)
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment