src/MEDCoupling/MEDCouplingCartesianAMRMesh.cxx

   1 // Copyright (C) 2007-2014  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 // Author : Anthony Geay
  20
  21 #include "MEDCouplingCartesianAMRMesh.hxx"
  22 #include "MEDCoupling1GTUMesh.hxx"
  23 #include "MEDCouplingIMesh.hxx"
  24 #include "MEDCouplingUMesh.hxx"
  25
  26 #include <limits>
  27 #include <sstream>
  28 #include <numeric>
  29
  30 using namespace ParaMEDMEM;
  31
  32 /// @cond INTERNAL
  33
  34 /*!
  35  * \param [in] mesh not null pointer of refined mesh replacing the cell range of \a father defined by the bottom left and top right just after.
  36  * \param [in] bottomLeftTopRight a vector equal to the space dimension of \a mesh that specifies for each dimension, the included cell start of the range for the first element of the pair,
  37  *                                a the end cell (\b excluded) of the range for the second element of the pair.
  38  */
  39 MEDCouplingCartesianAMRPatch::MEDCouplingCartesianAMRPatch(MEDCouplingCartesianAMRMesh *mesh, const std::vector< std::pair<int,int> >& bottomLeftTopRight)
  40 {
  41   if(!mesh)
  42     throw INTERP_KERNEL::Exception("EDCouplingCartesianAMRPatch constructor : input mesh is NULL !");
  43   _mesh=mesh; _mesh->incrRef();
  44   int dim((int)bottomLeftTopRight.size()),dimExp(_mesh->getSpaceDimension());
  45   if(dim!=dimExp)
  46     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch constructor : space dimension of father and input bottomLeft/topRight size mismatches !");
  47   _bl_tr=bottomLeftTopRight;
  48 }
  49
  50 int MEDCouplingCartesianAMRPatch::getNumberOfCellsRecursiveWithOverlap() const
  51 {
  52   return _mesh->getNumberOfCellsRecursiveWithOverlap();
  53 }
  54
  55 int MEDCouplingCartesianAMRPatch::getNumberOfCellsRecursiveWithoutOverlap() const
  56 {
  57   return _mesh->getNumberOfCellsRecursiveWithoutOverlap();
  58 }
  59
  60 int MEDCouplingCartesianAMRPatch::getMaxNumberOfLevelsRelativeToThis() const
  61 {
  62   return _mesh->getMaxNumberOfLevelsRelativeToThis();
  63 }
  64
  65 void MEDCouplingCartesianAMRPatch::addPatch(const std::vector< std::pair<int,int> >& bottomLeftTopRight, const std::vector<int>& factors)
  66 {
  67   return _mesh->addPatch(bottomLeftTopRight,factors);
  68 }
  69
  70 int MEDCouplingCartesianAMRPatch::getNumberOfOverlapedCellsForFather() const
  71 {
  72   return MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(_bl_tr);
  73 }
  74
  75 /*!
  76  * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
  77  * the must be >= 0.
  78  *
  79  * \param [in] other - The other patch
  80  * \param [in] ghostLev - The size of the neighborhood zone.
  81  *
  82  * \throw if \a this or \a other are invalid (end before start).
  83  * \throw if \a ghostLev is \b not >= 0 .
  84  * \throw if \a this and \a other have not the same space dimension.
  85  */
  86 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhood(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
  87 {
  88   if(ghostLev<0)
  89     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the size of the neighborhood must be >= 0 !");
  90   if(!other)
  91     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the input patch is NULL !");
  92   const std::vector< std::pair<int,int> >& thisp(getBLTRRange());
  93   const std::vector< std::pair<int,int> >& otherp(other->getBLTRRange());
  94   std::size_t thispsize(thisp.size());
  95   if(thispsize!=otherp.size())
  96     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the dimensions must be the same !");
  97   for(std::size_t i=0;i<thispsize;i++)
  98     {
  99       const std::pair<int,int>& thispp(thisp[i]);
 100       const std::pair<int,int>& otherpp(otherp[i]);
 101       if(thispp.second<thispp.first)
 102         throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
 103       if(otherpp.second<otherpp.first)
 104         throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
 105       if(otherpp.first==thispp.second+ghostLev-1)
 106         continue;
 107       if(otherpp.second+ghostLev-1==thispp.first)
 108         continue;
 109       int start(std::max(thispp.first,otherpp.first)),end(std::min(thispp.second,otherpp.second));
 110       if(end<start)
 111         return false;
 112     }
 113   return true;
 114 }
 115
 116 std::size_t MEDCouplingCartesianAMRPatch::getHeapMemorySizeWithoutChildren() const
 117 {
 118   std::size_t ret(sizeof(MEDCouplingCartesianAMRPatch));
 119   ret+=_bl_tr.capacity()*sizeof(std::pair<int,int>);
 120   return ret;
 121 }
 122
 123 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRPatch::getDirectChildren() const
 124 {
 125   std::vector<const BigMemoryObject *> ret;
 126   if((const MEDCouplingCartesianAMRMesh *)_mesh)
 127     ret.push_back((const MEDCouplingCartesianAMRMesh *)_mesh);
 128   return ret;
 129 }
 130
 131 /// @endcond
 132
 133
 134 MEDCouplingCartesianAMRMesh *MEDCouplingCartesianAMRMesh::New(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
 135                                                               const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop)
 136 {
 137   return new MEDCouplingCartesianAMRMesh(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop);
 138 }
 139
 140 int MEDCouplingCartesianAMRMesh::getSpaceDimension() const
 141 {
 142   return _mesh->getSpaceDimension();
 143 }
 144
 145 void MEDCouplingCartesianAMRMesh::setFactors(const std::vector<int>& newFactors)
 146 {
 147   if(getSpaceDimension()!=(int)newFactors.size())
 148     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::setFactors : size of input factors is not equal to the space dimension !");
 149   if(_factors.empty())
 150     {
 151       _factors=newFactors;
 152       return ;
 153     }
 154   if(_factors==newFactors)
 155     return ;
 156   if(!_patches.empty())
 157     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::setFactors : modification of factors is not allowed when presence of patches !");
 158   _factors=newFactors;
 159 }
 160
 161 int MEDCouplingCartesianAMRMesh::getMaxNumberOfLevelsRelativeToThis() const
 162 {
 163   int ret(1);
 164   for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
 165     ret=std::max(ret,(*it)->getMaxNumberOfLevelsRelativeToThis()+1);
 166   return ret;
 167 }
 168
 169 int MEDCouplingCartesianAMRMesh::getNumberOfCellsAtCurrentLevel() const
 170 {
 171   return _mesh->getNumberOfCells();
 172 }
 173
 174 int MEDCouplingCartesianAMRMesh::getNumberOfCellsRecursiveWithOverlap() const
 175 {
 176   int ret(_mesh->getNumberOfCells());
 177   for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
 178     {
 179       ret+=(*it)->getNumberOfCellsRecursiveWithOverlap();
 180     }
 181   return ret;
 182 }
 183
 184 int MEDCouplingCartesianAMRMesh::getNumberOfCellsRecursiveWithoutOverlap() const
 185 {
 186   int ret(_mesh->getNumberOfCells());
 187   for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
 188     {
 189       ret-=(*it)->getNumberOfOverlapedCellsForFather();
 190       ret+=(*it)->getNumberOfCellsRecursiveWithoutOverlap();
 191     }
 192   return ret;
 193 }
 194
 195 const MEDCouplingCartesianAMRMesh *MEDCouplingCartesianAMRMesh::getFather() const
 196 {
 197   return _father;
 198 }
 199
 200 const MEDCouplingCartesianAMRMesh *MEDCouplingCartesianAMRMesh::getGodFather() const
 201 {
 202   if(_father==0)
 203     return this;
 204   else
 205     return _father->getGodFather();
 206 }
 207
 208 void MEDCouplingCartesianAMRMesh::detachFromFather()
 209 {
 210   _father=0;
 211 }
 212
 213 /*!
 214  * \param [in] bottomLeftTopRight a vector equal to the space dimension of \a mesh that specifies for each dimension, the included cell start of the range for the first element of the pair,
 215  *                                a the end cell (\b excluded) of the range for the second element of the pair.
 216  * \param [in] factors The factor of refinement per axis (different from 0).
 217  */
 218 void MEDCouplingCartesianAMRMesh::addPatch(const std::vector< std::pair<int,int> >& bottomLeftTopRight, const std::vector<int>& factors)
 219 {
 220   checkFactorsAndIfNotSetAssign(factors);
 221   MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> mesh(static_cast<MEDCouplingIMesh *>(_mesh->buildStructuredSubPart(bottomLeftTopRight)));
 222   mesh->refineWithFactor(factors);
 223   MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRMesh> zeMesh(new MEDCouplingCartesianAMRMesh(this,mesh));
 224   MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> elt(new MEDCouplingCartesianAMRPatch(zeMesh,bottomLeftTopRight));
 225   _patches.push_back(elt);
 226 }
 227
 228 /// @cond INTERNAL
 229
 230 class InternalPatch : public RefCountObjectOnly
 231 {
 232 public:
 233   InternalPatch():_nb_of_true(0) { }
 234   int getDimension() const { return (int)_part.size(); }
 235   double getEfficiency() const { return (double)_nb_of_true/(double)_crit.size(); }
 236   int getNumberOfCells() const { return (int)_crit.size(); }
 237   void setNumberOfTrue(int nboft) { _nb_of_true=nboft; }
 238   std::vector<bool>& getCriterion() { return _crit; }
 239   const std::vector<bool>& getConstCriterion() const { return _crit; }
 240   void setPart(const std::vector< std::pair<int,int> >& part) { _part=part; }
 241   std::vector< std::pair<int,int> >& getPart() { return _part; }
 242   const std::vector< std::pair<int,int> >& getConstPart() const { return _part; }
 243   bool presenceOfTrue() const { return _nb_of_true>0; }
 244   std::vector<int> computeCGS() const { return MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(_part); }
 245   std::vector< std::vector<int> > computeSignature() const { return MEDCouplingStructuredMesh::ComputeSignaturePerAxisOf(computeCGS(),getConstCriterion()); }
 246   double getEfficiencyPerAxis(int axisId) const { return (double)_nb_of_true/((double)(_part[axisId].second-_part[axisId].first)); }
 247   void zipToFitOnCriterion();
 248   void updateNumberOfTrue() const;
 249   MEDCouplingAutoRefCountObjectPtr<InternalPatch> extractPart(const std::vector< std::pair<int,int> >&partInGlobal) const;
 250   MEDCouplingAutoRefCountObjectPtr<InternalPatch> deepCpy() const;
 251 protected:
 252   ~InternalPatch() { }
 253 private:
 254   mutable int _nb_of_true;
 255   std::vector<bool> _crit;
 256   //! _part is global
 257   std::vector< std::pair<int,int> > _part;
 258 };
 259
 260 void InternalPatch::zipToFitOnCriterion()
 261 {
 262   std::vector<int> cgs(computeCGS());
 263   std::vector<bool> newCrit;
 264   std::vector< std::pair<int,int> > newPart,newPart2;
 265   int newNbOfTrue(MEDCouplingStructuredMesh::FindMinimalPartOf(cgs,_crit,newCrit,newPart));
 266   MEDCouplingStructuredMesh::ChangeReferenceToGlobalOfCompactFrmt(_part,newPart,newPart2);
 267   if(newNbOfTrue!=_nb_of_true)
 268     throw INTERP_KERNEL::Exception("InternalPatch::zipToFitOnCrit : internal error !");
 269   _crit=newCrit; _part=newPart2;
 270 }
 271
 272 void InternalPatch::updateNumberOfTrue() const
 273 {
 274   _nb_of_true=(int)std::count(_crit.begin(),_crit.end(),true);
 275 }
 276
 277 MEDCouplingAutoRefCountObjectPtr<InternalPatch> InternalPatch::extractPart(const std::vector< std::pair<int,int> >&partInGlobal) const
 278 {
 279   MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(new InternalPatch);
 280   std::vector<int> cgs(computeCGS());
 281   std::vector< std::pair<int,int> > newPart;
 282   MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(_part,partInGlobal,newPart);
 283   MEDCouplingStructuredMesh::ExtractFieldOfBoolFrom(cgs,_crit,newPart,ret->getCriterion());
 284   ret->setPart(partInGlobal);
 285   ret->updateNumberOfTrue();
 286   return ret;
 287 }
 288
 289 MEDCouplingAutoRefCountObjectPtr<InternalPatch> InternalPatch::deepCpy() const
 290 {
 291   MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(new InternalPatch);
 292   (*ret)=*this;
 293   return ret;
 294 }
 295
 296 void DissectBigPatch(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int rangeOfAxisId, bool& cutFound, int& cutPlace)
 297 {
 298   cutFound=false; cutPlace=-1;
 299   std::vector<double> ratio(rangeOfAxisId-1);
 300   for(int id=0;id<rangeOfAxisId-1;id++)
 301     {
 302       double efficiency[2];
 303       for(int h=0;h<2;h++)
 304         {
 305           std::vector< std::pair<int,int> > rectH(patchToBeSplit->getConstPart());
 306           if(h==0)
 307             rectH[axisId].second=patchToBeSplit->getConstPart()[axisId].first+id;
 308           else
 309             rectH[axisId].first=patchToBeSplit->getConstPart()[axisId].first+id;
 310           MEDCouplingAutoRefCountObjectPtr<InternalPatch> p(patchToBeSplit->deepCpy());
 311           p->zipToFitOnCriterion();
 312           //anouar rectH ?
 313           efficiency[h]=p->getEfficiencyPerAxis(axisId);
 314         }
 315       ratio[id]=std::max(efficiency[0],efficiency[1])/std::min(efficiency[0],efficiency[1]);
 316     }
 317   int minCellDirection(bso.getMinCellDirection()),indexMin(-1);
 318   int dimRatioInner(rangeOfAxisId-1-2*(minCellDirection-1));
 319   std::vector<double> ratio_inner(dimRatioInner);
 320   double minRatio(1.e10);
 321   for(int i=0; i<dimRatioInner; i++)
 322     {
 323       if(ratio[minCellDirection-1+i]<minRatio)
 324         {
 325           minRatio=ratio[minCellDirection-1+i];
 326           indexMin=i+minCellDirection;
 327         }
 328     }
 329   cutFound=true; cutPlace=indexMin+patchToBeSplit->getConstPart()[axisId].first-1;
 330 }
 331
 332 void FindHole(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int& axisId, bool& cutFound, int& cutPlace)
 333 {
 334   cutPlace=-1; cutFound=false;
 335   int minCellDirection(bso.getMinCellDirection());
 336   const int dim(patchToBeSplit->getDimension());
 337   std::vector< std::vector<int> > signatures(patchToBeSplit->computeSignature());
 338   for(int id=0;id<dim;id++)
 339     {
 340       const std::vector<int>& signature(signatures[id]);
 341       std::vector<int> hole;
 342       std::vector<double> distance;
 343       int len((int)signature.size());
 344       for(int i=0;i<len;i++)
 345         if(signature[i]==0)
 346           if(len>= 2*minCellDirection && i >= minCellDirection-1 && i <= len-minCellDirection-1)
 347             hole.push_back(i);
 348       if(!hole.empty())
 349         {
 350           double center(((double)len/2.));
 351           for(std::size_t i=0;i<hole.size();i++)
 352             distance.push_back(fabs(hole[i]+1.-center));
 353
 354           std::size_t posDistanceMin(std::distance(distance.begin(),std::min_element(distance.begin(),distance.end())));
 355           cutFound=true;
 356           axisId=id;
 357           cutPlace=hole[posDistanceMin]+patchToBeSplit->getConstPart()[axisId].first+1;
 358           return ;
 359         }
 360     }
 361 }
 362
 363 void FindInflection(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, bool& cutFound, int& cutPlace, int& axisId)
 364 {
 365   cutFound=false; cutPlace=-1;// do not set axisId before to be sure that cutFound was set to true
 366
 367   const std::vector< std::pair<int,int> >& part(patchToBeSplit->getConstPart());
 368   int sign,minCellDirection(bso.getMinCellDirection());
 369   const int dim(patchToBeSplit->getDimension());
 370
 371   std::vector<int> zeroCrossDims(dim,-1);
 372   std::vector<int> zeroCrossVals(dim,-1);
 373   std::vector< std::vector<int> > signatures(patchToBeSplit->computeSignature());
 374   for (int id=0;id<dim;id++)
 375     {
 376       const std::vector<int>& signature(signatures[id]);
 377
 378       std::vector<int> derivate_second_order,gradient_absolute,signe_change,zero_cross,edge,max_cross_list ;
 379       std::vector<double> distance ;
 380
 381       for (unsigned int i=1;i<signature.size()-1;i++)
 382         derivate_second_order.push_back(signature[i-1]-2*signature[i]+signature[i+1]) ;
 383
 384       // Gradient absolute value
 385       for ( unsigned int i=1;i<derivate_second_order.size();i++)
 386         gradient_absolute.push_back(fabs(derivate_second_order[i]-derivate_second_order[i-1])) ;
 387       if(derivate_second_order.empty())
 388         continue;
 389       for (unsigned int i=0;i<derivate_second_order.size()-1;i++)
 390         {
 391           if (derivate_second_order[i]*derivate_second_order[i+1] < 0 )
 392             sign = -1 ;
 393           if (derivate_second_order[i]*derivate_second_order[i+1] > 0 )
 394             sign = 1 ;
 395           if (derivate_second_order[i]*derivate_second_order[i+1] == 0 )
 396             sign = 0 ;
 397           if ( sign==0 || sign==-1 )
 398             if ( i >= (unsigned int)minCellDirection-2 && i <= signature.size()-minCellDirection-2 )
 399               {
 400                 zero_cross.push_back(i) ;
 401                 edge.push_back(gradient_absolute[i]) ;
 402               }
 403           signe_change.push_back(sign) ;
 404         }
 405       if ( zero_cross.size() > 0 )
 406         {
 407           int max_cross=*max_element(edge.begin(),edge.end()) ;
 408           for (unsigned int i=0;i<edge.size();i++)
 409             if (edge[i]==max_cross)
 410               max_cross_list.push_back(zero_cross[i]+1) ;
 411
 412           double center((signature.size()/2.0));
 413           for (unsigned int i=0;i<max_cross_list.size();i++)
 414             distance.push_back(fabs(max_cross_list[i]+1-center));
 415
 416           float distance_min=*min_element(distance.begin(),distance.end()) ;
 417           int pos_distance_min=find(distance.begin(),distance.end(),distance_min)-distance.begin();
 418           int best_place = max_cross_list[pos_distance_min] + part[id].first ;
 419           if ( max_cross >=0 )
 420             {
 421               zeroCrossDims[id] = best_place ;
 422               zeroCrossVals[id] = max_cross ;
 423             }
 424         }
 425       derivate_second_order.clear() ;
 426       gradient_absolute.clear() ;
 427       signe_change.clear() ;
 428       zero_cross.clear() ;
 429       edge.clear() ;
 430       max_cross_list.clear() ;
 431       distance.clear() ;
 432     }
 433
 434   if ( zeroCrossDims[0]!=-1 || zeroCrossDims[1]!=-1  )
 435     {
 436       int max_cross_dims = *max_element(zeroCrossVals.begin(),zeroCrossVals.end()) ;
 437
 438       if (zeroCrossVals[0]==max_cross_dims &&  zeroCrossVals[1]==max_cross_dims )
 439         {
 440           int nl_left(part[0].second-part[0].first);
 441           int nc_left(part[1].second-part[1].first);
 442           if ( nl_left >=  nc_left )
 443             max_cross_dims = 0 ;
 444           else
 445             max_cross_dims = 1 ;
 446         }
 447       else
 448         max_cross_dims=std::find(zeroCrossVals.begin(),zeroCrossVals.end(),max_cross_dims)-zeroCrossVals.begin();
 449       cutFound=true;
 450       cutPlace=zeroCrossDims[max_cross_dims];
 451       axisId=max_cross_dims ;
 452     }
 453 }
 454
 455 void TryAction4(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int rangeOfAxisId, bool& cutFound, int& cutPlace)
 456 {
 457   cutFound=false;
 458   if(patchToBeSplit->getEfficiency()<=bso.getEffeciencySnd())
 459     {
 460       if(rangeOfAxisId>=2*bso.getMinCellDirection())
 461         {
 462           cutFound=true;
 463           cutPlace=rangeOfAxisId/2+patchToBeSplit->getConstPart()[axisId].first-1;
 464         }
 465     }
 466   else
 467     {
 468       if(patchToBeSplit->getNumberOfCells()>bso.getMaxCells())
 469         {
 470           DissectBigPatch(bso,patchToBeSplit,axisId,rangeOfAxisId,cutFound,cutPlace);
 471         }
 472     }
 473 }
 474
 475 MEDCouplingAutoRefCountObjectPtr<InternalPatch> DealWithNoCut(const InternalPatch *patch)
 476 {
 477   MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(const_cast<InternalPatch *>(patch));
 478   ret->incrRef();
 479   return ret;
 480 }
 481
 482 void DealWithCut(const InternalPatch *patchToBeSplit, int axisId, int cutPlace, std::vector<MEDCouplingAutoRefCountObjectPtr<InternalPatch> >& listOfPatches)
 483 {
 484   MEDCouplingAutoRefCountObjectPtr<InternalPatch> leftPart,rightPart;
 485   std::vector< std::pair<int,int> > rect(patchToBeSplit->getConstPart());
 486   std::vector< std::pair<int,int> > leftRect(rect),rightRect(rect);
 487   leftRect[axisId].second=cutPlace+1;
 488   rightRect[axisId].first=cutPlace+1;
 489   leftPart=patchToBeSplit->extractPart(leftRect);
 490   rightPart=patchToBeSplit->extractPart(rightRect);
 491   leftPart->zipToFitOnCriterion(); rightPart->zipToFitOnCriterion();
 492   listOfPatches.push_back(leftPart);
 493   listOfPatches.push_back(rightPart);
 494 }
 495
 496 /// @endcond
 497
 498 /*!
 499  * This method creates patches in \a this (by destroying the patches if any). This method uses \a criterion array as a field on cells on this level.
 500  */
 501 void MEDCouplingCartesianAMRMesh::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const std::vector<bool>& criterion, const std::vector<int>& factors)
 502 {
 503   int nbCells(getNumberOfCellsAtCurrentLevel());
 504   if(nbCells!=(int)criterion.size())
 505     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createPatchesFromCriterion : the number of tuples of criterion array must be equal to the number of cells at the current level !");
 506   _patches.clear();
 507   std::vector<int> cgs(_mesh->getCellGridStructure());
 508   std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> > listOfPatches,listOfPatchesOK;
 509   //
 510   MEDCouplingAutoRefCountObjectPtr<InternalPatch> p(new InternalPatch);
 511   p->setNumberOfTrue(MEDCouplingStructuredMesh::FindMinimalPartOf(cgs,criterion,p->getCriterion(),p->getPart()));
 512   if(p->presenceOfTrue())
 513     listOfPatches.push_back(p);
 514   while(!listOfPatches.empty())
 515     {
 516       std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> > listOfPatchesTmp;
 517       for(std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> >::iterator it=listOfPatches.begin();it!=listOfPatches.end();it++)
 518         {
 519           //
 520           int axisId,rangeOfAxisId,cutPlace;
 521           bool cutFound;
 522           MEDCouplingStructuredMesh::FindTheWidestAxisOfGivenRangeInCompactFrmt((*it)->getConstPart(),axisId,rangeOfAxisId);
 523           if((*it)->getEfficiency()>=bso.getEffeciency() && (*it)->getNumberOfCells()<bso.getMaxCells())
 524             { listOfPatchesOK.push_back(DealWithNoCut(*it)); continue; }//action 1
 525           FindHole(bso,*it,axisId,cutFound,cutPlace);//axisId overwritten here if FindHole equal to true !
 526           if(cutFound)
 527             { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 2
 528           FindInflection(bso,*it,cutFound,cutPlace,axisId);//axisId overwritten here if cutFound equal to true !
 529           if(cutFound)
 530             { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 3
 531           TryAction4(bso,*it,axisId,rangeOfAxisId,cutFound,cutPlace);
 532           if(cutFound)
 533             { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 4
 534           listOfPatchesOK.push_back(DealWithNoCut(*it));
 535         }
 536       listOfPatches=listOfPatchesTmp;
 537     }
 538   for(std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> >::const_iterator it=listOfPatchesOK.begin();it!=listOfPatchesOK.end();it++)
 539     addPatch((*it)->getConstPart(),factors);
 540 }
 541
 542 /*!
 543  * This method creates patches in \a this (by destroying the patches if any). This method uses \a criterion array as a field on cells on this level.
 544  */
 545 void MEDCouplingCartesianAMRMesh::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const DataArrayByte *criterion, const std::vector<int>& factors)
 546 {
 547   if(!criterion || !criterion->isAllocated())
 548     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createPatchesFromCriterion : the criterion DataArrayByte instance must be allocated and not NULL !");
 549   std::vector<bool> crit(criterion->toVectorOfBool());//check that criterion has one component.
 550   createPatchesFromCriterion(bso,crit,factors);
 551 }
 552
 553 void MEDCouplingCartesianAMRMesh::removeAllPatches()
 554 {
 555   _patches.clear();
 556   declareAsNew();
 557 }
 558
 559 void MEDCouplingCartesianAMRMesh::removePatch(int patchId)
 560 {
 561   checkPatchId(patchId);
 562   int sz((int)_patches.size()),j(0);
 563   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> > patches(sz-1);
 564   for(int i=0;i<sz;i++)
 565     if(i!=patchId)
 566       patches[j++]=_patches[i];
 567   (const_cast<MEDCouplingCartesianAMRMesh *>(_patches[patchId]->getMesh()))->detachFromFather();
 568   _patches=patches;
 569   declareAsNew();
 570 }
 571
 572 int MEDCouplingCartesianAMRMesh::getNumberOfPatches() const
 573 {
 574   return (int)_patches.size();
 575 }
 576
 577 const MEDCouplingCartesianAMRPatch *MEDCouplingCartesianAMRMesh::getPatch(int patchId) const
 578 {
 579   checkPatchId(patchId);
 580   return _patches[patchId];
 581 }
 582
 583 /*!
 584  * This method states if patch2 (with id \a patchId2) is in the neighborhood of patch1 (with id \a patchId1).
 585  * The neighborhood size is defined by \a ghostLev in the reference of \a this ( \b not in the reference of patches !).
 586  */
 587 bool MEDCouplingCartesianAMRMesh::isPatchInNeighborhoodOf(int patchId1, int patchId2, int ghostLev) const
 588 {
 589   if(ghostLev<0)
 590     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::isPatchInNeighborhoodOf : the ghost size must be >=0 !");
 591   const MEDCouplingCartesianAMRPatch *p1(getPatch(patchId1)),*p2(getPatch(patchId2));
 592   if(_factors.empty())
 593     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::isPatchInNeighborhoodOf : no factors defined !");
 594   int ghostLevInPatchRef;
 595   if(ghostLev==0)
 596     ghostLevInPatchRef=0;
 597   else
 598     {
 599       ghostLevInPatchRef=(ghostLev-1)/_factors[0]+1;
 600       for(std::size_t i=0;i<_factors.size();i++)
 601         ghostLevInPatchRef=std::max(ghostLevInPatchRef,(ghostLev-1)/_factors[i]+1);
 602     }
 603   return p1->isInMyNeighborhood(p2,ghostLevInPatchRef);
 604 }
 605
 606 /*!
 607  * This method creates a new cell field array on given \a patchId patch in \a this starting from a coarse cell field on \a this \a cellFieldOnThis.
 608  * This method can be seen as a fast projection from the cell field \a cellFieldOnThis on \c this->getImageMesh() to a refined part of \a this
 609  * defined by the patch with id \a patchId.
 610  *
 611  * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
 612  * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
 613  * \return DataArrayDouble * - The array of the cell field on the requested patch
 614  *
 615  * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
 616  * \throw if \a cellFieldOnThis is NULL or not allocated
 617  * \sa fillCellFieldOnPatch, MEDCouplingIMesh::SpreadCoarseToFine
 618  */
 619 DataArrayDouble *MEDCouplingCartesianAMRMesh::createCellFieldOnPatch(int patchId, const DataArrayDouble *cellFieldOnThis) const
 620 {
 621   if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
 622     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatch : the input cell field array is NULL or not allocated !");
 623   const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
 624   const MEDCouplingIMesh *fine(patch->getMesh()->getImageMesh());
 625   MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(fine->getNumberOfCells(),cellFieldOnThis->getNumberOfComponents());
 626   ret->copyStringInfoFrom(*cellFieldOnThis);
 627   MEDCouplingIMesh::SpreadCoarseToFine(cellFieldOnThis,_mesh->getCellGridStructure(),ret,patch->getBLTRRange(),getFactors());
 628   return ret.retn();
 629 }
 630
 631 /*!
 632  * This method is equivalent to MEDCouplingCartesianAMRMesh::createCellFieldOnPatch except that here instead of creating a new instance
 633  * it fills the value into the \a cellFieldOnPatch data.
 634  *
 635  * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
 636  * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
 637  * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
 638  *
 639  * \sa createCellFieldOnPatch, fillCellFieldComingFromPatch
 640  */
 641 void MEDCouplingCartesianAMRMesh::fillCellFieldOnPatch(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch) const
 642 {
 643   if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
 644     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatch : the input cell field array is NULL or not allocated !");
 645   const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
 646   MEDCouplingIMesh::SpreadCoarseToFine(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors());
 647 }
 648
 649 /*!
 650  * This method is the generalization of fillCellFieldOnPatch method. This method only projects coarse to fine without considering the
 651  * potential neighbor patches covered by the ghost cells of patch with id \a patchId.
 652  *
 653  * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
 654  * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
 655  * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
 656  * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
 657  *
 658  * \sa fillCellFieldOnPatch, fillCellFieldOnPatchGhostAdv
 659  */
 660 void MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchGhost(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch, int ghostLev) const
 661 {
 662   if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
 663     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatchGhost : the input cell field array is NULL or not allocated !");
 664   const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
 665   MEDCouplingIMesh::SpreadCoarseToFineGhost(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),ghostLev);
 666 }
 667
 668 /*!
 669  * This method is a refinement of fillCellFieldOnPatchGhost. fillCellFieldOnPatchGhost is first called.
 670  * Then for all other patches than those pointed by \a patchId that overlap the ghost zone of the patch impact the ghost zone adequately.
 671  *
 672  * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
 673  * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
 674  * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
 675  * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
 676  * \param [in] arrsOnPatches - \b WARNING arrsOnPatches[patchId] is \b NOT \b const. All others are const.
 677  */
 678 void MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchGhostAdv(int patchId, const DataArrayDouble *cellFieldOnThis, int ghostLev, const std::vector<const DataArrayDouble *>& arrsOnPatches) const
 679 {
 680   int nbp(getNumberOfPatches()),dim(getSpaceDimension());
 681   if(nbp!=(int)arrsOnPatches.size())
 682     {
 683       std::ostringstream oss; oss << "MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchGhostAdv : there are " << nbp << " patches in this and " << arrsOnPatches.size() << " arrays in the last parameter !";
 684       throw INTERP_KERNEL::Exception(oss.str().c_str());
 685     }
 686   DataArrayDouble *theFieldToFill(const_cast<DataArrayDouble *>(arrsOnPatches[patchId]));
 687   // first, do as usual
 688   fillCellFieldOnPatchGhost(patchId,cellFieldOnThis,theFieldToFill,ghostLev);
 689   // all reference patch stuff
 690   const MEDCouplingCartesianAMRPatch *refP(getPatch(patchId));
 691   const std::vector< std::pair<int,int> >& refBLTR(refP->getBLTRRange());//[(1,4),(2,4)]
 692   std::vector<int> dimsCoarse(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(refBLTR));//[3,2]
 693   std::transform(dimsCoarse.begin(),dimsCoarse.end(),_factors.begin(),dimsCoarse.begin(),std::multiplies<int>());//[12,8]
 694   std::transform(dimsCoarse.begin(),dimsCoarse.end(),dimsCoarse.begin(),std::bind2nd(std::plus<int>(),2*ghostLev));//[14,10]
 695   std::vector< std::pair<int,int> > rangeCoarse(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(dimsCoarse));//[(0,14),(0,10)]
 696   std::vector<int> fakeFactors(dim,1);
 697   //
 698   for(int i=0;i<nbp;i++)
 699     {
 700       if(i!=patchId)
 701         if(isPatchInNeighborhoodOf(i,patchId,ghostLev))
 702           {
 703             const MEDCouplingCartesianAMRPatch *otherP(getPatch(i));
 704             const std::vector< std::pair<int,int> >& otherBLTR(otherP->getBLTRRange());//[(4,5),(3,4)]
 705             std::vector< std::pair<int,int> > tmp0,tmp1,tmp2;
 706             MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(refBLTR,otherBLTR,tmp0,false);//tmp0=[(3,4),(1,2)]
 707             ApplyFactorsOnCompactFrmt(tmp0,_factors);//tmp0=[(12,16),(4,8)]
 708             ApplyGhostOnCompactFrmt(tmp0,ghostLev);//tmp0=[(13,17),(5,9)]
 709             std::vector< std::pair<int,int> > interstRange(MEDCouplingStructuredMesh::IntersectRanges(tmp0,rangeCoarse));//interstRange=[(13,14),(5,9)]
 710             MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(otherBLTR,refBLTR,tmp1,false);//tmp1=[(-3,0),(-1,1)]
 711             ApplyFactorsOnCompactFrmt(tmp1,_factors);//tmp1=[(-12,-4),(-4,0)]
 712             MEDCouplingStructuredMesh::ChangeReferenceToGlobalOfCompactFrmt(tmp1,interstRange,tmp2,false);//tmp2=[(1,2),(1,5)]
 713             //
 714             std::vector< std::pair<int,int> > dimsFine(otherBLTR);
 715             ApplyFactorsOnCompactFrmt(dimsFine,_factors);
 716             ApplyAllGhostOnCompactFrmt(dimsFine,ghostLev);
 717             //
 718             MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ghostVals(MEDCouplingStructuredMesh::ExtractFieldOfDoubleFrom(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(dimsFine),arrsOnPatches[i],tmp2));
 719             MEDCouplingIMesh::CondenseFineToCoarse(dimsCoarse,ghostVals,interstRange,fakeFactors,theFieldToFill);
 720           }
 721     }
 722 }
 723
 724 /*!
 725  * This method updates \a cellFieldOnThis part of values coming from the cell field \a cellFieldOnPatch lying on patch having id \a patchId.
 726  *
 727  * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
 728  * \param [in] cellFieldOnPatch - The array of the cell field on patch with id \a patchId.
 729  * \param [in,out] cellFieldOnThis The array of the cell field on \a this to be updated only on the part concerning the patch with id \a patchId.
 730  *
 731  * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
 732  * \throw if \a cellFieldOnPatch is NULL or not allocated
 733  * \sa createCellFieldOnPatch, MEDCouplingIMesh::CondenseFineToCoarse,fillCellFieldComingFromPatchGhost
 734  */
 735 void MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch(int patchId, const DataArrayDouble *cellFieldOnPatch, DataArrayDouble *cellFieldOnThis) const
 736 {
 737   if(!cellFieldOnPatch || !cellFieldOnPatch->isAllocated())
 738       throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch : the input cell field array is NULL or not allocated !");
 739   const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
 740   MEDCouplingIMesh::CondenseFineToCoarse(_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),cellFieldOnThis);
 741 }
 742
 743 /*!
 744  * This method is the extension of MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch managing the ghost cells. If this
 745  * method is called with \a ghostLev equal to 0 it behaves exactly as MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch.
 746  *
 747  * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
 748  * \param [in] cellFieldOnPatch - The array of the cell field on patch with id \a patchId.
 749  * \param [in,out] cellFieldOnThis The array of the cell field on \a this to be updated only on the part concerning the patch with id \a patchId.
 750  * \param [in] ghostLev The size of ghost zone (must be >= 0 !)
 751  *
 752  * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
 753  * \throw if \a cellFieldOnPatch is NULL or not allocated
 754  * \sa fillCellFieldComingFromPatch
 755  */
 756 void MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatchGhost(int patchId, const DataArrayDouble *cellFieldOnPatch, DataArrayDouble *cellFieldOnThis, int ghostLev) const
 757 {
 758   if(!cellFieldOnPatch || !cellFieldOnPatch->isAllocated())
 759     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatchGhost : the input cell field array is NULL or not allocated !");
 760   const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
 761   MEDCouplingIMesh::CondenseFineToCoarseGhost(_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),cellFieldOnThis,ghostLev);
 762 }
 763
 764 /*!
 765  * This method finds all patches (located by their ids) that are in the neighborhood of patch with id \a patchId. The neighborhood size is
 766  * defined by ghostLev.
 767  *
 768  * \param [in] patchId - the id of the considered patch.
 769  * \param [in] ghostLev - the size of the neighborhood.
 770  * \return DataArrayInt * - the newly allocated array containing the list of patches in the neighborhood of the considered patch. This array is to be deallocated by the caller.
 771  */
 772 DataArrayInt *MEDCouplingCartesianAMRMesh::findPatchesInTheNeighborhoodOf(int patchId, int ghostLev) const
 773 {
 774   int nbp(getNumberOfPatches());
 775   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(0,1);
 776   for(int i=0;i<nbp;i++)
 777     {
 778       if(i!=patchId)
 779         if(isPatchInNeighborhoodOf(i,patchId,ghostLev))
 780           ret->pushBackSilent(i);
 781     }
 782   return ret.retn();
 783 }
 784
 785 MEDCouplingUMesh *MEDCouplingCartesianAMRMesh::buildUnstructured() const
 786 {
 787   MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> part(_mesh->buildUnstructured());
 788   std::vector<bool> bs(_mesh->getNumberOfCells(),false);
 789   std::vector<int> cgs(_mesh->getCellGridStructure());
 790   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> > msSafe(_patches.size()+1);
 791   std::size_t ii(0);
 792   for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ii++)
 793     {
 794       MEDCouplingStructuredMesh::SwitchOnIdsFrom(cgs,(*it)->getBLTRRange(),bs);
 795       msSafe[ii+1]=(*it)->getMesh()->buildUnstructured();
 796     }
 797   MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsOff(DataArrayInt::BuildListOfSwitchedOff(bs));
 798   msSafe[0]=static_cast<MEDCouplingUMesh *>(part->buildPartOfMySelf(eltsOff->begin(),eltsOff->end(),false));
 799   std::vector< const MEDCouplingUMesh * > ms(msSafe.size());
 800   for(std::size_t i=0;i<msSafe.size();i++)
 801     ms[i]=msSafe[i];
 802   return MEDCouplingUMesh::MergeUMeshes(ms);
 803 }
 804
 805 /*!
 806  * This method returns a mesh containing as cells that there is patches at the current level.
 807  * The patches are seen like 'boxes' that is too say the refinement will not appear here.
 808  *
 809  * \return MEDCoupling1SGTUMesh * - A new object to be managed by the caller containing as cells as there are patches in \a this.
 810  */
 811 MEDCoupling1SGTUMesh *MEDCouplingCartesianAMRMesh::buildMeshFromPatchEnvelop() const
 812 {
 813   std::vector<const MEDCoupling1SGTUMesh *> cells;
 814   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > cellsSafe;
 815   for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
 816     {
 817       const MEDCouplingCartesianAMRPatch *patch(*it);
 818       if(patch)
 819         {
 820           MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> cell(patch->getMesh()->getImageMesh()->asSingleCell());
 821           MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> cell1SGT(cell->build1SGTUnstructured());
 822           cellsSafe.push_back(cell1SGT); cells.push_back(cell1SGT);
 823         }
 824     }
 825   return MEDCoupling1SGTUMesh::Merge1SGTUMeshes(cells);
 826 }
 827
 828 MEDCoupling1SGTUMesh *MEDCouplingCartesianAMRMesh::buildMeshOfDirectChildrenOnly() const
 829 {
 830   std::vector<const MEDCoupling1SGTUMesh *> patches;
 831   std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > patchesSafe;
 832   for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
 833       {
 834         const MEDCouplingCartesianAMRPatch *patch(*it);
 835         if(patch)
 836           {
 837             MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> patchMesh(patch->getMesh()->getImageMesh()->build1SGTUnstructured());
 838             patchesSafe.push_back(patchMesh); patches.push_back(patchMesh);
 839           }
 840       }
 841     return MEDCoupling1SGTUMesh::Merge1SGTUMeshes(patches);
 842 }
 843
 844 MEDCouplingCartesianAMRMesh::MEDCouplingCartesianAMRMesh(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
 845                                                          const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop):_father(0)
 846 {
 847   _mesh=MEDCouplingIMesh::New(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop);
 848 }
 849
 850 MEDCouplingCartesianAMRMesh::MEDCouplingCartesianAMRMesh(MEDCouplingCartesianAMRMesh *father, MEDCouplingIMesh *mesh):_father(father)
 851 {
 852   if(!_father)
 853     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh(MEDCouplingIMesh *mesh) constructor : empty father !");
 854   if(!mesh)
 855     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh(MEDCouplingIMesh *mesh) constructor : The input mesh is null !");
 856   mesh->checkCoherency();
 857   _mesh=mesh; _mesh->incrRef();
 858 }
 859
 860 void MEDCouplingCartesianAMRMesh::checkPatchId(int patchId) const
 861 {
 862   int sz(getNumberOfPatches());
 863   if(patchId<0 || patchId>=sz)
 864     {
 865       std::ostringstream oss; oss << "MEDCouplingCartesianAMRMesh::checkPatchId : invalid patchId (" << patchId << ") ! Must be in [0," << sz << ") !";
 866       throw INTERP_KERNEL::Exception(oss.str().c_str());
 867     }
 868 }
 869
 870 void MEDCouplingCartesianAMRMesh::checkFactorsAndIfNotSetAssign(const std::vector<int>& factors)
 871 {
 872   if(getSpaceDimension()!=(int)factors.size())
 873     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::checkFactorsAndIfNotSetAssign : invalid size of factors ! size must be equal to the spaceDimension !");
 874   if(_factors.empty())
 875     {
 876       _factors=factors;
 877     }
 878   else
 879     {
 880       if(_factors!=factors)
 881         throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::checkFactorsAndIfNotSetAssign : the factors ");
 882     }
 883 }
 884
 885 /*!
 886  * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in refined reference.
 887  * \param [in] factors - the factors per axis.
 888  */
 889 void MEDCouplingCartesianAMRMesh::ApplyFactorsOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, const std::vector<int>& factors)
 890 {
 891   std::size_t sz(factors.size());
 892   if(sz!=partBeforeFact.size())
 893     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::ApplyFactorsOnCompactFrmt : size of input vectors must be the same !");
 894   for(std::size_t i=0;i<sz;i++)
 895     {
 896       partBeforeFact[i].first*=factors[i];
 897       partBeforeFact[i].second*=factors[i];
 898     }
 899 }
 900
 901 /*!
 902  * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in ghost reference.
 903  * \param [in] ghostSize - the ghost size of zone for all axis.
 904  */
 905 void MEDCouplingCartesianAMRMesh::ApplyGhostOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, int ghostSize)
 906 {
 907   if(ghostSize<0)
 908     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::ApplyGhostOnCompactFrmt : ghost size must be >= 0 !");
 909   std::size_t sz(partBeforeFact.size());
 910   for(std::size_t i=0;i<sz;i++)
 911     {
 912       partBeforeFact[i].first+=ghostSize;
 913       partBeforeFact[i].second+=ghostSize;
 914     }
 915 }
 916
 917 /*!
 918  * This method is different than ApplyGhostOnCompactFrmt
 919  *
 920  * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in ghost reference.
 921  * \param [in] ghostSize - the ghost size of zone for all axis.
 922  */
 923 void MEDCouplingCartesianAMRMesh::ApplyAllGhostOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, int ghostSize)
 924 {
 925   if(ghostSize<0)
 926     throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::ApplyAllGhostOnCompactFrmt : ghost size must be >= 0 !");
 927   std::size_t sz(partBeforeFact.size());
 928   for(std::size_t i=0;i<sz;i++)
 929     {
 930       partBeforeFact[i].first-=ghostSize;
 931       partBeforeFact[i].second+=ghostSize;
 932     }
 933 }
 934
 935 std::size_t MEDCouplingCartesianAMRMesh::getHeapMemorySizeWithoutChildren() const
 936 {
 937   return sizeof(MEDCouplingCartesianAMRMesh);
 938 }
 939
 940 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRMesh::getDirectChildren() const
 941 {
 942   std::vector<const BigMemoryObject *> ret;
 943   if((const MEDCouplingIMesh *)_mesh)
 944     ret.push_back((const MEDCouplingIMesh *)_mesh);
 945   for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
 946     {
 947       if((const MEDCouplingCartesianAMRPatch*)*it)
 948         ret.push_back((const MEDCouplingCartesianAMRPatch*)*it);
 949     }
 950   return ret;
 951 }
 952
 953 void MEDCouplingCartesianAMRMesh::updateTime() const
 954 {
 955   if((const MEDCouplingIMesh *)_mesh)
 956     updateTimeWith(*_mesh);
 957   for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
 958     {
 959       const MEDCouplingCartesianAMRPatch *elt(*it);
 960       if(!elt)
 961         continue;
 962       const MEDCouplingCartesianAMRMesh *mesh(elt->getMesh());
 963       if(mesh)
 964         updateTimeWith(*mesh);
 965     }
 966 }