From: abn Date: Fri, 12 May 2023 15:46:03 +0000 (+0200) Subject: second try without mpi4py.h ... X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=014b609a87c18bdd2489008942d6cf019ad2fe5e;p=tools%2Fmedcoupling.git second try without mpi4py.h ... --- diff --git a/CMakeLists.txt b/CMakeLists.txt index 181e98cb3..952609497 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -220,14 +220,6 @@ IF(MEDCOUPLING_ENABLE_PYTHON) IF(MEDCOUPLING_SWIG4_COMPAT) LIST(APPEND CMAKE_SWIG_FLAGS "-DMEDCOUPLING_SWIG4_COMPAT") ENDIF() - - # Possibility to specify a mpi4py installation - SET(MEDCOUPLING_MPI4PY_DIR "" CACHE PATH "Path to the mpi4py installation") - IF(MEDCOUPLING_MPI4PY_DIR) - MESSAGE(STATUS "Using mpi4py directory: ${MEDCOUPLING_MPI4PY_DIR}") - SET(MEDCOUPLING_MPI4PY_INCLUDE_DIR "${MEDCOUPLING_MPI4PY_DIR}/include") - SET(MEDCOUPLING_MPI4PY_I_FILE "${MEDCOUPLING_MPI4PY_DIR}/include/mpi4py/mpi4py.i") - ENDIF() ENDIF(MEDCOUPLING_ENABLE_PYTHON) IF(MEDCOUPLING_BUILD_DOC) diff --git a/doc/developer/doxygen/doxfiles/reference/misc/swig_mpi.dox b/doc/developer/doxygen/doxfiles/reference/misc/swig_mpi.dox index 84004333e..dc8898b45 100644 --- a/doc/developer/doxygen/doxfiles/reference/misc/swig_mpi.dox +++ b/doc/developer/doxygen/doxfiles/reference/misc/swig_mpi.dox @@ -21,4 +21,11 @@ void method_taking_a_MPI_Comm(MPI_Comm mpicomm); where the path to mpi4py.i needs to be adapted to your configuration. +In MEDCoupling DEC APIs (InterpKernelDEC and OverlapDEC constructors), the Python interface deals with that using the following trick to get the C++ adress of the MPI communicator: + +\code{.py} +from mpi4py import MPI +MPI._addressof(mpicomm) # returns the C++ adress of the MPI communicator +\endcode + */ diff --git a/src/ParaMEDMEM/InterpKernelDEC.hxx b/src/ParaMEDMEM/InterpKernelDEC.hxx index 95c40f705..f4b156329 100644 --- a/src/ParaMEDMEM/InterpKernelDEC.hxx +++ b/src/ParaMEDMEM/InterpKernelDEC.hxx @@ -130,8 +130,7 @@ namespace MEDCoupling public: InterpKernelDEC(); InterpKernelDEC(ProcessorGroup& source_group, ProcessorGroup& target_group); - InterpKernelDEC(const std::set& src_ids, const std::set& trg_ids, - const MPI_Comm& world_comm=MPI_COMM_WORLD); + InterpKernelDEC(const std::set& src_ids, const std::set& trg_ids, const MPI_Comm& world_comm=MPI_COMM_WORLD); virtual ~InterpKernelDEC(); void release(); diff --git a/src/ParaMEDMEM_Swig/CMakeLists.txt b/src/ParaMEDMEM_Swig/CMakeLists.txt index 28a11ad5e..47320923d 100644 --- a/src/ParaMEDMEM_Swig/CMakeLists.txt +++ b/src/ParaMEDMEM_Swig/CMakeLists.txt @@ -36,10 +36,6 @@ SET(SWIG_MODULE_ParaMEDMEM_EXTRA_FLAGS "${NUMPY_DEFINITIONS};${SCIPY_DEFINITIONS IF(MEDCOUPLING_USE_64BIT_IDS) STRING(APPEND SWIG_MODULE_ParaMEDMEM_EXTRA_FLAGS ";-DMEDCOUPLING_USE_64BIT_IDS") ENDIF(MEDCOUPLING_USE_64BIT_IDS) -IF(MEDCOUPLING_MPI4PY_DIR) - INCLUDE_DIRECTORIES(${MEDCOUPLING_MPI4PY_INCLUDE_DIR}) - STRING(APPEND SWIG_MODULE_ParaMEDMEM_EXTRA_FLAGS ";-DMEDCOUPLING_MPI4PY_I_FILE=${MEDCOUPLING_MPI4PY_I_FILE}") -ENDIF() INCLUDE_DIRECTORIES( ${PYTHON_INCLUDE_DIRS} diff --git a/src/ParaMEDMEM_Swig/ParaMEDMEMCommon.i b/src/ParaMEDMEM_Swig/ParaMEDMEMCommon.i index e8d1cf114..50c42b136 100644 --- a/src/ParaMEDMEM_Swig/ParaMEDMEMCommon.i +++ b/src/ParaMEDMEM_Swig/ParaMEDMEMCommon.i @@ -60,9 +60,7 @@ using namespace ICoCo; %include "ComponentTopology.hxx" %include "DEC.hxx" %include "DisjointDEC.hxx" -%include "InterpKernelDEC.hxx" %include "StructuredCoincidentDEC.hxx" -%include "OverlapDEC.hxx" %newobject MEDCoupling::ParaUMesh::New; %newobject MEDCoupling::ParaUMesh::getMesh; @@ -81,6 +79,9 @@ using namespace ICoCo; %newobject MEDCoupling::ParaSkyLineArray::getSkyLineArray; %newobject MEDCoupling::ParaSkyLineArray::getGlobalIdsArray; +%newobject MEDCoupling::InterpKernelDEC::_NewWithComm_internal; +%newobject MEDCoupling::OverlapDEC::_NewWithComm_internal; + %feature("unref") ParaSkyLineArray "$this->decrRef();" %feature("unref") ParaUMesh "$this->decrRef();" %feature("unref") ParaDataArrayInt32 "$this->decrRef();" @@ -280,16 +281,80 @@ namespace MEDCoupling } } }; -} -/* This object can be used only if MED_ENABLE_FVM is defined*/ -#ifdef MED_ENABLE_FVM -class NonCoincidentDEC : public DEC -{ -public: - NonCoincidentDEC(ProcessorGroup& source, ProcessorGroup& target); -}; -#endif + /* This object can be used only if MED_ENABLE_FVM is defined*/ + #ifdef MED_ENABLE_FVM + class NonCoincidentDEC : public DEC + { + public: + NonCoincidentDEC(ProcessorGroup& source, ProcessorGroup& target); + }; + #endif + + class InterpKernelDEC : public DisjointDEC, public INTERP_KERNEL::InterpolationOptions + { + public: + InterpKernelDEC(); + InterpKernelDEC(ProcessorGroup& source_group, ProcessorGroup& target_group); + virtual ~InterpKernelDEC(); + void release(); + + void synchronize(); + void recvData(); + void recvData(double time); + void sendData(); + void sendData(double time , double deltatime); + void prepareSourceDE(); + void prepareTargetDE(); + + %extend { + // This one should really not be called directly by the user since it still has an interface with a pointer to MPI_Comm + // which Swig doesn't handle nicely. + // It is just here to provide a constructor taking a **pointer** to a comm - See pythoncode below. + static InterpKernelDEC* _NewWithComm_internal(const std::set& src_ids, const std::set& trg_ids, long another_comm) + { + return new InterpKernelDEC(src_ids,trg_ids, *(MPI_Comm*)another_comm); // I know, ugly cast ... + } + } + }; + + class OverlapDEC : public DEC, public INTERP_KERNEL::InterpolationOptions + { + public: + OverlapDEC(const std::set& procIds); + virtual ~OverlapDEC(); + void release(); + + void sendRecvData(bool way=true); + void sendData(); + void recvData(); + void synchronize(); + void attachSourceLocalField(ParaFIELD *field, bool ownPt=false); + void attachTargetLocalField(ParaFIELD *field, bool ownPt=false); + void attachSourceLocalField(MEDCouplingFieldDouble *field); + void attachTargetLocalField(MEDCouplingFieldDouble *field); + void attachSourceLocalField(ICoCo::MEDDoubleField *field); + void attachTargetLocalField(ICoCo::MEDDoubleField *field); + ProcessorGroup *getGroup(); + bool isInGroup() const; + + void setDefaultValue(double val); + void setWorkSharingAlgo(int method); + + void debugPrintWorkSharing(std::ostream & ostr) const; + + %extend { + // This one should really not be called directly by the user since it still has an interface with a pointer to MPI_Comm + // which Swig doesn't handle nicely. + // It is just here to provide a constructor taking a **pointer** to a comm - See pythoncode below. + static OverlapDEC* _NewWithComm_internal(const std::set& ids, long another_comm) + { + return new OverlapDEC(ids, *(MPI_Comm*)another_comm); // I know, ugly cast ... + } + } + }; + +} // end namespace MEDCoupling %extend MEDCoupling::ParaMESH { @@ -330,3 +395,19 @@ if MEDCouplingUse64BitIDs(): else: ParaDataArrayInt = ParaDataArrayInt32 %} + +%pythoncode %{ + +# And here we use mpi4py ability to provide its internal (C++) pointer to the communicator: +def _InterpKernelDEC_WithComm_internal(src_procs, tgt_procs, mpicomm): + from mpi4py import MPI + return InterpKernelDEC._NewWithComm_internal(src_procs, tgt_procs, MPI._addressof(mpicomm)) + +def _OverlapDEC_WithComm_internal(procs, mpicomm): + from mpi4py import MPI + return OverlapDEC._NewWithComm_internal(procs, MPI._addressof(mpicomm)) + +InterpKernelDEC.NewWithCustomComm = _InterpKernelDEC_WithComm_internal +OverlapDEC.NewWithCustomComm = _OverlapDEC_WithComm_internal + +%} diff --git a/src/ParaMEDMEM_Swig/test_InterpKernelDEC.py b/src/ParaMEDMEM_Swig/test_InterpKernelDEC.py index ef68d76a0..79f8686b4 100755 --- a/src/ParaMEDMEM_Swig/test_InterpKernelDEC.py +++ b/src/ParaMEDMEM_Swig/test_InterpKernelDEC.py @@ -101,7 +101,8 @@ class ParaMEDMEM_IK_DEC_Tests(unittest.TestCase): interface = CommInterface() source_group = MPIProcessorGroup(interface, procs_source) target_group = MPIProcessorGroup(interface, procs_target) - idec = InterpKernelDEC(source_group, target_group) + #idec = InterpKernelDEC(source_group, target_group) + idec = InterpKernelDEC.NewWithCustomComm(procs_source, procs_target, MPI.COMM_WORLD) # Write out full size meshes/fields for inspection if rank == 0: diff --git a/src/ParaMEDMEM_Swig/test_OverlapDEC.py b/src/ParaMEDMEM_Swig/test_OverlapDEC.py index 48df8f4dc..3cbc2ffc3 100755 --- a/src/ParaMEDMEM_Swig/test_OverlapDEC.py +++ b/src/ParaMEDMEM_Swig/test_OverlapDEC.py @@ -98,7 +98,8 @@ class ParaMEDMEM_O_DEC_Tests(unittest.TestCase): # Define (single) processor group - note the difference with InterpKernelDEC which needs two groups. proc_group = list(range(size)) # No need for ProcessorGroup object here. - odec = OverlapDEC(proc_group) + #odec = OverlapDEC(proc_group) + odec = OverlapDEC.NewWithCustomComm(proc_group, MPI.COMM_WORLD) # Write out full size meshes/fields for inspection if rank == 0: diff --git a/src/PyWrapping/CMakeLists.txt b/src/PyWrapping/CMakeLists.txt index 94ab53ebd..5bf9f8181 100644 --- a/src/PyWrapping/CMakeLists.txt +++ b/src/PyWrapping/CMakeLists.txt @@ -36,10 +36,6 @@ SET(SWIG_MODULE_medcoupling_EXTRA_FLAGS "${NUMPY_DEFINITIONS};${SCIPY_DEFINITION IF(MEDCOUPLING_USE_64BIT_IDS) STRING(APPEND SWIG_MODULE_medcoupling_EXTRA_FLAGS ";-DMEDCOUPLING_USE_64BIT_IDS") ENDIF(MEDCOUPLING_USE_64BIT_IDS) -IF(MEDCOUPLING_MPI4PY_DIR) - INCLUDE_DIRECTORIES(${MEDCOUPLING_MPI4PY_INCLUDE_DIR}) - STRING(APPEND SWIG_MODULE_medcoupling_EXTRA_FLAGS ";-DMEDCOUPLING_MPI4PY_I_FILE=${MEDCOUPLING_MPI4PY_I_FILE}") -ENDIF() SET(medcoupling_SWIG_DPYS_FILES medcoupling.i)