1 // Copyright (C) 2007-2016 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License, or (at your option) any later version.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
22 %include "ParaMEDMEM.typemap"
23 %include "MEDLoaderCommon.i"
26 #include "CommInterface.hxx"
27 #include "ProcessorGroup.hxx"
28 #include "Topology.hxx"
29 #include "MPIProcessorGroup.hxx"
31 #include "InterpKernelDEC.hxx"
32 #include "NonCoincidentDEC.hxx"
33 #include "StructuredCoincidentDEC.hxx"
34 #include "ParaMESH.hxx"
35 #include "ParaFIELD.hxx"
36 #include "ICoCoMEDField.hxx"
37 #include "ComponentTopology.hxx"
41 using namespace MEDCoupling;
42 using namespace ICoCo;
44 enum mpi_constants { mpi_comm_world, mpi_comm_self, mpi_double, mpi_int };
47 %include "CommInterface.hxx"
48 %include "ProcessorGroup.hxx"
49 %include "DECOptions.hxx"
50 %include "ParaMESH.hxx"
51 %include "ParaFIELD.hxx"
52 %include "MPIProcessorGroup.hxx"
53 %include "ComponentTopology.hxx"
55 %include "InterpKernelDEC.hxx"
56 %include "StructuredCoincidentDEC.hxx"
58 %rename(ICoCoMEDField) ICoCo::MEDField;
59 %include "ICoCoMEDField.hxx"
63 /* This object can be used only if MED_ENABLE_FVM is defined*/
65 class NonCoincidentDEC : public DEC
68 NonCoincidentDEC(ProcessorGroup& source, ProcessorGroup& target);
72 %extend MEDCoupling::ParaMESH
74 PyObject *getGlobalNumberingCell2() const
76 const int *tmp=self->getGlobalNumberingCell();
77 int size=self->getCellMesh()->getNumberOfCells();
78 PyObject *ret=PyList_New(size);
79 for(int i=0;i<size;i++)
80 PyList_SetItem(ret,i,PyInt_FromLong(tmp[i]));
84 PyObject *getGlobalNumberingFace2() const
86 const int *tmp=self->getGlobalNumberingFace();
87 int size=self->getFaceMesh()->getNumberOfCells();
88 PyObject *ret=PyList_New(size);
89 for(int i=0;i<size;i++)
90 PyList_SetItem(ret,i,PyInt_FromLong(tmp[i]));
94 PyObject *getGlobalNumberingNode2() const
96 const int *tmp=self->getGlobalNumberingNode();
97 int size=self->getCellMesh()->getNumberOfNodes();
98 PyObject *ret=PyList_New(size);
99 for(int i=0;i<size;i++)
100 PyList_SetItem(ret,i,PyInt_FromLong(tmp[i]));
105 //=============================================================================================
106 // Interface for MPI-realization-specific constants like MPI_COMM_WORLD.
108 // Type and values of constants like MPI_COMM_WORLD depends on MPI realization
109 // and usually such constants actually are macros. To have such symbols in python
110 // and translate them into correct values we use the following technique.
111 // We define some constants (enum mpi_constants) and map them into real MPI values
112 // using typemaps, and we create needed python symbols equal to 'mpi_constants'
113 // via %pythoncode directive.
115 // Constants corresponding to similar MPI definitions
116 enum mpi_constants { mpi_comm_world, mpi_comm_self, mpi_double, mpi_int };
118 // Map mpi_comm_world and mpi_comm_self -> MPI_COMM_WORLD and MPI_COMM_SELF
119 %typemap(in) MPI_Comm
121 switch (PyInt_AsLong($input))
123 case mpi_comm_world: $1 = MPI_COMM_WORLD; break;
124 case mpi_comm_self: $1 = MPI_COMM_SELF; break;
126 PyErr_SetString(PyExc_TypeError,"unexpected value of MPI_Comm");
130 // Map mpi_double and mpi_int -> MPI_DOUBLE and MPI_INT
131 %typemap(in) MPI_Datatype
133 switch (PyInt_AsLong($input))
135 case mpi_double: $1 = MPI_DOUBLE; break;
136 case mpi_int: $1 = MPI_INT; break;
138 PyErr_SetString(PyExc_TypeError,"unexpected value of MPI_Datatype");
142 // The following code gets inserted into the result python file:
143 // create needed python symbols
145 MPI_COMM_WORLD = mpi_comm_world
146 MPI_COMM_SELF = mpi_comm_self
147 MPI_DOUBLE = mpi_double
150 //=============================================================================================
155 %inline %{ PyObject* MPI_Comm_size(MPI_Comm comm)
158 int err = MPI_Comm_size(comm, &res);
159 if ( err != MPI_SUCCESS )
161 PyErr_SetString(PyExc_RuntimeError,"Erorr in MPI_Comm_size()");
164 return PyInt_FromLong( res );
170 %inline %{ PyObject* MPI_Comm_rank(MPI_Comm comm)
173 int err = MPI_Comm_rank(comm, &res);
174 if ( err != MPI_SUCCESS )
176 PyErr_SetString(PyExc_RuntimeError,"Erorr in MPI_Comm_rank()");
179 return PyInt_FromLong( res );
183 int MPI_Init(int *argc, char ***argv );
184 int MPI_Barrier(MPI_Comm comm);
191 %inline %{ PyObject* MPI_Bcast(PyObject* buffer, int nb, MPI_Datatype type, int root, MPI_Comm c)
193 // buffer must be a list
194 if (!PyList_Check(buffer))
196 PyErr_SetString(PyExc_TypeError, "buffer is expected to be a list");
200 int aSize = PyList_Size(buffer);
203 std::ostringstream stream; stream << "buffer is expected to be of size " << nb;
204 PyErr_SetString(PyExc_ValueError, stream.str().c_str());
207 // allocate and fill a buffer
211 if ( type == MPI_DOUBLE )
213 aBuf = (void*) ( dblBuf = new double[ nb ] );
214 for ( int i = 0; i < aSize; ++i )
215 dblBuf[i] = PyFloat_AS_DOUBLE( PyList_GetItem( buffer, i ));
217 else if ( type == MPI_INT )
219 aBuf = (void*) ( intBuf = new int[ nb ] );
220 for ( int i = 0; i < aSize; ++i )
221 intBuf[i] = int( PyInt_AS_LONG( PyList_GetItem( buffer, i )));
225 PyErr_SetString(PyExc_TypeError, "Only MPI_DOUBLE and MPI_INT supported");
229 int err = MPI_Bcast(aBuf, nb, type, root, c);
231 if ( err != MPI_SUCCESS )
233 PyErr_SetString(PyExc_RuntimeError,"Erorr in MPI_Bcast()");
234 delete [] intBuf; delete [] dblBuf;
237 // put recieved data into the list
239 if ( type == MPI_DOUBLE )
241 for ( int i = 0; i < aSize && !pyerr; ++i )
242 pyerr = PyList_SetItem(buffer, i, PyFloat_FromDouble( dblBuf[i] ));
247 for ( int i = 0; i < aSize && !pyerr; ++i )
248 pyerr = PyList_SetItem(buffer, i, PyInt_FromLong( intBuf[i] ));
253 PyErr_SetString(PyExc_RuntimeError, "Error of PyList_SetItem()");
256 return PyInt_FromLong( err );
262 def ParaMEDMEMDataArrayDoublenew(cls,*args):
264 return _ParaMEDMEM.DataArrayDouble____new___(cls,args)
265 def ParaMEDMEMDataArrayDoubleIadd(self,*args):
267 return _ParaMEDMEM.DataArrayDouble____iadd___(self, self, *args)
268 def ParaMEDMEMDataArrayDoubleIsub(self,*args):
270 return _ParaMEDMEM.DataArrayDouble____isub___(self, self, *args)
271 def ParaMEDMEMDataArrayDoubleImul(self,*args):
273 return _ParaMEDMEM.DataArrayDouble____imul___(self, self, *args)
274 def ParaMEDMEMDataArrayDoubleIdiv(self,*args):
276 return _ParaMEDMEM.DataArrayDouble____idiv___(self, self, *args)
277 def ParaMEDMEMDataArrayDoubleIpow(self,*args):
279 return _ParaMEDMEM.DataArrayDouble____ipow___(self, self, *args)
280 def ParaMEDMEMDataArrayDoubleTupleIadd(self,*args):
282 return _ParaMEDMEM.DataArrayDoubleTuple____iadd___(self, self, *args)
283 def ParaMEDMEMDataArrayDoubleTupleIsub(self,*args):
285 return _ParaMEDMEM.DataArrayDoubleTuple____isub___(self, self, *args)
286 def ParaMEDMEMDataArrayDoubleTupleImul(self,*args):
288 return _ParaMEDMEM.DataArrayDoubleTuple____imul___(self, self, *args)
289 def ParaMEDMEMDataArrayDoubleTupleIdiv(self,*args):
291 return _ParaMEDMEM.DataArrayDoubleTuple____idiv___(self, self, *args)
292 def MEDCouplingFieldDoublenew(cls,*args):
294 return _ParaMEDMEM.MEDCouplingFieldDouble____new___(cls,args)
295 def MEDCouplingFieldDoubleIadd(self,*args):
297 return _ParaMEDMEM.MEDCouplingFieldDouble____iadd___(self, self, *args)
298 def MEDCouplingFieldDoubleIsub(self,*args):
300 return _ParaMEDMEM.MEDCouplingFieldDouble____isub___(self, self, *args)
301 def MEDCouplingFieldDoubleImul(self,*args):
303 return _ParaMEDMEM.MEDCouplingFieldDouble____imul___(self, self, *args)
304 def MEDCouplingFieldDoubleIdiv(self,*args):
306 return _ParaMEDMEM.MEDCouplingFieldDouble____idiv___(self, self, *args)
307 def MEDCouplingFieldDoubleIpow(self,*args):
309 return _ParaMEDMEM.MEDCouplingFieldDouble____ipow___(self, self, *args)
310 def ParaMEDMEMDataArrayIntnew(cls,*args):
312 return _ParaMEDMEM.DataArrayInt____new___(cls,args)
313 def ParaMEDMEMDataArrayIntIadd(self,*args):
315 return _ParaMEDMEM.DataArrayInt____iadd___(self, self, *args)
316 def ParaMEDMEMDataArrayIntIsub(self,*args):
318 return _ParaMEDMEM.DataArrayInt____isub___(self, self, *args)
319 def ParaMEDMEMDataArrayIntImul(self,*args):
321 return _ParaMEDMEM.DataArrayInt____imul___(self, self, *args)
322 def ParaMEDMEMDataArrayIntIdiv(self,*args):
324 return _ParaMEDMEM.DataArrayInt____idiv___(self, self, *args)
325 def ParaMEDMEMDataArrayIntImod(self,*args):
327 return _ParaMEDMEM.DataArrayInt____imod___(self, self, *args)
328 def ParaMEDMEMDataArrayIntIpow(self,*args):
330 return _ParaMEDMEM.DataArrayInt____ipow___(self, self, *args)
331 def ParaMEDMEMDataArrayIntTupleIadd(self,*args):
333 return _ParaMEDMEM.DataArrayIntTuple____iadd___(self, self, *args)
334 def ParaMEDMEMDataArrayIntTupleIsub(self,*args):
336 return _ParaMEDMEM.DataArrayIntTuple____isub___(self, self, *args)
337 def ParaMEDMEMDataArrayIntTupleImul(self,*args):
339 return _ParaMEDMEM.DataArrayIntTuple____imul___(self, self, *args)
340 def ParaMEDMEMDataArrayIntTupleIdiv(self,*args):
342 return _ParaMEDMEM.DataArrayIntTuple____idiv___(self, self, *args)
343 def ParaMEDMEMDataArrayIntTupleImod(self,*args):
345 return _ParaMEDMEM.DataArrayIntTuple____imod___(self, self, *args)
348 %include "MEDCouplingFinalize.i"