]> SALOME platform Git repositories - tools/medcoupling.git/commitdiff
Salome HOME
Let mpi4py do the MPI wrapping part job
authorAnthony Geay <anthony.geay@edf.fr>
Wed, 13 Dec 2017 13:46:45 +0000 (14:46 +0100)
committerAnthony Geay <anthony.geay@edf.fr>
Wed, 13 Dec 2017 13:46:45 +0000 (14:46 +0100)
src/ParaMEDMEM_Swig/CMakeLists.txt
src/ParaMEDMEM_Swig/ParaMEDMEM.i
src/ParaMEDMEM_Swig/ParaMEDMEM.typemap [deleted file]
src/ParaMEDMEM_Swig/test_InterpKernelDEC.py

index 17c86d796e8715fc3bcadd290c9ad72b5962f843..93403995c9b2eca31ce4bd61ef80ae41c400ee66 100644 (file)
@@ -30,9 +30,6 @@ ELSE()
 ENDIF()
 SET(SWIG_MODULE_ParaMEDMEM_EXTRA_FLAGS "${NUMPY_DEFINITIONS};${SCIPY_DEFINITIONS}")
 
-SET (ParaMEDMEM_SWIG_DPYS_FILES
-    ParaMEDMEM.typemap)
-
 INCLUDE_DIRECTORIES(
   ${PYTHON_INCLUDE_DIRS}
   ${NUMPY_INCLUDE_DIR}
@@ -46,7 +43,7 @@ INCLUDE_DIRECTORIES(
   ${CMAKE_CURRENT_SOURCE_DIR}/../INTERP_KERNEL/Bases
   )
 
-SET (SWIG_MODULE_ParaMEDMEM_EXTRA_DEPS ${ParaMEDMEM_SWIG_DPYS_FILES}
+SET (SWIG_MODULE_ParaMEDMEM_EXTRA_DEPS
     ${paramedmem_HEADERS_HXX}
     ${medloader_HEADERS_HXX}
     ${medcoupling_HEADERS_HXX} ${medcoupling_HEADERS_TXX}
index da9810e8a1a5edd965329c7585c423a28339f12d..d56de82a14d2bb67e50727ff1ccbc2bfb22bca69 100644 (file)
 #define MEDCOUPLING_EXPORT
 #define INTERPKERNEL_EXPORT
 
-%include "ParaMEDMEM.typemap"
 %include "MEDCouplingCommon.i"
 
+%include std_set.i
+%include std_string.i
+
+%template() std::set<int>;
+
 %{
 #include "CommInterface.hxx"
 #include "ProcessorGroup.hxx"
 #include "ICoCoMEDField.hxx"
 #include "ComponentTopology.hxx"
 
-#include <mpi.h>
-
 using namespace INTERP_KERNEL;
 using namespace MEDCoupling;
 using namespace ICoCo;
-      
-enum mpi_constants { mpi_comm_world, mpi_comm_self, mpi_double, mpi_int };
 %}
 
 %include "InterpolationOptions.hxx"
@@ -109,162 +109,6 @@ public:
   }
 }
 
-//=============================================================================================
-// Interface for MPI-realization-specific constants like MPI_COMM_WORLD.
-//
-// Type and values of constants like MPI_COMM_WORLD depends on MPI realization
-// and usually such constants actually are macros. To have such symbols in python
-// and translate them into correct values we use the following technique.
-// We define some constants (enum mpi_constants) and map them into real MPI values
-// using typemaps, and we create needed python symbols equal to 'mpi_constants'
-// via %pythoncode directive.
-
-// Constants corresponding to similar MPI definitions
-enum mpi_constants { mpi_comm_world, mpi_comm_self, mpi_double, mpi_int };
-
-// Map mpi_comm_world and mpi_comm_self -> MPI_COMM_WORLD and MPI_COMM_SELF
-%typemap(in) MPI_Comm
-{ 
-  switch (PyInt_AsLong($input))
-    {
-    case mpi_comm_world: $1 = MPI_COMM_WORLD; break;
-    case mpi_comm_self:  $1 = MPI_COMM_SELF;  break;
-    default:
-      PyErr_SetString(PyExc_TypeError,"unexpected value of MPI_Comm");
-      return NULL;
-    }
-}
-// Map mpi_double and mpi_int -> MPI_DOUBLE and MPI_INT
-%typemap(in) MPI_Datatype
-{ 
-  switch (PyInt_AsLong($input))
-    {
-    case mpi_double:     $1 = MPI_DOUBLE;     break;
-    case mpi_int:        $1 = MPI_INT;        break;
-    default:
-      PyErr_SetString(PyExc_TypeError,"unexpected value of MPI_Datatype");
-      return NULL;
-    }
-}
-// The following code gets inserted into the result python file:
-// create needed python symbols
-%pythoncode %{
-MPI_COMM_WORLD = mpi_comm_world
-MPI_COMM_SELF  = mpi_comm_self
-MPI_DOUBLE     = mpi_double
-MPI_INT        = mpi_int
-%}
-//=============================================================================================
-
-// ==============
-// MPI_Comm_size
-// ==============
-%inline %{ PyObject* MPI_Comm_size(MPI_Comm comm)
-  {
-    int res = 0;
-    int err = MPI_Comm_size(comm, &res);
-    if ( err != MPI_SUCCESS )
-      {
-        PyErr_SetString(PyExc_RuntimeError,"Erorr in MPI_Comm_size()");
-        return NULL;
-      }
-    return PyInt_FromLong( res );
-  } %}
-
-// ==============
-// MPI_Comm_rank
-// ==============
-%inline %{ PyObject* MPI_Comm_rank(MPI_Comm comm)
-  {
-    int res = 0;
-    int err = MPI_Comm_rank(comm, &res);
-    if ( err != MPI_SUCCESS )
-      {
-        PyErr_SetString(PyExc_RuntimeError,"Erorr in MPI_Comm_rank()");
-        return NULL;
-      }
-    return PyInt_FromLong( res );
-  } 
-  %}
-
-int MPI_Init(int *argc, char ***argv );
-int MPI_Barrier(MPI_Comm comm);
-int MPI_Finalize();
-
-// ==========
-// MPI_Bcast
-// ==========
-
-%inline %{ PyObject* MPI_Bcast(PyObject* buffer, int nb, MPI_Datatype type, int root, MPI_Comm c)
-  {
-    // buffer must be a list
-    if (!PyList_Check(buffer))
-      {
-        PyErr_SetString(PyExc_TypeError, "buffer is expected to be a list");
-        return NULL;
-      }
-    // check list size
-    int aSize = PyList_Size(buffer);
-    if ( aSize != nb )
-      {
-        std::ostringstream stream; stream << "buffer is expected to be of size " << nb;
-        PyErr_SetString(PyExc_ValueError, stream.str().c_str());
-        return NULL;
-      }
-    // allocate and fill a buffer
-    void* aBuf = 0;
-    int* intBuf = 0;
-    double* dblBuf = 0;
-    if ( type == MPI_DOUBLE )
-      {
-        aBuf = (void*) ( dblBuf = new double[ nb ] );
-        for ( int i = 0; i < aSize; ++i )
-          dblBuf[i] = PyFloat_AS_DOUBLE( PyList_GetItem( buffer, i ));
-      }
-    else if ( type == MPI_INT )
-      {
-        aBuf = (void*) ( intBuf = new int[ nb ] );
-        for ( int i = 0; i < aSize; ++i )
-          intBuf[i] = int( PyInt_AS_LONG( PyList_GetItem( buffer, i )));
-      }
-    else
-      {
-        PyErr_SetString(PyExc_TypeError, "Only MPI_DOUBLE and MPI_INT supported");
-        return NULL;
-      }
-    // call MPI_Bcast
-    int err = MPI_Bcast(aBuf, nb, type, root, c);
-    // treat error
-    if ( err != MPI_SUCCESS )
-      {
-        PyErr_SetString(PyExc_RuntimeError,"Erorr in MPI_Bcast()");
-        delete [] intBuf; delete [] dblBuf;
-        return NULL;
-      }
-    // put received data into the list
-    int pyerr = 0;
-    if ( type == MPI_DOUBLE )
-      {
-        for ( int i = 0; i < aSize && !pyerr; ++i )
-          pyerr = PyList_SetItem(buffer, i, PyFloat_FromDouble( dblBuf[i] ));
-        delete [] dblBuf;
-      }
-    else
-      {
-        for ( int i = 0; i < aSize && !pyerr; ++i )
-          pyerr = PyList_SetItem(buffer, i, PyInt_FromLong( intBuf[i] ));
-        delete [] intBuf;
-      }
-    if ( pyerr )
-      {
-        PyErr_SetString(PyExc_RuntimeError, "Error of PyList_SetItem()");
-        return NULL;
-      }
-    return PyInt_FromLong( err );
-
-  }
-  %}
-
 %pythoncode %{
 def MEDCouplingDataArrayDoubleIadd(self,*args):
     import _ParaMEDMEM
diff --git a/src/ParaMEDMEM_Swig/ParaMEDMEM.typemap b/src/ParaMEDMEM_Swig/ParaMEDMEM.typemap
deleted file mode 100644 (file)
index fa2a0ed..0000000
+++ /dev/null
@@ -1,89 +0,0 @@
-// Copyright (C) 2007-2016  CEA/DEN, EDF R&D
-//
-// This library is free software; you can redistribute it and/or
-// modify it under the terms of the GNU Lesser General Public
-// License as published by the Free Software Foundation; either
-// version 2.1 of the License, or (at your option) any later version.
-//
-// This library is distributed in the hope that it will be useful,
-// but WITHOUT ANY WARRANTY; without even the implied warranty of
-// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
-// Lesser General Public License for more details.
-//
-// You should have received a copy of the GNU Lesser General Public
-// License along with this library; if not, write to the Free Software
-// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
-//
-// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
-//
-
-%include std_set.i
-%include std_string.i
-
-%template() std::set<int>;
-
-// Creates "int *argc, char ***argv" parameters from input list
-%typemap(in) (int *argc, char ***argv) {
-  int i;
-  if (!PyList_Check($input)) {
-    PyErr_SetString(PyExc_ValueError, "Expecting a list");
-    return NULL;
-  }
-  int aSize = PyList_Size($input);
-  $1 = &aSize;
-  char** aStrs = (char **) malloc((aSize+1)*sizeof(char *));
-  for (i = 0; i < aSize; i++) {
-    PyObject *s = PyList_GetItem($input,i);
-    if (PyString_Check(s))
-      aStrs[i] = PyString_AsString(s);
-%#if PY_VERSION_HEX >= 0x03000000
-    else if (PyUnicode_Check(s))
-      aStrs[i] = PyUnicode_AsUTF8(s);
-%#endif
-    else {
-        free(aStrs);
-        PyErr_SetString(PyExc_ValueError, "List items must be strings");
-        return NULL;
-    }
-  }
-  aStrs[i] = 0;
-  $2 = &aStrs;
-}
-
-%typemap(freearg) (int *argc, char ***argv) {
-   if ($2) free(*($2));
-}
-
-/*  MACRO: IN typemap for std::set<TYPE> C++ object */
-%define TYPEMAP_INPUT_SET_BY_VALUE( TYPE )
-{
-  /* typemap in for set<TYPE> */
-  /* Check if is a list */
-  if (PyList_Check($input))
-  {
-    int size = PyList_Size($input);
-    std::set< TYPE > tmpSet;
-
-    for (int i=0; i < size; i++)
-    {
-      PyObject * tmp = PyList_GetItem($input,i);
-      TYPE elem = PyInt_AsLong(tmp);
-      tmpSet.insert(elem);
-    }
-    $1 = tmpSet;
-  }
-  else
-  {
-    PyErr_SetString(PyExc_TypeError,"not a list");
-    return NULL;
-  }
-}
-%enddef
-
-%typemap(in) std::set<int>
-{ 
-  TYPEMAP_INPUT_SET_BY_VALUE( int ) 
-}
-%typecheck(SWIG_TYPECHECK_POINTER) std::set<int> {
-  $1 = PyList_Check($input) ? 1 : 0;
-}
index 8a42cc3572658d21857319f565e1a43bb773932e..539d452c7972044994d963326b50d8eb1df0f810 100755 (executable)
@@ -58,7 +58,7 @@ class ParaMEDMEMBasicsTest(unittest.TestCase):
         filename_xml1 = os.path.join(data_dir, "share/resources/med/square1_split")
         filename_xml2 = os.path.join(data_dir, "share/resources/med/square2_split")
 
-        MPI_Barrier(MPI_COMM_WORLD)
+        MPI.COMM_WORLD.Barrier()
         if source_group.containsMyRank():
             filename = filename_xml1 + str(rank+1) + ".med"
             meshname = "Mesh_2_" + str(rank+1)
@@ -113,8 +113,8 @@ class ParaMEDMEMBasicsTest(unittest.TestCase):
         paramesh   =0
         parafield  =0
         icocofield =0
-        MPI_Barrier(MPI_COMM_WORLD)
-        MPI_Finalize()
+        MPI.COMM_WORLD.Barrier()
+        MPI.Finalize()
         pass
     pass