Salome HOME
ParaUMesh.redistributeCells implementation.
[tools/medcoupling.git] / src / ParaMEDMEM_Swig / ParaMEDMEMCommon.i
1 // Copyright (C) 2017-2020  CEA/DEN, EDF R&D
2 //
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.
7 //
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.
12 //
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
16 //
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
18 //
19
20 %include std_set.i
21
22 %template() std::set<int>;
23
24 %{
25 #include "CommInterface.hxx"
26 #include "ProcessorGroup.hxx"
27 #include "Topology.hxx"
28 #include "MPIProcessorGroup.hxx"
29 #include "DEC.hxx"
30 #include "InterpKernelDEC.hxx"
31 #include "NonCoincidentDEC.hxx"
32 #include "StructuredCoincidentDEC.hxx"
33 #include "ParaMESH.hxx"
34 #include "ParaFIELD.hxx"
35 #include "ICoCoMEDField.hxx"
36 #include "ComponentTopology.hxx"
37 #include "ParaUMesh.hxx"
38 #include "ParaSkyLineArray.hxx"
39
40 using namespace INTERP_KERNEL;
41 using namespace MEDCoupling;
42 using namespace ICoCo;
43 %}
44
45 %include "InterpolationOptions.hxx"
46 %include "ProcessorGroup.hxx"
47 %include "DECOptions.hxx"
48 %include "ParaMESH.hxx"
49 %include "ParaFIELD.hxx"
50 %include "MPIProcessorGroup.hxx"
51 %include "ComponentTopology.hxx"
52 %include "DEC.hxx"
53 %include "DisjointDEC.hxx"
54 %include "InterpKernelDEC.hxx"
55 %include "StructuredCoincidentDEC.hxx"
56
57 %include "ICoCoField.hxx"
58 %rename(ICoCoMEDField) ICoCo::MEDField;
59 %include "ICoCoMEDField.hxx"
60
61 %newobject MEDCoupling::ParaUMesh::New;
62 %newobject MEDCoupling::ParaUMesh::getMesh;
63 %newobject MEDCoupling::ParaUMesh::getGlobalCellIds;
64 %newobject MEDCoupling::ParaUMesh::getGlobalNodeIds;
65 %newobject MEDCoupling::ParaUMesh::getCellIdsLyingOnNodes;
66 %newobject MEDCoupling::ParaUMesh::redistributeCells;
67 %newobject MEDCoupling::ParaSkyLineArray::New;
68 %newobject MEDCoupling::ParaSkyLineArray::equiRedistribute;
69 %newobject MEDCoupling::ParaSkyLineArray::getSkyLineArray;
70 %newobject MEDCoupling::ParaSkyLineArray::getGlobalIdsArray;
71
72 %feature("unref") ParaSkyLineArray "$this->decrRef();"
73 %feature("unref") ParaUMesh "$this->decrRef();"
74
75 %nodefaultctor;
76
77 namespace MEDCoupling
78 {
79   class CommInterface
80   {
81   public:
82     CommInterface();
83     virtual ~CommInterface();
84     int worldSize() const;
85     int commSize(MPI_Comm comm, int* size) const;
86     int commRank(MPI_Comm comm, int* rank) const;
87     int commGroup(MPI_Comm comm, MPI_Group* group) const;
88     int groupIncl(MPI_Group group, int size, int* ranks, MPI_Group* group_output) const;
89     int commCreate(MPI_Comm comm, MPI_Group group, MPI_Comm* comm_output) const;
90     int groupFree(MPI_Group* group) const;
91     int commFree(MPI_Comm* comm) const;
92
93     int send(void* buffer, int count, MPI_Datatype datatype, int target, int tag, MPI_Comm comm) const;
94     int recv(void* buffer, int count, MPI_Datatype datatype, int source, int tag, MPI_Comm comm, MPI_Status* status) const;
95     int sendRecv(void* sendbuf, int sendcount, MPI_Datatype sendtype, 
96                  int dest, int sendtag, void* recvbuf, int recvcount, 
97                  MPI_Datatype recvtype, int source, int recvtag, MPI_Comm comm,
98                  MPI_Status* status);
99
100     int Isend(void* buffer, int count, MPI_Datatype datatype, int target,
101               int tag, MPI_Comm comm, MPI_Request *request) const;
102     int Irecv(void* buffer, int count, MPI_Datatype datatype, int source,
103               int tag, MPI_Comm comm, MPI_Request* request) const;
104
105     int wait(MPI_Request *request, MPI_Status *status) const;
106     int test(MPI_Request *request, int *flag, MPI_Status *status) const;
107     int requestFree(MPI_Request *request) const;
108     int waitany(int count, MPI_Request *array_of_requests, int *index, MPI_Status *status) const;
109     int testany(int count, MPI_Request *array_of_requests, int *index, int *flag, MPI_Status *status) const;
110     int waitall(int count, MPI_Request *array_of_requests, MPI_Status *array_of_status) const { return MPI_Waitall(count, array_of_requests, array_of_status); }
111     int testall(int count, MPI_Request *array_of_requests, int *flag, MPI_Status *array_of_status) const;
112     int waitsome(int incount, MPI_Request *array_of_requests,int *outcount, int *array_of_indices, MPI_Status *array_of_status) const;
113     int testsome(int incount, MPI_Request *array_of_requests, int *outcount,
114                  int *array_of_indices, MPI_Status *array_of_status) const;
115     int probe(int source, int tag, MPI_Comm comm, MPI_Status *status) const;
116     int Iprobe(int source, int tag, MPI_Comm comm, int *flag, MPI_Status *status) const;
117     int cancel(MPI_Request *request) const;
118     int testCancelled(MPI_Status *status, int *flag) const;
119     int barrier(MPI_Comm comm) const;
120     int errorString(int errorcode, char *string, int *resultlen) const;
121     int getCount(MPI_Status *status, MPI_Datatype datatype, int *count) const;
122
123     int broadcast(void* buffer, int count, MPI_Datatype datatype, int root, MPI_Comm comm) const;
124     int allGather(void* sendbuf, int sendcount, MPI_Datatype sendtype,
125                   void* recvbuf, int recvcount, MPI_Datatype recvtype,
126                   MPI_Comm comm) const;
127     int allToAll(void* sendbuf, int sendcount, MPI_Datatype sendtype,
128                  void* recvbuf, int recvcount, MPI_Datatype recvtype,
129                  MPI_Comm comm) const;
130     int allToAllV(void* sendbuf, int* sendcounts, int* senddispls,
131                   MPI_Datatype sendtype, void* recvbuf, int* recvcounts,
132                   int* recvdispls, MPI_Datatype recvtype, 
133                   MPI_Comm comm) const;
134
135     int reduce(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm) const;
136     int allReduce(void* sendbuf, void* recvbuf, int count, MPI_Datatype datatype, MPI_Op op, MPI_Comm comm) const;
137     %extend
138     {
139       PyObject *allGatherArrays(const DataArrayIdType *array) const
140       {
141         std::vector< MCAuto<DataArrayIdType> > ret;
142         self->allGatherArrays(MPI_COMM_WORLD,array,ret);
143         return convertFromVectorAutoObjToPyObj<DataArrayIdType>(ret,SWIGTITraits<mcIdType>::TI);
144       }
145
146       PyObject *allToAllArrays(PyObject *arrays) const
147       {
148         std::vector< DataArrayIdType * > arraysIn;
149         std::vector< MCAuto<DataArrayIdType> > arrayOut;
150         convertFromPyObjVectorOfObj<MEDCoupling::DataArrayIdType*>(arrays,SWIGTITraits<mcIdType>::TI,"DataArrayIdType",arraysIn);
151         std::vector< MCAuto<DataArrayIdType> > arraysIn2(FromVecToVecAuto<DataArrayIdType>(arraysIn));
152         self->allToAllArrays(MPI_COMM_WORLD,arraysIn2,arrayOut);
153         return convertFromVectorAutoObjToPyObj<DataArrayIdType>(arrayOut,SWIGTITraits<mcIdType>::TI);
154       }
155     }
156   };
157
158   class ParaUMesh : public RefCountObject
159   {
160   public:
161     static ParaUMesh *New(MEDCouplingUMesh *mesh, DataArrayIdType *globalCellIds, DataArrayIdType *globalNodeIds);
162     ParaUMesh *redistributeCells(const DataArrayIdType *globalCellIds) const;
163     %extend
164     {
165       ParaUMesh(MEDCouplingUMesh *mesh, DataArrayIdType *globalCellIds, DataArrayIdType *globalNodeIds)
166       {
167         return ParaUMesh::New(mesh,globalCellIds,globalNodeIds);
168       }
169
170       MEDCouplingUMesh *getMesh()
171       {
172         MEDCouplingUMesh *ret(self->getMesh());
173         if(ret) ret->incrRef();
174         return ret;
175       }
176
177       DataArrayIdType *getGlobalCellIds()
178       {
179         DataArrayIdType *ret(self->getGlobalCellIds());
180         if(ret) ret->incrRef();
181         return ret;
182       }
183
184       DataArrayIdType *getGlobalNodeIds()
185       {
186         DataArrayIdType *ret(self->getGlobalNodeIds());
187         if(ret) ret->incrRef();
188         return ret;
189       }
190
191       DataArrayIdType *getCellIdsLyingOnNodes(const DataArrayIdType *globalNodeIds, bool fullyIn) const
192       { 
193         MCAuto<DataArrayIdType> ret(self->getCellIdsLyingOnNodes(globalNodeIds,fullyIn));
194         return ret.retn();
195       }
196     }
197   };
198
199   class ParaSkyLineArray : public RefCountObject
200   {
201   public:
202     static ParaSkyLineArray *New(MEDCouplingSkyLineArray *ska, DataArrayIdType *globalIds);
203     %extend
204     {
205       ParaSkyLineArray(MEDCouplingSkyLineArray *ska, DataArrayIdType *globalIds)
206       {
207         return ParaSkyLineArray::New(ska,globalIds);
208       }
209
210       ParaSkyLineArray *equiRedistribute(mcIdType nbOfEntities) const
211       {
212         MCAuto<ParaSkyLineArray> ret(self->equiRedistribute(nbOfEntities));
213         return ret.retn();
214       }
215
216       MEDCouplingSkyLineArray *getSkyLineArray() const
217       {
218         MEDCouplingSkyLineArray *ret(self->getSkyLineArray());
219         if(ret)
220           ret->incrRef();
221         return ret;
222       }
223       
224       DataArrayIdType *getGlobalIdsArray() const
225       {
226         DataArrayIdType *ret(self->getGlobalIdsArray());
227         if(ret)
228           ret->incrRef();
229         return ret;
230       }
231     }
232   };
233 }
234
235 /* This object can be used only if MED_ENABLE_FVM is defined*/
236 #ifdef MED_ENABLE_FVM
237 class NonCoincidentDEC : public DEC
238 {
239 public:
240   NonCoincidentDEC(ProcessorGroup& source, ProcessorGroup& target);
241 };
242 #endif
243
244 %extend MEDCoupling::ParaMESH
245 {
246   PyObject *getGlobalNumberingCell2() const
247   {
248     const mcIdType *tmp=self->getGlobalNumberingCell();
249     mcIdType size=self->getCellMesh()->getNumberOfCells();
250     PyObject *ret=PyList_New(size);
251     for(mcIdType i=0;i<size;i++)
252       PyList_SetItem(ret,i,PyInt_FromLong(tmp[i])); 
253     return ret;
254   }
255
256   PyObject *getGlobalNumberingFace2() const
257   {
258     const mcIdType *tmp=self->getGlobalNumberingFace();
259     mcIdType size=self->getFaceMesh()->getNumberOfCells();
260     PyObject *ret=PyList_New(size);
261     for(mcIdType i=0;i<size;i++)
262       PyList_SetItem(ret,i,PyInt_FromLong(tmp[i])); 
263     return ret;
264   }
265
266   PyObject *getGlobalNumberingNode2() const
267   {
268     const mcIdType *tmp=self->getGlobalNumberingNode();
269     mcIdType size=self->getCellMesh()->getNumberOfNodes();
270     PyObject *ret=PyList_New(size);
271     for(mcIdType i=0;i<size;i++)
272       PyList_SetItem(ret,i,PyInt_FromLong(tmp[i])); 
273     return ret;
274   }
275 }