src/ParaMEDMEM/InterpKernelDEC.cxx

   1 // Copyright (C) 2007-2024  CEA, EDF
   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 <mpi.h>
  21 #include "CommInterface.hxx"
  22 #include "Topology.hxx"
  23 #include "BlockTopology.hxx"
  24 #include "ComponentTopology.hxx"
  25 #include "ParaFIELD.hxx"
  26 #include "MPIProcessorGroup.hxx"
  27 #include "ParaMESH.hxx"
  28 #include "DEC.hxx"
  29 #include "InterpolationMatrix.hxx"
  30 #include "InterpKernelDEC.hxx"
  31 #include "ElementLocator.hxx"
  32
  33 namespace MEDCoupling
  34 {
  35   InterpKernelDEC::InterpKernelDEC():
  36     DisjointDEC(),
  37     _interpolation_matrix(0)
  38   {
  39   }
  40
  41   /*!
  42     This constructor creates an InterpKernelDEC which has \a source_group as a working side
  43     and  \a target_group as an idle side. All the processors will actually participate, but intersection computations will be performed on the working side during the \a synchronize() phase.
  44     The constructor must be called synchronously on all processors of both processor groups.
  45     The source group and target group MUST form a partition of all the procs within the communicator passed as 'world_comm'
  46     when building the group.
  47
  48     \param source_group working side ProcessorGroup
  49     \param target_group lazy side ProcessorGroup
  50
  51   */
  52   InterpKernelDEC::InterpKernelDEC(ProcessorGroup& source_group, ProcessorGroup& target_group):
  53     DisjointDEC(source_group, target_group),
  54     _interpolation_matrix(0)
  55   {
  56
  57   }
  58
  59
  60   /*!
  61    * Creates an InterpKernelDEC from a set of source procs IDs and target group IDs.
  62    * The difference with the ctor using groups is that the set of procs might not cover entirely MPI_COMM_WORLD
  63    * (a sub-communicator holding the union of source and target procs is recreated internally).
  64    */
  65   InterpKernelDEC::InterpKernelDEC(const std::set<int>& src_ids, const std::set<int>& trg_ids,
  66                                    const MPI_Comm& world_comm):
  67     DisjointDEC(src_ids,trg_ids,world_comm),
  68     _interpolation_matrix(0)
  69   {
  70   }
  71
  72   InterpKernelDEC::~InterpKernelDEC()
  73   {
  74     release();
  75   }
  76
  77   void InterpKernelDEC::release()
  78   {
  79     if (_interpolation_matrix != nullptr)
  80       delete _interpolation_matrix;
  81     _interpolation_matrix = nullptr;
  82     DisjointDEC::cleanInstance();
  83   }
  84
  85
  86   /*!
  87     \brief Synchronization process for exchanging topologies.
  88
  89     This method prepares all the structures necessary for sending data from a processor group to the other. It uses the mesh
  90     underlying the fields that have been set with attachLocalField method.
  91     It works in four steps :
  92     -# Bounding boxes are computed for each sub-domain,
  93     -# The lazy side mesh parts that are likely to intersect the working side local processor are sent to the working side,
  94     -# The working side calls the interpolation kernel to compute the intersection between local and imported mesh.
  95     -# The lazy side is updated so that it knows the structure of the data that will be sent by
  96     the working side during a \a sendData() call.
  97
  98   */
  99   void InterpKernelDEC::synchronize()
 100   {
 101     if(!isInUnion())
 102       return ;
 103     delete _interpolation_matrix;
 104     _interpolation_matrix = new InterpolationMatrix (_local_field, *_source_group,*_target_group,*this,*this);
 105
 106     //setting up the communication DEC on both sides
 107     if (_source_group->containsMyRank())
 108       {
 109         //locate the distant meshes
 110         ElementLocator locator(*_local_field, *_target_group, *_source_group);
 111         //transferring option from InterpKernelDEC to ElementLocator
 112         locator.copyOptions(*this);
 113         MEDCouplingPointSet* distant_mesh=0;
 114         mcIdType* distant_ids=0;
 115         std::string distantMeth;
 116         for (int i=0; i<_target_group->size(); i++)
 117           {
 118             //        int idistant_proc = (i+_source_group->myRank())%_target_group->size();
 119             int idistant_proc=i;
 120
 121             //gathers pieces of the target meshes that can intersect the local mesh
 122             locator.exchangeMesh(idistant_proc,distant_mesh,distant_ids);
 123             if (distant_mesh !=0)
 124               {
 125                 locator.exchangeMethod(_method,idistant_proc,distantMeth);
 126                 //adds the contribution of the distant mesh on the local one
 127                 int idistant_proc_in_union=_union_group->translateRank(_target_group,idistant_proc);
 128                 //std::cout <<"add contribution from proc "<<idistant_proc_in_union<<" to proc "<<_union_group->myRank()<<std::endl;
 129                 _interpolation_matrix->addContribution(*distant_mesh,idistant_proc_in_union,distant_ids,_method,distantMeth);
 130                 distant_mesh->decrRef();
 131                 delete [] distant_ids;
 132                 distant_mesh=0;
 133                 distant_ids=0;
 134               }
 135           }
 136        _interpolation_matrix->finishContributionW(locator);
 137       }
 138
 139     if (_target_group->containsMyRank())
 140       {
 141         ElementLocator locator(*_local_field, *_source_group, *_target_group);
 142         //transferring option from InterpKernelDEC to ElementLocator
 143         locator.copyOptions(*this);
 144         MEDCouplingPointSet* distant_mesh=0;
 145         mcIdType* distant_ids=0;
 146         for (int i=0; i<_source_group->size(); i++)
 147           {
 148             //        int idistant_proc = (i+_target_group->myRank())%_source_group->size();
 149             int  idistant_proc=i;
 150             //gathers pieces of the target meshes that can intersect the local mesh
 151             locator.exchangeMesh(idistant_proc,distant_mesh,distant_ids);
 152             //std::cout << " Data sent from "<<_union_group->myRank()<<" to source proc "<< idistant_proc<<std::endl;
 153             if (distant_mesh!=0)
 154               {
 155                 std::string distantMeth;
 156                 locator.exchangeMethod(_method,idistant_proc,distantMeth);
 157                 distant_mesh->decrRef();
 158                 delete [] distant_ids;
 159                 distant_mesh=0;
 160                 distant_ids=0;
 161               }
 162           }
 163         _interpolation_matrix->finishContributionL(locator);
 164       }
 165     _interpolation_matrix->prepare();
 166   }
 167
 168   /*!
 169    * Set a default value for non fetched entities
 170    */
 171   void InterpKernelDEC::synchronizeWithDefaultValue(double val)
 172   {
 173     this->synchronize();
 174     if(_interpolation_matrix )
 175       _interpolation_matrix->setDefaultValue(val);
 176   }
 177
 178   /*!
 179     Receives the data whether the processor is on the working side or on the lazy side. It must match a \a sendData() call on the other side.
 180   */
 181   void InterpKernelDEC::recvData()
 182   {
 183     if (_source_group->containsMyRank())
 184       _interpolation_matrix->transposeMultiply(*_local_field->getField());
 185     else if (_target_group->containsMyRank())
 186       {
 187         _interpolation_matrix->multiply(*_local_field->getField());
 188         if (getForcedRenormalization())
 189           renormalizeTargetField(getMeasureAbsStatus());
 190       }
 191   }
 192
 193   MCAuto<DataArrayIdType> InterpKernelDEC::retrieveNonFetchedIds() const
 194   {
 195     if( _source_group->containsMyRank() )
 196     {
 197       return this->retrieveNonFetchedIdsSource();
 198     }
 199     if( _target_group->containsMyRank() )
 200     {
 201       return this->retrieveNonFetchedIdsTarget();
 202     }
 203     THROW_IK_EXCEPTION("Not detected side of rank !");
 204   }
 205
 206   MCAuto<DataArrayIdType> InterpKernelDEC::retrieveNonFetchedIdsSource() const
 207   {
 208     return _interpolation_matrix->retrieveNonFetchedIdsSource();
 209   }
 210
 211   MCAuto<DataArrayIdType> InterpKernelDEC::retrieveNonFetchedIdsTarget() const
 212   {
 213     mcIdType nbTuples = _local_field->getField()->getNumberOfTuplesExpected();
 214     return _interpolation_matrix->retrieveNonFetchedIdsTarget(nbTuples);
 215   }
 216
 217   /*!
 218     Receives the data at time \a time in asynchronous mode. The value of the field
 219     will be time-interpolated from the field values received.
 220     \param time time at which the value is desired
 221   */
 222   void InterpKernelDEC::recvData( double time )
 223   {
 224     _interpolation_matrix->getAccessDEC()->setTime(time);
 225     recvData() ;
 226   }
 227
 228   /*!
 229     Sends the data whether the processor is on the working side or on the lazy side.
 230     It must match a recvData() call on the other side.
 231   */
 232   void InterpKernelDEC::sendData()
 233   {
 234     if (_source_group->containsMyRank())
 235       {
 236
 237         _interpolation_matrix->multiply(*_local_field->getField());
 238         if (getForcedRenormalization())
 239           renormalizeTargetField(getMeasureAbsStatus());
 240
 241       }
 242     else if (_target_group->containsMyRank())
 243       _interpolation_matrix->transposeMultiply(*_local_field->getField());
 244   }
 245
 246   /*!
 247     Sends the data available at time \a time in asynchronous mode.
 248     \param time time at which the value is available
 249     \param deltatime time interval between the value presently sent and the next one.
 250   */
 251   void InterpKernelDEC::sendData( double time , double deltatime )
 252   {
 253     _interpolation_matrix->getAccessDEC()->setTime(time,deltatime);
 254     sendData() ;
 255   }
 256
 257 }