From: abn Date: Mon, 15 May 2023 12:24:30 +0000 (+0200) Subject: never say never ... X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=refs%2Fheads%2Fabn%2Fmpi4py;p=tools%2Fmedcoupling.git never say never ... --- diff --git a/src/ParaMEDMEM_Swig/ParaMEDMEMCommon.i b/src/ParaMEDMEM_Swig/ParaMEDMEMCommon.i index 87fe9941e..72e59b989 100644 --- a/src/ParaMEDMEM_Swig/ParaMEDMEMCommon.i +++ b/src/ParaMEDMEM_Swig/ParaMEDMEMCommon.i @@ -297,6 +297,7 @@ namespace MEDCoupling public: InterpKernelDEC(); InterpKernelDEC(ProcessorGroup& source_group, ProcessorGroup& target_group); + InterpKernelDEC(const std::set& src_ids, const std::set& trg_ids); // hide last optional parameter! virtual ~InterpKernelDEC(); void release(); @@ -309,6 +310,12 @@ namespace MEDCoupling void prepareTargetDE(); %extend { + // Provides a direct ctor for which the communicator can be passed with "MPI._addressof(the_com)": + InterpKernelDEC(const std::set& src_ids, const std::set& trg_ids, long long comm_ptr) + { + return new InterpKernelDEC(src_ids, trg_ids, *((MPI_Comm*)comm_ptr)); + } + // 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. @@ -316,8 +323,8 @@ namespace MEDCoupling { return new InterpKernelDEC(source_group,target_group); } - - static InterpKernelDEC* _NewWithComm_internal(const std::set& src_ids, const std::set& trg_ids, long another_comm) + + static InterpKernelDEC* _NewWithComm_internal(const std::set& src_ids, const std::set& trg_ids, long long another_comm) { return new InterpKernelDEC(src_ids,trg_ids, *(MPI_Comm*)another_comm); // I know, ugly cast ... } @@ -327,7 +334,7 @@ namespace MEDCoupling class OverlapDEC : public DEC, public INTERP_KERNEL::InterpolationOptions { public: - OverlapDEC(const std::set& procIds); + OverlapDEC(const std::set& procIds); // hide optional param comm virtual ~OverlapDEC(); void release(); @@ -350,10 +357,15 @@ namespace MEDCoupling void debugPrintWorkSharing(std::ostream & ostr) const; %extend { + OverlapDEC(const std::set& ids, long long comm_ptr) + { + return new OverlapDEC(ids, *((MPI_Comm*)comm_ptr)); + } + // 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) + static OverlapDEC* _NewWithComm_internal(const std::set& ids, long long another_comm) { return new OverlapDEC(ids, *(MPI_Comm*)another_comm); // I know, ugly cast ... } @@ -405,6 +417,7 @@ else: %pythoncode %{ # And here we use mpi4py ability to provide its internal (C++) pointer to the communicator: +# NB: doing a proper typemap from MPI_Comm from Python to C++ requires the inclusion of mpi4py headers and .i file ... an extra dependency ... def _IKDEC_WithComm_internal(src_procs, tgt_procs, mpicomm=None): from mpi4py import MPI # Check iterable: @@ -413,6 +426,9 @@ def _IKDEC_WithComm_internal(src_procs, tgt_procs, mpicomm=None): except: s, t = None, None msg = "InterpKernelDEC: invalid type in ctor arguments! Possible signatures are:\n" + msg += " - InterpKernelDEC(ProcessorGroup, ProcessorGroup)\n" + msg += " - InterpKernelDEC(, )\n" + msg += " - InterpKernelDEC(, , MPI_Comm*) : WARNING here the address of the communicator should be passed with MPI._addressof(the_com)\n" msg += " - InterpKernelDEC.New(ProcessorGroup, ProcessorGroup)\n" msg += " - InterpKernelDEC.New(, )\n" msg += " - InterpKernelDEC.New(, , MPI_Comm)\n" @@ -436,6 +452,8 @@ def _ODEC_WithComm_internal(procs, mpicomm=None): msg = "OverlapDEC: invalid type in ctor arguments! Possible signatures are:\n" msg += " - OverlapDEC.New()\n" msg += " - OverlapDEC.New(, MPI_Comm)\n" + msg += " - OverlapDEC()\n" + msg += " - OverlapDEC(, MPI_Comm*) : WARNING here the address of the communicator should be passed with MPI._addressof(the_com)\n" raise InterpKernelException(msg) if mpicomm is None: return OverlapDEC(g) diff --git a/src/ParaMEDMEM_Swig/test_InterpKernelDEC.py b/src/ParaMEDMEM_Swig/test_InterpKernelDEC.py index e9c5a0d21..fe30d14cb 100755 --- a/src/ParaMEDMEM_Swig/test_InterpKernelDEC.py +++ b/src/ParaMEDMEM_Swig/test_InterpKernelDEC.py @@ -93,18 +93,23 @@ class ParaMEDMEM_IK_DEC_Tests(unittest.TestCase): l1, l2 = range(nproc_source), range(size - nproc_source, size) # With 2 iterables: i1 = InterpKernelDEC.New(l1, l2) + # Should also work directly: + i2 = InterpKernelDEC(l1, l2) # With 2 proc groups: interface = CommInterface() source_group = MPIProcessorGroup(interface, list(l1)) target_group = MPIProcessorGroup(interface, list(l2)) - i2 = InterpKernelDEC.New(source_group, target_group) + i3 = InterpKernelDEC.New(source_group, target_group) # Should also work directly: - i3 = InterpKernelDEC(source_group, target_group) + i4 = InterpKernelDEC(source_group, target_group) # With 2 iterables and a custom comm: - i4 = InterpKernelDEC.New(l1, l2, MPI.COMM_WORLD) + i5 = InterpKernelDEC.New(l1, l2, MPI.COMM_WORLD) + # Work directly with the **hack** + i6 = InterpKernelDEC(l1, l2, MPI._addressof(MPI.COMM_WORLD)) # Should fail with 2 proc groups **and** a communicator self.assertRaises(InterpKernelException, InterpKernelDEC.New, source_group, target_group, MPI.COMM_WORLD) - i4.release(); i3.release(); i2.release(); i1.release() + self.assertRaises(NotImplementedError, InterpKernelDEC, source_group, target_group, MPI.COMM_WORLD) + i6.release(); i5.release(); i4.release(); i3.release(); i2.release(); i1.release() source_group.release() target_group.release() diff --git a/src/ParaMEDMEM_Swig/test_OverlapDEC.py b/src/ParaMEDMEM_Swig/test_OverlapDEC.py index 80f075686..778f40b03 100755 --- a/src/ParaMEDMEM_Swig/test_OverlapDEC.py +++ b/src/ParaMEDMEM_Swig/test_OverlapDEC.py @@ -102,7 +102,10 @@ class ParaMEDMEM_O_DEC_Tests(unittest.TestCase): o2 = OverlapDEC(proc_group) # With an iterable and a custom comm: o3 = OverlapDEC.New(proc_group, MPI.COMM_WORLD) - o3.release(); o2.release(); o1.release() + # Also work directly with the **hack** on the comm: + o4 = OverlapDEC(proc_group, MPI._addressof(MPI.COMM_WORLD)) + self.assertRaises(NotImplementedError, OverlapDEC, proc_group, MPI.COMM_WORLD) + o4.release(); o3.release(); o2.release(); o1.release() @WriteInTmpDir def testOverlapDEC_2D_py_1(self):