Salome HOME
A forgotten C++ test
[tools/medcoupling.git] / src / ParaMEDMEM / OverlapDEC.cxx
index c799c23aaeff315b73718fea1fc82766e7873943..af842fb9dd08f600d6886b8afa3d109c564d4660 100644 (file)
@@ -20,6 +20,7 @@
 
 #include "OverlapDEC.hxx"
 #include "CommInterface.hxx"
+#include "ParaMESH.hxx"
 #include "ParaFIELD.hxx"
 #include "MPIProcessorGroup.hxx"
 #include "OverlapElementLocator.hxx"
@@ -31,25 +32,29 @@ namespace ParaMEDMEM
     \anchor OverlapDEC-det
     \class OverlapDEC
 
+    \section OverlapDEC-over Overview
+
     The \c OverlapDEC enables the \ref InterpKerRemapGlobal "conservative remapping" of fields between
     two parallel codes. This remapping is based on the computation of intersection volumes on
-    a \b same \b processor \b group. On this processor group are defined two field-templates called A
+    a \b single \b processor \b group. On this processor group are defined two field-templates called A
     and B. The computation is possible for 3D meshes, 2D meshes, 3D-surface meshes, 1D meshes and
     2D-curve meshes. Dimensions must be similar for the distribution templates A and B.
-    The main difference with \ref InterpKernelDEC-det is that this \ref para-dec "DEC" manages 2 field templates
-    on each processor of the processor group (A and B) called source and target.
-    Furthermore all processors in processor group cooperates in global interpolation matrix
-    computation. In this respect \ref InterpKernelDEC is a specialization of \c OverlapDEC.
 
-    \section ParaMEDMEMOverlapDECAlgorithmDescription Algorithm Description
+    The main difference with \ref InterpKernelDEC-det "InterpKernelDEC" is that this
+    \ref para-dec "DEC" works with a *single* processor group, in which processors will share the work.
+    Consequently each processor manages two \ref MEDCouplingFieldTemplatesPage "field templates" (A and B)
+    called source and target.
+    Furthermore all processors in the processor group cooperate in the global interpolation matrix
+    computation. In this respect \c InterpKernelDEC is a specialization of \c OverlapDEC.
+
+    \section ParaMEDMEMOverlapDECAlgorithmDescription Algorithm description
 
     Let's consider the following use case that is ran in ParaMEDMEMTest_OverlapDEC.cxx to describes
     the different steps of the computation. The processor group contains 3 processors.
     \anchor ParaMEDMEMOverlapDECImgTest1
-    \image html OverlapDEC1.png "Example showing the use case in order to explain the different steps."
+    \image html OverlapDEC1.png "Example split of the source and target mesh among the 3 procs"
 
-    \subsection ParaMEDMEMOverlapDECAlgoStep1 Step 1 : Bounding box exchange and global interaction
-    between procs computation.
+    \subsection ParaMEDMEMOverlapDECAlgoStep1 Step 1 : Bounding box exchange and global interaction between procs computation.
 
     In order to reduce as much as possible the amount of communications between distant processors,
     every processor computes a bounding box for A and B. Then a AllToAll communication is performed
@@ -161,7 +166,7 @@ namespace ParaMEDMEM
     the \b local TODO list per proc is expected to
     be as well balanced as possible.
 
-    The interpolation is performed as \ref ParaMEDMEM::MEDCouplingRemapper "Remapper" does.
+    The interpolation is performed as the \ref ParaMEDMEM::MEDCouplingRemapper "remapper" does.
 
     This operation is performed by OverlapInterpolationMatrix::addContribution method.
 
@@ -203,9 +208,13 @@ namespace ParaMEDMEM
 
     The method in charge to perform this is : ParaMEDMEM::OverlapMapping::prepare.
 */
-  OverlapDEC::OverlapDEC(const std::set<int>& procIds, const MPI_Comm& world_comm):_own_group(true),_interpolation_matrix(0),
-                                                                                   _source_field(0),_own_source_field(false),
-                                                                                   _target_field(0),_own_target_field(false)
+  OverlapDEC::OverlapDEC(const std::set<int>& procIds, const MPI_Comm& world_comm):
+      _load_balancing_algo(1),
+      _own_group(true),_interpolation_matrix(0), _locator(0),
+      _source_field(0),_own_source_field(false),
+      _target_field(0),_own_target_field(false),
+      _default_field_value(0.0),
+      _comm(MPI_COMM_NULL)
   {
     ParaMEDMEM::CommInterface comm;
     int *ranks_world=new int[procIds.size()]; // ranks of sources and targets in world_comm
@@ -214,10 +223,10 @@ namespace ParaMEDMEM
     comm.commGroup(world_comm,&world_group);
     comm.groupIncl(world_group,procIds.size(),ranks_world,&group);
     delete [] ranks_world;
-    MPI_Comm theComm;
-    comm.commCreate(world_comm,group,&theComm);
+    comm.commCreate(world_comm,group,&_comm);
     comm.groupFree(&group);
-    if(theComm==MPI_COMM_NULL)
+    comm.groupFree(&world_group);
+    if(_comm==MPI_COMM_NULL)
       {
         _group=0;
         return ;
@@ -225,7 +234,7 @@ namespace ParaMEDMEM
     std::set<int> idsUnion;
     for(std::size_t i=0;i<procIds.size();i++)
       idsUnion.insert(i);
-    _group=new MPIProcessorGroup(comm,idsUnion,theComm);
+    _group=new MPIProcessorGroup(comm,idsUnion,_comm);
   }
 
   OverlapDEC::~OverlapDEC()
@@ -237,6 +246,12 @@ namespace ParaMEDMEM
     if(_own_target_field)
       delete _target_field;
     delete _interpolation_matrix;
+    delete _locator;
+    if (_comm != MPI_COMM_NULL)
+      {
+        ParaMEDMEM::CommInterface comm;
+        comm.commFree(&_comm);
+      }
   }
 
   void OverlapDEC::sendRecvData(bool way)
@@ -249,7 +264,7 @@ namespace ParaMEDMEM
 
   void OverlapDEC::sendData()
   {
-    _interpolation_matrix->multiply();
+    _interpolation_matrix->multiply(_default_field_value);
   }
 
   void OverlapDEC::recvData()
@@ -262,24 +277,31 @@ namespace ParaMEDMEM
   {
     if(!isInGroup())
       return ;
+    // Check number of components of field on both side (for now allowing void field/mesh on one proc is not allowed)
+    if (!_source_field || !_source_field->getField())
+      throw INTERP_KERNEL::Exception("OverlapDEC::synchronize(): currently, having a void source field on a proc is not allowed!");
+    if (!_target_field || !_target_field->getField())
+      throw INTERP_KERNEL::Exception("OverlapDEC::synchronize(): currently, having a void target field on a proc is not allowed!");
+    if (_target_field->getField()->getNumberOfComponents() != _source_field->getField()->getNumberOfComponents())
+      throw INTERP_KERNEL::Exception("OverlapDEC::synchronize(): source and target field have different number of components!");
     delete _interpolation_matrix;
-    _interpolation_matrix=new OverlapInterpolationMatrix(_source_field,_target_field,*_group,*this,*this);
-    OverlapElementLocator locator(_source_field,_target_field,*_group);
-    locator.copyOptions(*this);
-    locator.exchangeMeshes(*_interpolation_matrix);
-    std::vector< std::pair<int,int> > jobs=locator.getToDoList();
-    std::string srcMeth=locator.getSourceMethod();
-    std::string trgMeth=locator.getTargetMethod();
+    _locator = new OverlapElementLocator(_source_field,_target_field,*_group, getBoundingBoxAdjustmentAbs(), _load_balancing_algo);
+    _interpolation_matrix=new OverlapInterpolationMatrix(_source_field,_target_field,*_group,*this,*this, *_locator);
+    _locator->copyOptions(*this);
+    _locator->exchangeMeshes(*_interpolation_matrix);
+    std::vector< std::pair<int,int> > jobs=_locator->getToDoList();
+    std::string srcMeth=_locator->getSourceMethod();
+    std::string trgMeth=_locator->getTargetMethod();
     for(std::vector< std::pair<int,int> >::const_iterator it=jobs.begin();it!=jobs.end();it++)
       {
-        const MEDCouplingPointSet *src=locator.getSourceMesh((*it).first);
-        const DataArrayInt *srcIds=locator.getSourceIds((*it).first);
-        const MEDCouplingPointSet *trg=locator.getTargetMesh((*it).second);
-        const DataArrayInt *trgIds=locator.getTargetIds((*it).second);
-        _interpolation_matrix->addContribution(src,srcIds,srcMeth,(*it).first,trg,trgIds,trgMeth,(*it).second);
+        const MEDCouplingPointSet *src=_locator->getSourceMesh((*it).first);
+        const DataArrayInt *srcIds=_locator->getSourceIds((*it).first);
+        const MEDCouplingPointSet *trg=_locator->getTargetMesh((*it).second);
+        const DataArrayInt *trgIds=_locator->getTargetIds((*it).second);
+        _interpolation_matrix->computeLocalIntersection(src,srcIds,srcMeth,(*it).first,trg,trgIds,trgMeth,(*it).second);
       }
-    _interpolation_matrix->prepare(locator.getProcsInInteraction());
-    _interpolation_matrix->computeDeno();
+    _interpolation_matrix->prepare(_locator->getProcsToSendFieldData());
+    _interpolation_matrix->computeSurfacesAndDeno();
   }
 
   void OverlapDEC::attachSourceLocalField(ParaFIELD *field, bool ownPt)
@@ -302,10 +324,39 @@ namespace ParaMEDMEM
     _own_target_field=ownPt;
   }
 
+  void OverlapDEC::attachSourceLocalField(MEDCouplingFieldDouble *field)
+  {
+    if(!isInGroup())
+      return ;
+
+    ParaMESH *paramesh = new ParaMESH(static_cast<MEDCouplingPointSet *>(const_cast<MEDCouplingMesh *>(field->getMesh())),
+                                      *_group,field->getMesh()->getName());
+    ParaFIELD *tmpField=new ParaFIELD(field, paramesh, *_group);
+    tmpField->setOwnSupport(true);
+    attachSourceLocalField(tmpField,true);
+  }
+
+  void OverlapDEC::attachTargetLocalField(MEDCouplingFieldDouble *field)
+  {
+    if(!isInGroup())
+      return ;
+
+    ParaMESH *paramesh = new ParaMESH(static_cast<MEDCouplingPointSet *>(const_cast<MEDCouplingMesh *>(field->getMesh())),
+                                      *_group,field->getMesh()->getName());
+    ParaFIELD *tmpField=new ParaFIELD(field, paramesh, *_group);
+    tmpField->setOwnSupport(true);
+    attachTargetLocalField(tmpField,true);
+  }
+
   bool OverlapDEC::isInGroup() const
   {
     if(!_group)
       return false;
     return _group->containsMyRank();
   }
+
+  void OverlapDEC::debugPrintWorkSharing(std::ostream & ostr) const
+  {
+    _locator->debugPrintWorkSharing(ostr);
+  }
 }