1 // Copyright (C) 2007-2014 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
19 // Author : Anthony Geay
21 #include "MEDCouplingCartesianAMRMesh.hxx"
22 #include "MEDCouplingAMRAttribute.hxx"
23 #include "MEDCouplingFieldDouble.hxx"
24 #include "MEDCoupling1GTUMesh.hxx"
25 #include "MEDCouplingIMesh.hxx"
26 #include "MEDCouplingUMesh.hxx"
32 using namespace ParaMEDMEM;
36 int MEDCouplingCartesianAMRPatchGen::getNumberOfCellsRecursiveWithOverlap() const
38 return _mesh->getNumberOfCellsRecursiveWithOverlap();
41 int MEDCouplingCartesianAMRPatchGen::getNumberOfCellsRecursiveWithoutOverlap() const
43 return _mesh->getNumberOfCellsRecursiveWithoutOverlap();
46 int MEDCouplingCartesianAMRPatchGen::getMaxNumberOfLevelsRelativeToThis() const
48 return _mesh->getMaxNumberOfLevelsRelativeToThis();
51 MEDCouplingCartesianAMRPatchGen::MEDCouplingCartesianAMRPatchGen(MEDCouplingCartesianAMRMeshGen *mesh)
54 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatchGen constructor : input mesh is NULL !");
55 _mesh=mesh; _mesh->incrRef();
58 MEDCouplingCartesianAMRPatchGen::MEDCouplingCartesianAMRPatchGen(const MEDCouplingCartesianAMRPatchGen& other, MEDCouplingCartesianAMRMeshGen *father):RefCountObject(other),_mesh(other._mesh)
60 const MEDCouplingCartesianAMRMeshGen *mesh(other._mesh);
62 _mesh=mesh->deepCpy(father);
65 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatchGen::getMeshSafe() const
67 const MEDCouplingCartesianAMRMeshGen *mesh(_mesh);
69 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatchGen::getMeshSafe const : the mesh is NULL !");
73 MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatchGen::getMeshSafe()
75 MEDCouplingCartesianAMRMeshGen *mesh(_mesh);
77 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatchGen::getMeshSafe : the mesh is NULL !");
81 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRPatchGen::getDirectChildren() const
83 std::vector<const BigMemoryObject *> ret;
84 if((const MEDCouplingCartesianAMRMeshGen *)_mesh)
85 ret.push_back((const MEDCouplingCartesianAMRMeshGen *)_mesh);
90 * \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.
91 * \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,
92 * a the end cell (\b excluded) of the range for the second element of the pair.
94 MEDCouplingCartesianAMRPatch::MEDCouplingCartesianAMRPatch(MEDCouplingCartesianAMRMeshGen *mesh, const std::vector< std::pair<int,int> >& bottomLeftTopRight):MEDCouplingCartesianAMRPatchGen(mesh),_bl_tr(bottomLeftTopRight)
96 int dim((int)bottomLeftTopRight.size()),dimExp(_mesh->getSpaceDimension());
98 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch constructor : space dimension of father and input bottomLeft/topRight size mismatches !");
101 MEDCouplingCartesianAMRPatch *MEDCouplingCartesianAMRPatch::deepCpy(MEDCouplingCartesianAMRMeshGen *father) const
103 return new MEDCouplingCartesianAMRPatch(*this,father);
106 void MEDCouplingCartesianAMRPatch::addPatch(const std::vector< std::pair<int,int> >& bottomLeftTopRight, const std::vector<int>& factors)
108 return getMeshSafe()->addPatch(bottomLeftTopRight,factors);
111 int MEDCouplingCartesianAMRPatch::getNumberOfOverlapedCellsForFather() const
113 return MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(_bl_tr);
117 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
118 * the must be >= 0. \b WARNING this method only works if \a this and \a other share the same father (no check of this will be done !).
119 * Call isInMyNeighborhoodExt to deal with 2 patches not sharing directly the same father.
121 * \param [in] other - The other patch
122 * \param [in] ghostLev - The size of the neighborhood zone.
124 * \throw if \a this or \a other are invalid (end before start).
125 * \throw if \a ghostLev is \b not >= 0 .
126 * \throw if \a this and \a other have not the same space dimension.
128 * \sa isInMyNeighborhoodExt
130 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhood(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
133 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the size of the neighborhood must be >= 0 !");
135 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the input patch is NULL !");
136 const std::vector< std::pair<int,int> >& thisp(getBLTRRange());
137 const std::vector< std::pair<int,int> >& otherp(other->getBLTRRange());
138 return IsInMyNeighborhood(ghostLev==0?0:1,thisp,otherp);//make hypothesis that nb this->_mesh->getFather->getFactors() is >= ghostLev
142 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
143 * the must be >= 0. This method works even if \a this and \a other does not share the same father. But the level between their common
144 * ancestor must be the same. If they don't have the same ancestor an exception will be thrown.
146 * \param [in] other - The other patch
147 * \param [in] ghostLev - The size of the neighborhood zone.
149 * \throw if \a this or \a other are invalid (end before start).
150 * \throw if \a ghostLev is \b not >= 0 .
151 * \throw if \a this and \a other have not the same space dimension.
152 * \throw if there is not common ancestor of \a this and \a other.
154 * \sa isInMyNeighborhood, isInMyNeighborhoodDiffLev
156 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
159 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt : the size of the neighborhood must be >= 0 !");
161 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt : the input patch is NULL !");
163 const MEDCouplingCartesianAMRMeshGen *com(FindCommonAncestor(this,other,lev));//check that factors are OK
165 return isInMyNeighborhood(other,ghostLev);
166 std::vector<int> offset(ComputeOffsetFromTwoToOne(com,lev,this,other));
167 const std::vector< std::pair<int,int> >& thisp(getBLTRRange());
168 std::vector< std::pair<int,int> > otherp(other->getBLTRRange());
169 otherp=MEDCouplingStructuredMesh::TranslateCompactFrmt(otherp,offset);
170 return IsInMyNeighborhood(ghostLev,thisp,otherp);
174 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
175 * the must be >= 0. This method works even if \a this and \a other does not share the same father.
176 * \a this is expected to be more refined than \a other. That is to say lev of \a this is greater than level of \a other.
178 * \param [in] other - The other patch
179 * \param [in] ghostLev - The size of the neighborhood zone.
181 * \throw if \a this or \a other are invalid (end before start).
182 * \throw if \a ghostLev is \b not >= 0 .
183 * \throw if \a this and \a other have not the same space dimension.
184 * \throw if there is not common ancestor of \a this and \a other.
186 * \sa isInMyNeighborhoodExt
188 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhoodDiffLev(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
190 std::vector< std::pair<int,int> > thispp,otherpp;
191 std::vector<int> factors;
192 ComputeZonesOfTwoRelativeToOneDiffLev(ghostLev,this,other,thispp,otherpp,factors);
193 return IsInMyNeighborhood(ghostLev>0?1:0,thispp,otherpp);//1 not ghostLev ! It is not a bug ( I hope :) ) ! Because as \a this is a refinement of \a other ghostLev is supposed to be <= factors
196 std::vector< std::pair<int,int> > MEDCouplingCartesianAMRPatch::getBLTRRangeRelativeToGF() const
198 std::vector< std::pair<int,int> > ret(_bl_tr);
199 const MEDCouplingCartesianAMRMeshGen *mesh(getMesh());
201 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::getBLTRRangeRelativeToGF : not valid !");
202 const MEDCouplingCartesianAMRMeshGen *fath(mesh->getFather());
204 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::getBLTRRangeRelativeToGF : not valid 2 !");
205 std::vector<int> factors(fath->getFactors());
206 std::size_t sz(ret.size());
207 for(std::size_t ii=0;ii<sz;ii++)
209 ret[ii].first*=factors[ii];
210 ret[ii].second*=factors[ii];
212 const MEDCouplingCartesianAMRMeshGen *oldFather(fath);
213 fath=oldFather->getFather();
216 int pos(fath->getPatchIdFromChildMesh(oldFather));
217 const MEDCouplingCartesianAMRPatch *p(fath->getPatch(pos));
218 const std::vector< std::pair<int,int> >& tmp(p->getBLTRRange());
219 const std::vector<int>& factors2(fath->getFactors());
220 std::transform(factors.begin(),factors.end(),factors2.begin(),factors.begin(),std::multiplies<int>());
221 for(std::size_t ii=0;ii<sz;ii++)
223 int delta(ret[ii].second-ret[ii].first);
224 ret[ii].first+=tmp[ii].first*factors[ii];
225 ret[ii].second=ret[ii].first+delta;
228 fath=oldFather->getFather();
233 std::vector<int> MEDCouplingCartesianAMRPatch::computeCellGridSt() const
235 const MEDCouplingCartesianAMRMeshGen *m(getMesh());
237 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::computeCellGridSt : no mesh held by this !");
238 const MEDCouplingCartesianAMRMeshGen *father(m->getFather());
240 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::computeCellGridSt : no father help by underlying mesh !");
241 const std::vector< std::pair<int,int> >& bltr(getBLTRRange());
242 const std::vector<int>& factors(father->getFactors());
243 std::vector<int> ret(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(bltr));
244 std::transform(ret.begin(),ret.end(),factors.begin(),ret.begin(),std::multiplies<int>());
248 bool MEDCouplingCartesianAMRPatch::IsInMyNeighborhood(int ghostLev, const std::vector< std::pair<int,int> >& p1, const std::vector< std::pair<int,int> >& p2)
250 std::size_t thispsize(p1.size());
251 if(thispsize!=p2.size())
252 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the dimensions must be the same !");
253 for(std::size_t i=0;i<thispsize;i++)
255 const std::pair<int,int>& thispp(p1[i]);
256 const std::pair<int,int>& otherpp(p2[i]);
257 if(thispp.second<thispp.first)
258 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
259 if(otherpp.second<otherpp.first)
260 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
261 if(otherpp.first==thispp.second+ghostLev-1)
263 if(otherpp.second+ghostLev-1==thispp.first)
265 int start(std::max(thispp.first,otherpp.first)),end(std::min(thispp.second,otherpp.second));
273 * \sa FindNeighborsOfSubPatchesOf
275 std::vector< std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > > MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOfSameLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2)
278 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOfSameLev : the input pointers must be not NULL !");
279 std::vector< std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > > ret;
280 std::vector< const MEDCouplingCartesianAMRPatch *> p1Work(p1->getMesh()->getPatches()),p2Work(p2->getMesh()->getPatches());
281 while(!p1Work.empty())
283 std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > retTmp;
284 std::vector<const MEDCouplingCartesianAMRPatch *> p1Work2,p2Work2;
285 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it1=p1Work.begin();it1!=p1Work.end();it1++)
287 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it2=p2Work.begin();it2!=p2Work.end();it2++)
289 if((*it1)->isInMyNeighborhoodExt(*it2,ghostLev>0?1:0))//1 not ghostLev ! It is not a bug ( I hope :) ) ! Because as \a this is a refinement of \a other ghostLev is supposed to be <= factors
290 retTmp.push_back(std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *>(*it1,*it2));
292 std::vector<const MEDCouplingCartesianAMRPatch *> tmp1((*it1)->getMesh()->getPatches());
293 p1Work2.insert(p1Work2.end(),tmp1.begin(),tmp1.end());
295 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it2=p2Work.begin();it2!=p2Work.end();it2++)
297 std::vector<const MEDCouplingCartesianAMRPatch *> tmp2((*it2)->getMesh()->getPatches());
298 p2Work2.insert(p2Work2.end(),tmp2.begin(),tmp2.end());
300 ret.push_back(retTmp);
308 * This method returns all pair of patches (pa,pb) so that pb is in the neighborhood of pa (size of neighborhood is \a ghostLev).
309 * pa is a refinement (a child) of \b p1 and pb is equal to \a p2. So the returned pair do not have the same level as it is the case for
310 * FindNeighborsOfSubPatchesOfSameLev.
312 * \sa FindNeighborsOfSubPatchesOfSameLev
314 void MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOf(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, std::vector< std::pair<const MEDCouplingCartesianAMRPatch *, const MEDCouplingCartesianAMRPatch *> >& ret)
317 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOf : the input pointers must be not NULL !");
318 std::vector< const MEDCouplingCartesianAMRPatch *> p1Work(p1->getMesh()->getPatches());
319 while(!p1Work.empty())
321 std::vector<const MEDCouplingCartesianAMRPatch *> p1Work2;
322 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it0=p1Work.begin();it0!=p1Work.end();it0++)
324 if((*it0)->isInMyNeighborhoodDiffLev(p2,ghostLev))
325 ret.push_back(std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *>(*it0,p2));
326 std::vector<const MEDCouplingCartesianAMRPatch *> tmp2((*it0)->getMesh()->getPatches());
327 p1Work2.insert(p1Work2.end(),tmp2.begin(),tmp2.end());
334 * \a p1 and \a p2 are expected to be neighbors (inside the \a ghostLev zone). This method updates \a dataOnP1 only in the ghost part using a part of \a dataOnP2.
336 * \saUpdateNeighborsOfOneWithTwoExt
338 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(int ghostLev, const std::vector<int>& factors, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
340 const std::vector< std::pair<int,int> >& p1BLTR(p1->getBLTRRange());
341 const std::vector< std::pair<int,int> >& p2BLTR(p2->getBLTRRange());
342 UpdateNeighborsOfOneWithTwoInternal(ghostLev,factors,p1BLTR,p2BLTR,dataOnP1,dataOnP2);
346 * Idem than UpdateNeighborsOfOneWithTwo, except that here \a p1 and \a p2 are not sharing the same direct father.
348 * \sa UpdateNeighborsOfOneWithTwo
350 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoExt(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
352 const std::vector< std::pair<int,int> >& p1BLTR(p1->getBLTRRange());//p1BLTR=[(10,12),(5,8)]
353 std::vector< std::pair<int,int> > p2BLTR(p2->getBLTRRange());//p2BLTR=[(0,1),(0,5)]
355 const MEDCouplingCartesianAMRMeshGen *ca(FindCommonAncestor(p1,p2,lev));
356 std::vector<int> offset(ComputeOffsetFromTwoToOne(ca,lev,p1,p2));//[12,4]
357 p2BLTR=MEDCouplingStructuredMesh::TranslateCompactFrmt(p2BLTR,offset);//p2BLTR=[(12,13),(4,9)]
358 UpdateNeighborsOfOneWithTwoInternal(ghostLev,p1->getMesh()->getFather()->getFactors(),p1BLTR,p2BLTR,dataOnP1,dataOnP2);
362 * \a p1 is expected to be more refined than \a p2. \a p1 and \a p2 have to share a common ancestor. Compared to UpdateNeighborsOfOneWithTwoExt here \a p1 and \a p2 are \b not at the same level !
364 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoMixedLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2, bool isConservative)
366 std::vector< std::pair<int,int> > p1pp,p2pp;
367 std::vector<int> factors;
368 ComputeZonesOfTwoRelativeToOneDiffLev(ghostLev,p1,p2,p1pp,p2pp,factors);
370 std::vector<int> dimsP2NotRefined(p2->computeCellGridSt());
371 std::vector<int> dimsP2Refined(dimsP2NotRefined);
372 std::transform(dimsP2NotRefined.begin(),dimsP2NotRefined.end(),factors.begin(),dimsP2Refined.begin(),std::multiplies<int>());
373 std::vector< std::pair<int,int> > p2RefinedAbs(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(dimsP2NotRefined));
374 std::vector<int> dimsP2RefinedGhost(dimsP2Refined.size());
375 std::transform(dimsP2Refined.begin(),dimsP2Refined.end(),dimsP2RefinedGhost.begin(),std::bind2nd(std::plus<int>(),2*ghostLev));
376 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> fineP2(DataArrayDouble::New()); fineP2->alloc(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(dimsP2RefinedGhost),dataOnP2->getNumberOfComponents());
377 MEDCouplingIMesh::SpreadCoarseToFineGhost(dataOnP2,dimsP2NotRefined,fineP2,p2RefinedAbs,factors,ghostLev);
380 int fact(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(factors));
381 std::transform(fineP2->begin(),fineP2->end(),fineP2->getPointer(),std::bind2nd(std::multiplies<double>(),1./((double)fact)));
384 UpdateNeighborsOfOneWithTwoInternal(ghostLev,p1->getMesh()->getFather()->getFactors(),p1pp,p2pp,dataOnP1,fineP2);
388 * \a p1 is expected to be more refined than \a p2. \a p1 and \a p2 have to share a common ancestor. Compared to UpdateNeighborsOfOneWithTwoExt here \a p1 and \a p2 are \b not at the same level !
389 * This method has 3 outputs. 2 two first are the resp the position of \a p1 and \a p2 relative to \a p1. And \a factToApplyOn2 is the coeff of refinement to be applied on \a p2 to be virtualy
390 * on the same level as \a p1.
392 void MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, std::vector< std::pair<int,int> >& p1Zone, std::vector< std::pair<int,int> >& p2Zone, std::vector<int>& factToApplyOn2)
394 std::vector<const MEDCouplingCartesianAMRMeshGen *> ancestorsOfThis;
395 const MEDCouplingCartesianAMRMeshGen *work(p1->getMesh()),*work2(0);
396 ancestorsOfThis.push_back(work);
399 work=work->getFather();
401 ancestorsOfThis.push_back(work);
406 std::size_t levThis(0),levOther(0);
407 while(work && !found)
410 work=work->getFather();
414 std::vector<const MEDCouplingCartesianAMRMeshGen *>::iterator it(std::find(ancestorsOfThis.begin(),ancestorsOfThis.end(),work));
415 if(it!=ancestorsOfThis.end())
417 levThis=std::distance(ancestorsOfThis.begin(),it);
423 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev : no common ancestor found !");
424 if(levThis<=levOther)
425 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev : this method is not called correctly !");
427 const MEDCouplingCartesianAMRMeshGen *comAncestor(ancestorsOfThis[levThis]);
428 int idThis(comAncestor->getPatchIdFromChildMesh(ancestorsOfThis[levThis-1])),idOther(comAncestor->getPatchIdFromChildMesh(work2));
429 const MEDCouplingCartesianAMRPatch *thisp(comAncestor->getPatch(idThis)),*otherp(comAncestor->getPatch(idOther));
430 std::vector<int> offset(ComputeOffsetFromTwoToOne(comAncestor,levOther,thisp,otherp));
431 p1Zone=thisp->getBLTRRange(); p2Zone=MEDCouplingStructuredMesh::TranslateCompactFrmt(otherp->getBLTRRange(),offset);
432 factToApplyOn2.resize(p1Zone.size()); std::fill(factToApplyOn2.begin(),factToApplyOn2.end(),1);
434 std::size_t nbOfTurn(levThis-levOther);
435 for(std::size_t i=0;i<nbOfTurn;i++)
437 std::vector< std::pair<int,int> > tmp0;
438 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p1Zone,p2Zone,tmp0,false);
440 const MEDCouplingCartesianAMRMeshGen *curAncestor(ancestorsOfThis[levThis-i]);
441 ApplyFactorsOnCompactFrmt(p2Zone,curAncestor->getFactors());
442 curAncestor=ancestorsOfThis[levThis-1-i];
443 const std::vector<int>& factors(curAncestor->getFactors());
444 std::transform(factToApplyOn2.begin(),factToApplyOn2.end(),factors.begin(),factToApplyOn2.begin(),std::multiplies<int>());
445 int tmpId(curAncestor->getPatchIdFromChildMesh(ancestorsOfThis[levThis-2-i]));
446 p1Zone=curAncestor->getPatch(tmpId)->getBLTRRange();
450 std::size_t MEDCouplingCartesianAMRPatch::getHeapMemorySizeWithoutChildren() const
452 std::size_t ret(sizeof(MEDCouplingCartesianAMRPatch));
453 ret+=_bl_tr.capacity()*sizeof(std::pair<int,int>);
457 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatch::FindCommonAncestor(const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, int& lev)
459 const MEDCouplingCartesianAMRMeshGen *f1(p1->_mesh),*f2(p2->_mesh);
461 while(f1!=f2 || f1==0 || f2==0)
463 f1=f1->getFather(); f2=f2->getFather();
464 if(f1->getFactors()!=f2->getFactors())
465 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindCommonAncestor : factors differ !");
469 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindCommonAncestor : no common ancestor between p1 and p2 !");
473 std::vector<int> MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne(const MEDCouplingCartesianAMRMeshGen *comAncestor, int lev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2)
476 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne : this method is useful only for lev > 0 !");
478 int dim(p1->getMesh()->getSpaceDimension());
479 if(p2->getMesh()->getSpaceDimension()!=dim)
480 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne : dimension must be the same !");
481 std::vector< int > ret(dim,0);
482 for(int i=0;i<zeLev;i++)
484 const MEDCouplingCartesianAMRMeshGen *f1(p1->_mesh),*f2(p2->_mesh);
485 const MEDCouplingCartesianAMRPatch *p1h(0),*p2h(0);
486 for(int j=0;j<lev-i;j++)
488 const MEDCouplingCartesianAMRMeshGen *f1tmp(f1->getFather()),*f2tmp(f2->getFather());
489 int pid1(f1tmp->getPatchIdFromChildMesh(f1)),pid2(f2tmp->getPatchIdFromChildMesh(f2));
490 p1h=f1tmp->getPatch(pid1); p2h=f2tmp->getPatch(pid2);
493 std::vector< std::pair<int,int> > p2c(p2h->getBLTRRange());
494 for(int k=0;k<dim;k++)
496 p2c[k].first+=ret[k];
497 p2c[k].second+=ret[k];
499 for(int k=0;k<dim;k++)
501 ret[k]=p2c[k].first-p1h->getBLTRRange()[k].first;
502 ret[k]*=f1->getFactors()[k];
508 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoInternal(int ghostLev, const std::vector<int>& factors, const std::vector< std::pair<int,int> >&p1 ,const std::vector< std::pair<int,int> >&p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
509 {//p1=[(1,4),(2,4)] p2=[(4,5),(3,4)]
510 int dim((int)factors.size());
511 std::vector<int> dimsCoarse(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(p1));//[3,2]
512 std::transform(dimsCoarse.begin(),dimsCoarse.end(),factors.begin(),dimsCoarse.begin(),std::multiplies<int>());//[12,8]
513 std::transform(dimsCoarse.begin(),dimsCoarse.end(),dimsCoarse.begin(),std::bind2nd(std::plus<int>(),2*ghostLev));//[14,10]
514 std::vector< std::pair<int,int> > rangeCoarse(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(dimsCoarse));//[(0,14),(0,10)]
515 std::vector<int> fakeFactors(dim,1);
517 std::vector< std::pair<int,int> > tmp0,tmp1,tmp2;
518 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p1,p2,tmp0,false);//tmp0=[(3,4),(1,2)]
519 ApplyFactorsOnCompactFrmt(tmp0,factors);//tmp0=[(12,16),(4,8)]
520 MEDCouplingStructuredMesh::ApplyGhostOnCompactFrmt(tmp0,ghostLev);//tmp0=[(13,17),(5,9)]
521 std::vector< std::pair<int,int> > interstRange(MEDCouplingStructuredMesh::IntersectRanges(tmp0,rangeCoarse));//interstRange=[(13,14),(5,9)]
522 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p2,p1,tmp1,false);//tmp1=[(-3,0),(-1,1)]
523 ApplyFactorsOnCompactFrmt(tmp1,factors);//tmp1=[(-12,-4),(-4,0)]
524 MEDCouplingStructuredMesh::ChangeReferenceToGlobalOfCompactFrmt(tmp1,interstRange,tmp2,false);//tmp2=[(1,2),(1,5)]
526 std::vector< std::pair<int,int> > dimsFine(p2);
527 ApplyFactorsOnCompactFrmt(dimsFine,factors);
528 ApplyAllGhostOnCompactFrmt(dimsFine,ghostLev);
530 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ghostVals(MEDCouplingStructuredMesh::ExtractFieldOfDoubleFrom(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(dimsFine),dataOnP2,tmp2));
531 MEDCouplingIMesh::CondenseFineToCoarse(dimsCoarse,ghostVals,interstRange,fakeFactors,dataOnP1);
534 MEDCouplingCartesianAMRPatch::MEDCouplingCartesianAMRPatch(const MEDCouplingCartesianAMRPatch& other, MEDCouplingCartesianAMRMeshGen *father):MEDCouplingCartesianAMRPatchGen(other,father),_bl_tr(other._bl_tr)
539 * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in refined reference.
540 * \param [in] factors - the factors per axis.
542 void MEDCouplingCartesianAMRPatch::ApplyFactorsOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, const std::vector<int>& factors)
544 std::size_t sz(factors.size());
545 if(sz!=partBeforeFact.size())
546 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ApplyFactorsOnCompactFrmt : size of input vectors must be the same !");
547 for(std::size_t i=0;i<sz;i++)
549 partBeforeFact[i].first*=factors[i];
550 partBeforeFact[i].second*=factors[i];
555 * This method is different than ApplyGhostOnCompactFrmt. The \a partBeforeFact parameter is enlarger contrary to ApplyGhostOnCompactFrmt.
557 * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in ghost reference.
558 * \param [in] ghostSize - the ghost size of zone for all axis.
560 void MEDCouplingCartesianAMRPatch::ApplyAllGhostOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, int ghostSize)
563 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ApplyAllGhostOnCompactFrmt : ghost size must be >= 0 !");
564 std::size_t sz(partBeforeFact.size());
565 for(std::size_t i=0;i<sz;i++)
567 partBeforeFact[i].first-=ghostSize;
568 partBeforeFact[i].second+=ghostSize;
572 MEDCouplingCartesianAMRPatchGF::MEDCouplingCartesianAMRPatchGF(MEDCouplingCartesianAMRMesh *mesh):MEDCouplingCartesianAMRPatchGen(mesh)
576 MEDCouplingCartesianAMRPatchGF *MEDCouplingCartesianAMRPatchGF::deepCpy(MEDCouplingCartesianAMRMeshGen *father) const
578 return new MEDCouplingCartesianAMRPatchGF(*this,father);
581 std::size_t MEDCouplingCartesianAMRPatchGF::getHeapMemorySizeWithoutChildren() const
583 return sizeof(MEDCouplingCartesianAMRPatchGF);
586 MEDCouplingCartesianAMRPatchGF::MEDCouplingCartesianAMRPatchGF(const MEDCouplingCartesianAMRPatchGF& other, MEDCouplingCartesianAMRMeshGen *father):MEDCouplingCartesianAMRPatchGen(other,father)
592 int MEDCouplingCartesianAMRMeshGen::getSpaceDimension() const
594 return _mesh->getSpaceDimension();
597 void MEDCouplingCartesianAMRMeshGen::setFactors(const std::vector<int>& newFactors)
599 if(getSpaceDimension()!=(int)newFactors.size())
600 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::setFactors : size of input factors is not equal to the space dimension !");
606 if(_factors==newFactors)
608 if(!_patches.empty())
609 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::setFactors : modification of factors is not allowed when presence of patches !");
614 int MEDCouplingCartesianAMRMeshGen::getMaxNumberOfLevelsRelativeToThis() const
617 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
618 ret=std::max(ret,(*it)->getMaxNumberOfLevelsRelativeToThis()+1);
623 * This method returns the number of cells of \a this with the help of the MEDCouplingIMesh instance representing \a this.
624 * The patches in \a this are ignored here.
625 * \sa getNumberOfCellsAtCurrentLevelGhost, getNumberOfCellsRecursiveWithOverlap
627 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsAtCurrentLevel() const
629 return _mesh->getNumberOfCells();
633 * This method returns the number of cells of \a this with the help of the MEDCouplingIMesh instance representing \a this enlarged by \a ghostLev size
634 * to take into account of the ghost cells for future computation.
635 * The patches in \a this are ignored here.
637 * \sa getNumberOfCellsAtCurrentLevel
639 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsAtCurrentLevelGhost(int ghostLev) const
641 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> tmp(_mesh->buildWithGhost(ghostLev));
642 return tmp->getNumberOfCells();
646 * This method returns the number of cells including the current level but \b also \b including recursively all cells of other levels
647 * starting from this. The set of cells which size is returned here are generally overlapping each other.
649 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsRecursiveWithOverlap() const
651 int ret(_mesh->getNumberOfCells());
652 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
654 ret+=(*it)->getNumberOfCellsRecursiveWithOverlap();
660 * This method returns the max number of cells covering all the space without overlapping.
661 * It returns the number of cells of the mesh with the highest resolution.
662 * The returned value is equal to the number of cells of mesh returned by buildUnstructured.
664 * \sa buildUnstructured
666 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsRecursiveWithoutOverlap() const
668 int ret(_mesh->getNumberOfCells());
669 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
671 ret-=(*it)->getNumberOfOverlapedCellsForFather();
672 ret+=(*it)->getNumberOfCellsRecursiveWithoutOverlap();
678 * This method returns a vector of size equal to getAbsoluteLevelRelativeTo. It allows to find position an absolute position of \a this
679 * relative to \a ref (that is typically the god father).
681 * \sa getPatchAtPosition
683 std::vector<int> MEDCouplingCartesianAMRMeshGen::getPositionRelativeTo(const MEDCouplingCartesianAMRMeshGen *ref) const
686 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::getPositionRelativeTo : input pointer is NULL !");
687 std::vector<int> ret;
688 getPositionRelativeToInternal(ref,ret);
689 std::reverse(ret.begin(),ret.end());
694 * \sa getPositionRelativeTo, getMeshAtPosition
696 const MEDCouplingCartesianAMRPatch *MEDCouplingCartesianAMRMeshGen::getPatchAtPosition(const std::vector<int>& pos) const
698 std::size_t sz(pos.size());
700 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::getPatchAtPosition : empty input -> no patch by definition !");
702 const MEDCouplingCartesianAMRPatch *elt(getPatch(patchId));
705 if(!elt || !elt->getMesh())
706 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::getPatchAtPosition : NULL element found during walk !");
707 std::vector<int> pos2(pos.begin()+1,pos.end());
708 return elt->getMesh()->getPatchAtPosition(pos2);
711 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMeshGen::getMeshAtPosition(const std::vector<int>& pos) const
713 std::size_t sz(pos.size());
717 const MEDCouplingCartesianAMRPatch *elt(getPatch(patchId));
721 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::getMeshAtPosition : NULL patch !");
722 return elt->getMesh();
724 if(!elt || !elt->getMesh())
725 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::getPatchAtPosition : NULL element found during walk !");
726 std::vector<int> pos2(pos.begin()+1,pos.end());
727 return elt->getMesh()->getMeshAtPosition(pos2);
731 * This method returns grids relative to god father to specified level \a absoluteLev.
733 * \return std::vector<MEDCouplingCartesianAMRPatchGen *> - objects in vector are to be managed (decrRef) by the caller.
735 std::vector<MEDCouplingCartesianAMRPatchGen *> MEDCouplingCartesianAMRMeshGen::retrieveGridsAt(int absoluteLev) const
738 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::retrieveGridsAt : absolute level must be >=0 !");
739 return getGodFather()->retrieveGridsAt(absoluteLev);
743 * \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,
744 * a the end cell (\b excluded) of the range for the second element of the pair.
745 * \param [in] factors The factor of refinement per axis (different from 0).
747 void MEDCouplingCartesianAMRMeshGen::addPatch(const std::vector< std::pair<int,int> >& bottomLeftTopRight, const std::vector<int>& factors)
749 checkFactorsAndIfNotSetAssign(factors);
750 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> mesh(static_cast<MEDCouplingIMesh *>(_mesh->buildStructuredSubPart(bottomLeftTopRight)));
751 mesh->refineWithFactor(factors);
752 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRMeshSub> zeMesh(new MEDCouplingCartesianAMRMeshSub(this,mesh));
753 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> elt(new MEDCouplingCartesianAMRPatch(zeMesh,bottomLeftTopRight));
754 _patches.push_back(elt);
760 class InternalPatch : public RefCountObjectOnly
763 InternalPatch():_nb_of_true(0) { }
764 int getDimension() const { return (int)_part.size(); }
765 double getEfficiency() const { return (double)_nb_of_true/(double)_crit.size(); }
766 int getNumberOfCells() const { return (int)_crit.size(); }
767 void setNumberOfTrue(int nboft) { _nb_of_true=nboft; }
768 std::vector<bool>& getCriterion() { return _crit; }
769 const std::vector<bool>& getConstCriterion() const { return _crit; }
770 void setPart(const std::vector< std::pair<int,int> >& part) { _part=part; }
771 std::vector< std::pair<int,int> >& getPart() { return _part; }
772 const std::vector< std::pair<int,int> >& getConstPart() const { return _part; }
773 bool presenceOfTrue() const { return _nb_of_true>0; }
774 std::vector<int> computeCGS() const { return MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(_part); }
775 std::vector< std::vector<int> > computeSignature() const { return MEDCouplingStructuredMesh::ComputeSignaturePerAxisOf(computeCGS(),getConstCriterion()); }
776 double getEfficiencyPerAxis(int axisId) const { return (double)_nb_of_true/((double)(_part[axisId].second-_part[axisId].first)); }
777 void zipToFitOnCriterion(int minPatchLgth);
778 void updateNumberOfTrue() const;
779 MEDCouplingAutoRefCountObjectPtr<InternalPatch> extractPart(const std::vector< std::pair<int,int> >&partInGlobal) const;
780 MEDCouplingAutoRefCountObjectPtr<InternalPatch> deepCpy() const;
784 mutable int _nb_of_true;
785 std::vector<bool> _crit;
787 std::vector< std::pair<int,int> > _part;
790 void InternalPatch::zipToFitOnCriterion(int minPatchLgth)
792 std::vector<int> cgs(computeCGS());
793 std::vector<bool> newCrit;
794 std::vector< std::pair<int,int> > newPart,newPart2;
795 int newNbOfTrue(MEDCouplingStructuredMesh::FindMinimalPartOf(minPatchLgth,cgs,_crit,newCrit,newPart));
796 MEDCouplingStructuredMesh::ChangeReferenceToGlobalOfCompactFrmt(_part,newPart,newPart2);
797 if(newNbOfTrue!=_nb_of_true)
798 throw INTERP_KERNEL::Exception("InternalPatch::zipToFitOnCrit : internal error !");
799 _crit=newCrit; _part=newPart2;
802 void InternalPatch::updateNumberOfTrue() const
804 _nb_of_true=(int)std::count(_crit.begin(),_crit.end(),true);
807 MEDCouplingAutoRefCountObjectPtr<InternalPatch> InternalPatch::extractPart(const std::vector< std::pair<int,int> >&partInGlobal) const
809 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(new InternalPatch);
810 std::vector<int> cgs(computeCGS());
811 std::vector< std::pair<int,int> > newPart;
812 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(_part,partInGlobal,newPart);
813 MEDCouplingStructuredMesh::ExtractFieldOfBoolFrom(cgs,_crit,newPart,ret->getCriterion());
814 ret->setPart(partInGlobal);
815 ret->updateNumberOfTrue();
819 MEDCouplingAutoRefCountObjectPtr<InternalPatch> InternalPatch::deepCpy() const
821 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(new InternalPatch);
826 void DissectBigPatch(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int largestLength, int& cutPlace)
828 int minimumPatchLength(bso.getMinimumPatchLength());
829 std::vector<double> ratio(largestLength-minimumPatchLength,std::numeric_limits<double>::max());
831 double minSemiEfficiencyRatio(std::numeric_limits<double>::max());
832 double efficiencyPerAxis[2];
834 for(int i=minimumPatchLength-1;i<largestLength-minimumPatchLength;i++)
838 std::vector< std::pair<int,int> > rectH(patchToBeSplit->getConstPart());
840 rectH[axisId].second=patchToBeSplit->getConstPart()[axisId].first+i;
842 rectH[axisId].first=patchToBeSplit->getConstPart()[axisId].first+i;
843 MEDCouplingAutoRefCountObjectPtr<InternalPatch> p(patchToBeSplit->deepCpy());
844 p->zipToFitOnCriterion(bso.getMinimumPatchLength());
845 efficiencyPerAxis[h]=p->getEfficiencyPerAxis(axisId);
847 ratio[i]=std::max(efficiencyPerAxis[0], efficiencyPerAxis[1]) / std::min(efficiencyPerAxis[0], efficiencyPerAxis[1]);
848 if(ratio[i]<minSemiEfficiencyRatio)
850 minSemiEfficiencyRatio = ratio[i];
856 throw INTERP_KERNEL::Exception("DissectBigPatch : just call to Arthur !");
858 cutPlace=index_min+patchToBeSplit->getConstPart()[axisId].first;
861 bool FindHole(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int& cutPlace)
864 int minimumPatchLength(bso.getMinimumPatchLength());
865 const int dim(patchToBeSplit->getDimension());
866 std::vector< std::vector<int> > signatures(patchToBeSplit->computeSignature());
867 for(int id=0;id<dim;id++)
869 const std::vector<int>& signature(signatures[id]);
870 std::vector<int> hole;
871 std::vector<double> distance;
872 int len((int)signature.size());
873 for(int i=minimumPatchLength-1;i<len-minimumPatchLength;i++)
878 int closestHoleToMiddle(hole[0]);
879 int oldDistanceToMiddle(std::abs(hole[0]-len/2));
880 int newDistanceToMiddle(oldDistanceToMiddle);
881 for(std::size_t i=0;i<hole.size();i++)
883 newDistanceToMiddle=std::abs(hole[i]-len/2);
884 if(newDistanceToMiddle < oldDistanceToMiddle)
886 oldDistanceToMiddle = newDistanceToMiddle;
887 closestHoleToMiddle = hole[i];
890 cutPlace=closestHoleToMiddle+patchToBeSplit->getConstPart()[axisId].first;
897 bool FindInflection(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int& cutPlace, int& axisId)
899 bool cutFound(false); cutPlace=-1;// do not set axisId before to be sure that cutFound was set to true
900 const std::vector< std::pair<int,int> >& part(patchToBeSplit->getConstPart());
901 int sign,minimumPatchLength(bso.getMinimumPatchLength());
902 const int dim(patchToBeSplit->getDimension());
904 std::vector<int> zeroCrossDims(dim,-1);
905 std::vector<int> zeroCrossVals(dim,-1);
906 std::vector< std::vector<int> > signatures(patchToBeSplit->computeSignature());
907 for (int id=0;id<dim;id++)
909 const std::vector<int>& signature(signatures[id]);
911 std::vector<int> derivate_second_order,gradient_absolute,zero_cross,edge,max_cross_list ;
912 std::vector<double> distance ;
914 for(std::size_t i=1;i<signature.size()-1;i++)
915 derivate_second_order.push_back(signature[i-1]-2*signature[i]+signature[i+1]) ;
917 // Gradient absolute value
918 for(std::size_t i=1;i<derivate_second_order.size();i++)
919 gradient_absolute.push_back(abs(derivate_second_order[i]-derivate_second_order[i-1])) ;
920 if(derivate_second_order.empty())
922 for(std::size_t i=1;i<derivate_second_order.size()-1;i++)
924 if (derivate_second_order[i]*derivate_second_order[i+1] < 0 )
926 if (derivate_second_order[i]*derivate_second_order[i+1] > 0 )
928 if (derivate_second_order[i]*derivate_second_order[i+1] == 0 )
930 if ( sign==0 || sign==-1 )
931 if ( i >= (std::size_t)minimumPatchLength-2 && i <= signature.size()-minimumPatchLength-2 )
933 zero_cross.push_back(i) ;
934 edge.push_back(gradient_absolute[i]) ;
937 if ( zero_cross.size() > 0 )
939 int max_cross=*max_element(edge.begin(),edge.end()) ;
940 for (unsigned int i=0;i<edge.size();i++)
941 if (edge[i]==max_cross)
942 max_cross_list.push_back(zero_cross[i]+1) ;
944 double center((signature.size()/2.0));
945 for (unsigned int i=0;i<max_cross_list.size();i++)
946 distance.push_back(fabs(max_cross_list[i]+1-center));
948 double distance_min=*min_element(distance.begin(),distance.end()) ;
949 int pos_distance_min=find(distance.begin(),distance.end(),distance_min)-distance.begin();
950 int best_place = max_cross_list[pos_distance_min] + part[id].first ;
953 zeroCrossDims[id] = best_place ;
954 zeroCrossVals[id] = max_cross ;
957 derivate_second_order.clear() ;
958 gradient_absolute.clear() ;
961 max_cross_list.clear() ;
965 if ( zeroCrossDims[0]!=-1 || zeroCrossDims[1]!=-1 )
967 int max_cross_dims = *max_element(zeroCrossVals.begin(),zeroCrossVals.end()) ;
969 if (zeroCrossVals[0]==max_cross_dims && zeroCrossVals[1]==max_cross_dims )
971 int nl_left(part[0].second-part[0].first);
972 int nc_left(part[1].second-part[1].first);
973 if ( nl_left >= nc_left )
979 max_cross_dims=std::find(zeroCrossVals.begin(),zeroCrossVals.end(),max_cross_dims)-zeroCrossVals.begin();
981 cutPlace=zeroCrossDims[max_cross_dims];
982 axisId=max_cross_dims ;
987 bool TryAction4(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int rangeOfAxisId, int& cutPlace)
989 if(patchToBeSplit->getEfficiency()<=bso.getEfficiencyGoal())
991 if(rangeOfAxisId>=2*bso.getMinimumPatchLength())
993 cutPlace=rangeOfAxisId/2+patchToBeSplit->getConstPart()[axisId].first-1;
1000 if(patchToBeSplit->getNumberOfCells()>bso.getMaximumNbOfCellsInPatch() || rangeOfAxisId>bso.getMaximumPatchLength())
1002 DissectBigPatch(bso,patchToBeSplit,axisId,rangeOfAxisId,cutPlace);
1010 MEDCouplingAutoRefCountObjectPtr<InternalPatch> DealWithNoCut(const InternalPatch *patch)
1012 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(const_cast<InternalPatch *>(patch));
1017 void DealWithCut(double minPatchLgth, const InternalPatch *patchToBeSplit, int axisId, int cutPlace, std::vector<MEDCouplingAutoRefCountObjectPtr<InternalPatch> >& listOfPatches)
1019 MEDCouplingAutoRefCountObjectPtr<InternalPatch> leftPart,rightPart;
1020 std::vector< std::pair<int,int> > rect(patchToBeSplit->getConstPart());
1021 std::vector< std::pair<int,int> > leftRect(rect),rightRect(rect);
1022 leftRect[axisId].second=cutPlace+1;
1023 rightRect[axisId].first=cutPlace+1;
1024 leftPart=patchToBeSplit->extractPart(leftRect);
1025 rightPart=patchToBeSplit->extractPart(rightRect);
1026 leftPart->zipToFitOnCriterion(minPatchLgth); rightPart->zipToFitOnCriterion(minPatchLgth);
1027 listOfPatches.push_back(leftPart);
1028 listOfPatches.push_back(rightPart);
1033 void MEDCouplingCartesianAMRMeshGen::removeAllPatches()
1039 void MEDCouplingCartesianAMRMeshGen::removePatch(int patchId)
1041 checkPatchId(patchId);
1042 int sz((int)_patches.size()),j(0);
1043 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> > patches(sz-1);
1044 for(int i=0;i<sz;i++)
1046 patches[j++]=_patches[i];
1047 (const_cast<MEDCouplingCartesianAMRMeshGen *>(_patches[patchId]->getMesh()))->detachFromFather();
1052 int MEDCouplingCartesianAMRMeshGen::getNumberOfPatches() const
1054 return (int)_patches.size();
1058 * This method is a generic algorithm to create patches in \a this (by destroying the patches if any).
1059 * This method uses \a criterion array as a field on cells on this level.
1060 * This method only create patches at level 0 relative to \a this.
1062 * This generic algorithm can be degenerated into three child ones, depending on the arguments given; in particular depending
1063 * on whether they are equal to 0 or not.
1064 * 1/ If minimumPatchLength = maximumPatchLength = maximumPatchVolume = 0, then we have the Berger-Rigoutsos algorithm.
1065 * This algorithm was developed in 1991 and seems appropriate for sequential programming.
1066 * 2/ If maximumPatchLength = 0, then we have the Livne algorithm.
1067 * This algorithm was developed in 2004 and is an improvement of the Berger-Rigoutsos algorithm.
1068 * 3/ If maximumPatchVolume = 0, the we have the lmin-lmax algorithm.
1069 * This algorithm was developed by Arthur TALPAERT in 2014 and is an improvement of the Livne algorithm. It is especially
1070 * appropriate for parallel computing, where one patch would be given to one CPU. See Arthur TALPAERT's 2014 CANUM poster
1071 * for more information.
1074 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const std::vector<bool>& criterion, const std::vector<int>& factors)
1076 int nbCells(getNumberOfCellsAtCurrentLevel());
1077 if(nbCells!=(int)criterion.size())
1078 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion : the number of tuples of criterion array must be equal to the number of cells at the current level !");
1080 std::vector<int> cgs(_mesh->getCellGridStructure());
1081 std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> > listOfPatches,listOfPatchesOK;
1083 MEDCouplingAutoRefCountObjectPtr<InternalPatch> p(new InternalPatch);
1084 p->setNumberOfTrue(MEDCouplingStructuredMesh::FindMinimalPartOf(bso.getMinimumPatchLength(),cgs,criterion,p->getCriterion(),p->getPart()));
1085 if(p->presenceOfTrue())
1086 listOfPatches.push_back(p);
1087 while(!listOfPatches.empty())
1089 std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> > listOfPatchesTmp;
1090 for(std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> >::iterator it=listOfPatches.begin();it!=listOfPatches.end();it++)
1093 int axisId,largestLength,cutPlace;
1094 MEDCouplingStructuredMesh::FindTheWidestAxisOfGivenRangeInCompactFrmt((*it)->getConstPart(),axisId,largestLength);
1095 if((*it)->getEfficiency()>=bso.getEfficiencyThreshold() && ((*it)->getNumberOfCells()>bso.getMaximumNbOfCellsInPatch() || largestLength>bso.getMaximumPatchLength()))
1097 DissectBigPatch(bso,*it,axisId,largestLength,cutPlace);
1098 DealWithCut(bso.getMinimumPatchLength(),*it,axisId,cutPlace,listOfPatchesTmp);
1101 if(FindHole(bso,*it,axisId,cutPlace))//axisId overwritten here if FindHole equal to true !
1102 { DealWithCut(bso.getMinimumPatchLength(),*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 2
1103 if(FindInflection(bso,*it,cutPlace,axisId))//axisId overwritten here if cutFound equal to true !
1104 { DealWithCut(bso.getMinimumPatchLength(),*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 3
1105 if(TryAction4(bso,*it,axisId,largestLength,cutPlace))
1106 { DealWithCut(bso.getMinimumPatchLength(),*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 4
1108 listOfPatchesOK.push_back(DealWithNoCut(*it));
1110 listOfPatches=listOfPatchesTmp;
1112 for(std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> >::const_iterator it=listOfPatchesOK.begin();it!=listOfPatchesOK.end();it++)
1113 addPatch((*it)->getConstPart(),factors);
1118 * 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.
1119 * This method only create patches at level 0 relative to \a this.
1121 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const DataArrayByte *criterion, const std::vector<int>& factors)
1123 if(!criterion || !criterion->isAllocated())
1124 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion : the criterion DataArrayByte instance must be allocated and not NULL !");
1125 std::vector<bool> crit(criterion->toVectorOfBool());//check that criterion has one component.
1126 createPatchesFromCriterion(bso,crit,factors);
1130 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const DataArrayDouble *criterion, const std::vector<int>& factors, double eps)
1133 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion : null criterion pointer !");
1134 std::vector<bool> inp(criterion->toVectorOfBool(eps));
1135 createPatchesFromCriterion(bso,inp,factors);
1138 int MEDCouplingCartesianAMRMeshGen::getPatchIdFromChildMesh(const MEDCouplingCartesianAMRMeshGen *mesh) const
1141 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ret++)
1143 if((*it)->getMesh()==mesh)
1146 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::getPatchIdFromChildMesh : no such a mesh in my direct progeny !");
1149 std::vector< const MEDCouplingCartesianAMRPatch *> MEDCouplingCartesianAMRMeshGen::getPatches() const
1151 std::size_t sz(_patches.size());
1152 std::vector< const MEDCouplingCartesianAMRPatch *> ret(sz);
1153 for(std::size_t i=0;i<sz;i++)
1158 const MEDCouplingCartesianAMRPatch *MEDCouplingCartesianAMRMeshGen::getPatch(int patchId) const
1160 checkPatchId(patchId);
1161 return _patches[patchId];
1165 * This method states if patch2 (with id \a patchId2) is in the neighborhood of patch1 (with id \a patchId1).
1166 * The neighborhood size is defined by \a ghostLev in the reference of \a this ( \b not in the reference of patches !).
1168 bool MEDCouplingCartesianAMRMeshGen::isPatchInNeighborhoodOf(int patchId1, int patchId2, int ghostLev) const
1170 const MEDCouplingCartesianAMRPatch *p1(getPatch(patchId1)),*p2(getPatch(patchId2));
1171 return p1->isInMyNeighborhood(p2,ghostLev);
1175 * 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.
1176 * 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
1177 * defined by the patch with id \a patchId.
1179 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1180 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1181 * \return DataArrayDouble * - The array of the cell field on the requested patch
1183 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1184 * \throw if \a cellFieldOnThis is NULL or not allocated
1185 * \sa fillCellFieldOnPatch, MEDCouplingIMesh::SpreadCoarseToFine
1187 DataArrayDouble *MEDCouplingCartesianAMRMeshGen::createCellFieldOnPatch(int patchId, const DataArrayDouble *cellFieldOnThis) const
1189 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1190 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatch : the input cell field array is NULL or not allocated !");
1191 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1192 const MEDCouplingIMesh *fine(patch->getMesh()->getImageMesh());
1193 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(fine->getNumberOfCells(),cellFieldOnThis->getNumberOfComponents());
1194 ret->copyStringInfoFrom(*cellFieldOnThis);
1195 MEDCouplingIMesh::SpreadCoarseToFine(cellFieldOnThis,_mesh->getCellGridStructure(),ret,patch->getBLTRRange(),getFactors());
1200 * This method is equivalent to MEDCouplingCartesianAMRMesh::createCellFieldOnPatch except that here instead of creating a new instance
1201 * it fills the value into the \a cellFieldOnPatch data.
1203 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1204 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1205 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1207 * \sa createCellFieldOnPatch, fillCellFieldComingFromPatch
1209 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatch(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch, bool isConservative) const
1211 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1212 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatch : the input cell field array is NULL or not allocated !");
1213 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1214 MEDCouplingIMesh::SpreadCoarseToFine(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors());
1217 int fact(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(getFactors()));
1218 std::transform(cellFieldOnPatch->begin(),cellFieldOnPatch->end(),cellFieldOnPatch->getPointer(),std::bind2nd(std::multiplies<double>(),1./((double)fact)));
1223 * This method is the generalization of fillCellFieldOnPatch method. This method only projects coarse to fine without considering the
1224 * potential neighbor patches covered by the ghost cells of patch with id \a patchId.
1226 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1227 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1228 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1229 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1231 * \sa fillCellFieldOnPatch, fillCellFieldOnPatchGhostAdv
1233 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchGhost(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch, int ghostLev, bool isConservative) const
1235 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1236 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatchGhost : the input cell field array is NULL or not allocated !");
1237 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1238 MEDCouplingIMesh::SpreadCoarseToFineGhost(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),ghostLev);
1241 int fact(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(getFactors()));
1242 std::transform(cellFieldOnPatch->begin(),cellFieldOnPatch->end(),cellFieldOnPatch->getPointer(),std::bind2nd(std::multiplies<double>(),1./((double)fact)));
1247 * This method is equivalent to fillCellFieldOnPatchGhost except that here \b ONLY \b the \b ghost \b zone will be updated
1248 * in \a cellFieldOnPatch.
1250 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1251 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1252 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled \b only \b in \b the \b ghost \b zone.
1253 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1255 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyOnGhostZone(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch, int ghostLev) const
1257 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1258 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchOnlyOnGhostZone : the input cell field array is NULL or not allocated !");
1259 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1260 MEDCouplingIMesh::SpreadCoarseToFineGhostZone(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),ghostLev);
1264 * This method is a refinement of fillCellFieldOnPatchGhost. fillCellFieldOnPatchGhost is first called.
1265 * Then for all other patches than those pointed by \a patchId that overlap the ghost zone of the patch impact the ghost zone adequately.
1267 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1268 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1269 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1270 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1271 * \param [in] arrsOnPatches - \b WARNING arrsOnPatches[patchId] is \b NOT \b const. All others are const.
1273 * \sa fillCellFieldOnPatchOnlyGhostAdv
1275 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchGhostAdv(int patchId, const DataArrayDouble *cellFieldOnThis, int ghostLev, const std::vector<const DataArrayDouble *>& arrsOnPatches, bool isConservative) const
1277 int nbp(getNumberOfPatches());
1278 if(nbp!=(int)arrsOnPatches.size())
1280 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchGhostAdv : there are " << nbp << " patches in this and " << arrsOnPatches.size() << " arrays in the last parameter !";
1281 throw INTERP_KERNEL::Exception(oss.str().c_str());
1283 DataArrayDouble *theFieldToFill(const_cast<DataArrayDouble *>(arrsOnPatches[patchId]));
1284 // first, do as usual
1285 fillCellFieldOnPatchGhost(patchId,cellFieldOnThis,theFieldToFill,ghostLev,isConservative);
1286 fillCellFieldOnPatchOnlyGhostAdv(patchId,ghostLev,arrsOnPatches);
1290 * This method updates the patch with id \a patchId considering the only the all the patches in \a this to fill ghost zone.
1291 * So \b warning, the DataArrayDouble instance \a arrsOnPatches[patchId] is non const.
1293 * \sa getPatchIdsInTheNeighborhoodOf
1295 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyGhostAdv(int patchId, int ghostLev, const std::vector<const DataArrayDouble *>& arrsOnPatches) const
1297 int nbp(getNumberOfPatches());
1298 if(nbp!=(int)arrsOnPatches.size())
1300 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchOnlyGhostAdv : there are " << nbp << " patches in this and " << arrsOnPatches.size() << " arrays in the last parameter !";
1301 throw INTERP_KERNEL::Exception(oss.str().c_str());
1303 const MEDCouplingCartesianAMRPatch *refP(getPatch(patchId));
1304 DataArrayDouble *theFieldToFill(const_cast<DataArrayDouble *>(arrsOnPatches[patchId]));
1305 std::vector<int> ids(getPatchIdsInTheNeighborhoodOf(patchId,ghostLev));
1306 for(std::vector<int>::const_iterator it=ids.begin();it!=ids.end();it++)
1308 const MEDCouplingCartesianAMRPatch *otherP(getPatch(*it));
1309 MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(ghostLev,_factors,refP,otherP,theFieldToFill,arrsOnPatches[*it]);
1313 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyOnGhostZoneWith(int ghostLev, const MEDCouplingCartesianAMRPatch *patchToBeModified, const MEDCouplingCartesianAMRPatch *neighborPatch, DataArrayDouble *cellFieldOnPatch, const DataArrayDouble *cellFieldNeighbor) const
1315 MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(ghostLev,_factors,patchToBeModified,neighborPatch,cellFieldOnPatch,cellFieldNeighbor);
1319 * This method updates \a cellFieldOnThis part of values coming from the cell field \a cellFieldOnPatch lying on patch having id \a patchId.
1321 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1322 * \param [in] cellFieldOnPatch - The array of the cell field on patch with id \a patchId.
1323 * \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.
1324 * \param [in] isConservative - true if the field needs to be conserved. false if maximum principle has to be applied.
1326 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1327 * \throw if \a cellFieldOnPatch is NULL or not allocated
1328 * \sa createCellFieldOnPatch, MEDCouplingIMesh::CondenseFineToCoarse,fillCellFieldComingFromPatchGhost
1330 void MEDCouplingCartesianAMRMeshGen::fillCellFieldComingFromPatch(int patchId, const DataArrayDouble *cellFieldOnPatch, DataArrayDouble *cellFieldOnThis, bool isConservative) const
1332 if(!cellFieldOnPatch || !cellFieldOnPatch->isAllocated())
1333 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch : the input cell field array is NULL or not allocated !");
1334 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1335 MEDCouplingIMesh::CondenseFineToCoarse(_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),cellFieldOnThis);
1338 int fact(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(getFactors()));
1339 MEDCouplingStructuredMesh::MultiplyPartOf(_mesh->getCellGridStructure(),patch->getBLTRRange(),1./((double)fact),cellFieldOnThis);
1344 * This method is the extension of MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch managing the ghost cells. If this
1345 * method is called with \a ghostLev equal to 0 it behaves exactly as MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch.
1347 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1348 * \param [in] cellFieldOnPatch - The array of the cell field on patch with id \a patchId.
1349 * \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.
1350 * \param [in] ghostLev The size of ghost zone (must be >= 0 !)
1351 * \param [in] isConservative - true if the field needs to be conserved. false if maximum principle has to be applied.
1353 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1354 * \throw if \a cellFieldOnPatch is NULL or not allocated
1355 * \sa fillCellFieldComingFromPatch
1357 void MEDCouplingCartesianAMRMeshGen::fillCellFieldComingFromPatchGhost(int patchId, const DataArrayDouble *cellFieldOnPatch, DataArrayDouble *cellFieldOnThis, int ghostLev, bool isConservative) const
1359 if(!cellFieldOnPatch || !cellFieldOnPatch->isAllocated())
1360 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatchGhost : the input cell field array is NULL or not allocated !");
1361 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1362 MEDCouplingIMesh::CondenseFineToCoarseGhost(_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),cellFieldOnThis,ghostLev);
1365 int fact(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(getFactors()));
1366 MEDCouplingStructuredMesh::MultiplyPartOfByGhost(_mesh->getCellGridStructure(),patch->getBLTRRange(),ghostLev,1./((double)fact),cellFieldOnThis);
1371 * This method finds all patches (located by their ids) that are in the neighborhood of patch with id \a patchId. The neighborhood size is
1372 * defined by ghostLev.
1374 * \param [in] patchId - the id of the considered patch.
1375 * \param [in] ghostLev - the size of the neighborhood.
1376 * \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.
1378 DataArrayInt *MEDCouplingCartesianAMRMeshGen::findPatchesInTheNeighborhoodOf(int patchId, int ghostLev) const
1380 int nbp(getNumberOfPatches());
1381 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(0,1);
1382 for(int i=0;i<nbp;i++)
1385 if(isPatchInNeighborhoodOf(i,patchId,ghostLev))
1386 ret->pushBackSilent(i);
1391 MEDCouplingUMesh *MEDCouplingCartesianAMRMeshGen::buildUnstructured() const
1393 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> part(_mesh->buildUnstructured());
1394 std::vector<bool> bs(_mesh->getNumberOfCells(),false);
1395 std::vector<int> cgs(_mesh->getCellGridStructure());
1396 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> > msSafe(_patches.size()+1);
1398 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ii++)
1400 MEDCouplingStructuredMesh::SwitchOnIdsFrom(cgs,(*it)->getBLTRRange(),bs);
1401 msSafe[ii+1]=(*it)->getMesh()->buildUnstructured();
1403 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsOff(DataArrayInt::BuildListOfSwitchedOff(bs));
1404 msSafe[0]=static_cast<MEDCouplingUMesh *>(part->buildPartOfMySelf(eltsOff->begin(),eltsOff->end(),false));
1405 std::vector< const MEDCouplingUMesh * > ms(msSafe.size());
1406 for(std::size_t i=0;i<msSafe.size();i++)
1408 return MEDCouplingUMesh::MergeUMeshes(ms);
1412 * This method returns a mesh containing as cells that there is patches at the current level.
1413 * The patches are seen like 'boxes' that is too say the refinement will not appear here.
1415 * \return MEDCoupling1SGTUMesh * - A new object to be managed by the caller containing as cells as there are patches in \a this.
1417 MEDCoupling1SGTUMesh *MEDCouplingCartesianAMRMeshGen::buildMeshFromPatchEnvelop() const
1419 std::vector<const MEDCoupling1SGTUMesh *> cells;
1420 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > cellsSafe;
1421 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1423 const MEDCouplingCartesianAMRPatch *patch(*it);
1426 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> cell(patch->getMesh()->getImageMesh()->asSingleCell());
1427 MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> cell1SGT(cell->build1SGTUnstructured());
1428 cellsSafe.push_back(cell1SGT); cells.push_back(cell1SGT);
1431 return MEDCoupling1SGTUMesh::Merge1SGTUMeshes(cells);
1434 MEDCoupling1SGTUMesh *MEDCouplingCartesianAMRMeshGen::buildMeshOfDirectChildrenOnly() const
1436 std::vector<const MEDCoupling1SGTUMesh *> patches;
1437 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > patchesSafe;
1438 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1440 const MEDCouplingCartesianAMRPatch *patch(*it);
1443 MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> patchMesh(patch->getMesh()->getImageMesh()->build1SGTUnstructured());
1444 patchesSafe.push_back(patchMesh); patches.push_back(patchMesh);
1447 return MEDCoupling1SGTUMesh::Merge1SGTUMeshes(patches);
1451 * This method works same as buildUnstructured except that arrays are given in input to build a field on cell in output.
1452 * \return MEDCouplingFieldDouble * - a newly created instance the caller has reponsability to deal with.
1453 * \sa buildUnstructured
1455 MEDCouplingFieldDouble *MEDCouplingCartesianAMRMeshGen::buildCellFieldOnRecurseWithoutOverlapWithoutGhost(int ghostSz, const std::vector<const DataArrayDouble *>& recurseArrs) const
1457 if(recurseArrs.empty())
1458 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::buildCellFieldOnRecurseWithoutOverlapWithoutGhost : array is empty ! Should never happen !");
1460 std::vector<bool> bs(_mesh->getNumberOfCells(),false);
1461 std::vector<int> cgs(_mesh->getCellGridStructure());
1462 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> > msSafe(_patches.size()+1);
1464 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ii++)
1466 MEDCouplingStructuredMesh::SwitchOnIdsFrom(cgs,(*it)->getBLTRRange(),bs);
1467 std::vector<const DataArrayDouble *> tmpArrs(extractSubTreeFromGlobalFlatten((*it)->getMesh(),recurseArrs));
1468 msSafe[ii+1]=(*it)->getMesh()->buildCellFieldOnRecurseWithoutOverlapWithoutGhost(ghostSz,tmpArrs);
1470 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsOff(DataArrayInt::BuildListOfSwitchedOff(bs));
1472 MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret(MEDCouplingFieldDouble::New(ON_CELLS));
1473 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr2(extractGhostFrom(ghostSz,recurseArrs[0]));
1474 arr2=arr2->selectByTupleIdSafe(eltsOff->begin(),eltsOff->end());
1475 ret->setArray(arr2);
1476 ret->setName(arr2->getName());
1477 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> part(_mesh->buildUnstructured());
1478 MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> mesh(part->buildPartOfMySelf(eltsOff->begin(),eltsOff->end(),false));
1482 std::vector< const MEDCouplingFieldDouble * > ms(msSafe.size());
1483 for(std::size_t i=0;i<msSafe.size();i++)
1486 return MEDCouplingFieldDouble::MergeFields(ms);
1490 * This method extracts from \arr arr the part inside \a arr by cutting the \a ghostSz external part.
1491 * \arr is expected to be an array having a number of tuples equal to \c getImageMesh()->buildWithGhost(ghostSz).
1493 DataArrayDouble *MEDCouplingCartesianAMRMeshGen::extractGhostFrom(int ghostSz, const DataArrayDouble *arr) const
1495 std::vector<int> st(_mesh->getCellGridStructure());
1496 std::vector< std::pair<int,int> > p(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(st));
1497 std::transform(st.begin(),st.end(),st.begin(),std::bind2nd(std::plus<int>(),2*ghostSz));
1498 MEDCouplingStructuredMesh::ApplyGhostOnCompactFrmt(p,ghostSz);
1499 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(MEDCouplingStructuredMesh::ExtractFieldOfDoubleFrom(st,arr,p));
1504 * This method returns all the patches in \a this not equal to \a patchId that are in neighborhood of patch with id \a patchId.
1506 * \sa fillCellFieldOnPatchOnlyGhostAdv
1508 std::vector<int> MEDCouplingCartesianAMRMeshGen::getPatchIdsInTheNeighborhoodOf(int patchId, int ghostLev) const
1510 std::vector<int> ret;
1511 int nbp(getNumberOfPatches());
1513 for(int i=0;i<nbp;i++)
1516 if(isPatchInNeighborhoodOf(i,patchId,ghostLev))
1523 * This method returns a dump python of \a this. It is useful for users of createPatchesFromCriterion method for debugging.
1525 * \sa dumpPatchesOf, createPatchesFromCriterion, createPatchesFromCriterionML
1527 std::string MEDCouplingCartesianAMRMeshGen::buildPythonDumpOfThis() const
1529 std::ostringstream oss;
1530 oss << "amr=MEDCouplingCartesianAMRMesh(\""<< getImageMesh()->getName() << "\"," << getSpaceDimension() << ",[";
1531 std::vector<int> ngs(getImageMesh()->getNodeGridStructure());
1532 std::vector<double> orig(getImageMesh()->getOrigin()),dxyz(getImageMesh()->getDXYZ());
1533 std::copy(ngs.begin(),ngs.end(),std::ostream_iterator<int>(oss,","));
1535 std::copy(orig.begin(),orig.end(),std::ostream_iterator<double>(oss,","));
1537 std::copy(dxyz.begin(),dxyz.end(),std::ostream_iterator<double>(oss,","));
1539 dumpPatchesOf("amr",oss);
1543 MEDCouplingCartesianAMRMeshGen::MEDCouplingCartesianAMRMeshGen(const MEDCouplingCartesianAMRMeshGen& other):RefCountObject(other),_mesh(other._mesh),_patches(other._patches),_factors(other._factors)
1545 const MEDCouplingIMesh *mesh(other._mesh);
1547 _mesh=static_cast<MEDCouplingIMesh *>(mesh->deepCpy());
1548 std::size_t sz(other._patches.size());
1549 for(std::size_t i=0;i<sz;i++)
1551 const MEDCouplingCartesianAMRPatch *patch(other._patches[i]);
1553 _patches[i]=patch->deepCpy(this);
1557 MEDCouplingCartesianAMRMeshGen::MEDCouplingCartesianAMRMeshGen(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1558 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop)
1560 _mesh=MEDCouplingIMesh::New(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop);
1563 MEDCouplingCartesianAMRMeshGen::MEDCouplingCartesianAMRMeshGen(MEDCouplingIMesh *mesh)
1566 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen(MEDCouplingIMesh *mesh) constructor : The input mesh is null !");
1567 mesh->checkCoherency();
1568 _mesh=mesh; _mesh->incrRef();
1571 void MEDCouplingCartesianAMRMeshGen::checkPatchId(int patchId) const
1573 int sz(getNumberOfPatches());
1574 if(patchId<0 || patchId>=sz)
1576 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMeshGen::checkPatchId : invalid patchId (" << patchId << ") ! Must be in [0," << sz << ") !";
1577 throw INTERP_KERNEL::Exception(oss.str().c_str());
1581 void MEDCouplingCartesianAMRMeshGen::checkFactorsAndIfNotSetAssign(const std::vector<int>& factors)
1583 if(getSpaceDimension()!=(int)factors.size())
1584 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::checkFactorsAndIfNotSetAssign : invalid size of factors ! size must be equal to the spaceDimension !");
1585 if(_factors.empty())
1591 if(_factors!=factors)
1592 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::checkFactorsAndIfNotSetAssign : the factors ");
1596 void MEDCouplingCartesianAMRMeshGen::retrieveGridsAtInternal(int lev, std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGen> >& grids) const
1600 const MEDCouplingCartesianAMRMesh *thisc(dynamic_cast<const MEDCouplingCartesianAMRMesh *>(this));//tony
1601 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGF> elt(new MEDCouplingCartesianAMRPatchGF(const_cast<MEDCouplingCartesianAMRMesh *>(thisc)));
1602 grids.push_back(DynamicCastSafe<MEDCouplingCartesianAMRPatchGF,MEDCouplingCartesianAMRPatchGen>(elt));
1606 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1608 const MEDCouplingCartesianAMRPatch *pt(*it);
1611 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> tmp1(*it);
1612 grids.push_back(DynamicCastSafe<MEDCouplingCartesianAMRPatch,MEDCouplingCartesianAMRPatchGen>(tmp1));
1618 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1620 const MEDCouplingCartesianAMRPatch *pt(*it);
1622 pt->getMesh()->retrieveGridsAtInternal(lev-1,grids);
1627 int MEDCouplingCartesianAMRMeshGen::GetGhostLevelInFineRef(int ghostLev, const std::vector<int>& factors)
1630 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::GetGhostLevelInFineRef : the ghost size must be >=0 !");
1632 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::GetGhostLevelInFineRef : no factors defined !");
1633 int ghostLevInPatchRef;
1635 ghostLevInPatchRef=0;
1638 ghostLevInPatchRef=(ghostLev-1)/factors[0]+1;
1639 for(std::size_t i=0;i<factors.size();i++)
1640 ghostLevInPatchRef=std::max(ghostLevInPatchRef,(ghostLev-1)/factors[i]+1);
1642 return ghostLevInPatchRef;
1646 * This method returns a sub set of \a all. The subset is defined by the \a head in the tree defined by \a this.
1647 * Elements in \a all are expected to be sorted from god father to most refined structure.
1649 std::vector<const DataArrayDouble *> MEDCouplingCartesianAMRMeshGen::extractSubTreeFromGlobalFlatten(const MEDCouplingCartesianAMRMeshGen *head, const std::vector<const DataArrayDouble *>& all) const
1651 int maxLev(getMaxNumberOfLevelsRelativeToThis());
1652 std::vector<const DataArrayDouble *> ret;
1653 std::vector<const MEDCouplingCartesianAMRMeshGen *> meshes(1,this);
1655 for(int i=0;i<maxLev;i++)
1657 std::vector<const MEDCouplingCartesianAMRMeshGen *> meshesTmp;
1658 for(std::vector<const MEDCouplingCartesianAMRMeshGen *>::const_iterator it=meshes.begin();it!=meshes.end();it++)
1660 if((*it)==head || head->isObjectInTheProgeny(*it))
1661 ret.push_back(all[kk]);
1663 std::vector< const MEDCouplingCartesianAMRPatch *> ps((*it)->getPatches());
1664 for(std::vector< const MEDCouplingCartesianAMRPatch *>::const_iterator it0=ps.begin();it0!=ps.end();it0++)
1666 const MEDCouplingCartesianAMRMeshGen *mesh((*it0)->getMesh());
1667 meshesTmp.push_back(mesh);
1673 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::extractSubTreeFromGlobalFlatten : the size of input vector is not compatible with number of leaves in this !");
1677 void MEDCouplingCartesianAMRMeshGen::dumpPatchesOf(const std::string& varName, std::ostream& oss) const
1680 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1682 const MEDCouplingCartesianAMRPatch *patch(*it);
1685 std::ostringstream oss2; oss2 << varName << ".addPatch([";
1686 const std::vector< std::pair<int,int> >& bltr(patch->getBLTRRange());
1687 std::size_t sz(bltr.size());
1688 for(std::size_t i=0;i<sz;i++)
1690 oss2 << "(" << bltr[i].first << "," << bltr[i].second << ")";
1695 std::copy(_factors.begin(),_factors.end(),std::ostream_iterator<int>(oss2,","));
1698 std::ostringstream oss3; oss3 << varName << "[" << j++ << "]";
1699 patch->getMesh()->dumpPatchesOf(oss3.str(),oss);
1704 std::size_t MEDCouplingCartesianAMRMeshGen::getHeapMemorySizeWithoutChildren() const
1706 return sizeof(MEDCouplingCartesianAMRMeshGen);
1709 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRMeshGen::getDirectChildren() const
1711 std::vector<const BigMemoryObject *> ret;
1712 if((const MEDCouplingIMesh *)_mesh)
1713 ret.push_back((const MEDCouplingIMesh *)_mesh);
1714 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1716 if((const MEDCouplingCartesianAMRPatch*)*it)
1717 ret.push_back((const MEDCouplingCartesianAMRPatch*)*it);
1722 void MEDCouplingCartesianAMRMeshGen::updateTime() const
1724 if((const MEDCouplingIMesh *)_mesh)
1725 updateTimeWith(*_mesh);
1726 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1728 const MEDCouplingCartesianAMRPatch *elt(*it);
1731 const MEDCouplingCartesianAMRMeshGen *mesh(elt->getMesh());
1733 updateTimeWith(*mesh);
1737 MEDCouplingCartesianAMRMeshSub::MEDCouplingCartesianAMRMeshSub(MEDCouplingCartesianAMRMeshGen *father, MEDCouplingIMesh *mesh):MEDCouplingCartesianAMRMeshGen(mesh),_father(father)
1740 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshSub(MEDCouplingCartesianAMRMeshGen *father, MEDCouplingIMesh *mesh) constructor : empty father !");
1743 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMeshSub::getFather() const
1748 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMeshSub::getGodFather() const
1751 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshSub::getGodFather : Impossible to find a god father because there is a hole in chain !");
1752 return _father->getGodFather();
1756 * This method returns the level of \a this. 0 for god father. 1 for children of god father ...
1758 int MEDCouplingCartesianAMRMeshSub::getAbsoluteLevel() const
1761 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshSub::getAbsoluteLevel : Impossible to find a god father because there is a hole in chain !");
1762 return _father->getAbsoluteLevel()+1;
1765 void MEDCouplingCartesianAMRMeshSub::detachFromFather()
1771 std::vector< std::pair<int,int> > MEDCouplingCartesianAMRMeshSub::positionRelativeToGodFather(std::vector<int>& st) const
1773 st=_father->getFactors();
1774 std::size_t dim(st.size());
1775 std::vector<int> prev(st);
1776 int id(_father->getPatchIdFromChildMesh(this));
1777 const MEDCouplingCartesianAMRPatch *p(_father->getPatch(id));
1778 std::vector< std::pair<int,int> > ret(p->getBLTRRange());
1779 std::vector<int> delta(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(ret)),start(dim);
1780 std::transform(delta.begin(),delta.end(),prev.begin(),delta.begin(),std::multiplies<int>());
1781 for(std::size_t i=0;i<dim;i++)
1782 start[i]=ret[i].first;
1783 std::transform(start.begin(),start.end(),prev.begin(),start.begin(),std::multiplies<int>());
1784 const MEDCouplingCartesianAMRMeshGen *it(_father);
1785 while(!dynamic_cast<const MEDCouplingCartesianAMRMesh *>(it))
1787 const MEDCouplingCartesianAMRMeshSub *itc(static_cast<const MEDCouplingCartesianAMRMeshSub *>(it));
1788 int id2(itc->_father->getPatchIdFromChildMesh(itc));
1789 const MEDCouplingCartesianAMRPatch *p2(itc->_father->getPatch(id2));
1790 const std::vector< std::pair<int,int> >& start2(p2->getBLTRRange());
1791 std::vector<int> tmp(dim);
1792 for(std::size_t i=0;i<dim;i++)
1793 tmp[i]=start2[i].first;
1795 prev=itc->_father->getFactors();
1796 std::transform(st.begin(),st.end(),prev.begin(),st.begin(),std::multiplies<int>());
1797 std::transform(st.begin(),st.end(),tmp.begin(),tmp.begin(),std::multiplies<int>());
1798 std::transform(start.begin(),start.end(),tmp.begin(),start.begin(),std::plus<int>());
1801 for(std::size_t i=0;i<dim;i++)
1803 ret[i].first=start[i];
1804 ret[i].second=start[i]+delta[i];
1809 int MEDCouplingCartesianAMRMeshSub::getAbsoluteLevelRelativeTo(const MEDCouplingCartesianAMRMeshGen *ref) const
1814 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshSub::getAbsoluteLevelRelativeTo : ref is not in the progeny of this !");
1816 return _father->getAbsoluteLevelRelativeTo(ref)+1;
1819 MEDCouplingCartesianAMRMeshSub::MEDCouplingCartesianAMRMeshSub(const MEDCouplingCartesianAMRMeshSub& other, MEDCouplingCartesianAMRMeshGen *father):MEDCouplingCartesianAMRMeshGen(other),_father(father)
1823 MEDCouplingCartesianAMRMeshSub *MEDCouplingCartesianAMRMeshSub::deepCpy(MEDCouplingCartesianAMRMeshGen *fath) const
1825 return new MEDCouplingCartesianAMRMeshSub(*this,fath);
1829 * \sa getPositionRelativeTo
1831 void MEDCouplingCartesianAMRMeshSub::getPositionRelativeToInternal(const MEDCouplingCartesianAMRMeshGen *ref, std::vector<int>& ret) const
1836 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshSub::getPositionRelativeToInternal : ref is not in the progeny of this !");
1837 int myId(_father->getPatchIdFromChildMesh(this));
1838 ret.push_back(myId);
1839 _father->getPositionRelativeToInternal(ref,ret);
1842 MEDCouplingCartesianAMRMesh *MEDCouplingCartesianAMRMesh::New(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1843 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop)
1845 return new MEDCouplingCartesianAMRMesh(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop);
1848 MEDCouplingCartesianAMRMesh *MEDCouplingCartesianAMRMesh::New(MEDCouplingIMesh *mesh)
1850 return new MEDCouplingCartesianAMRMesh(mesh);
1853 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMesh::getFather() const
1855 //I'm god father ! No father !
1859 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMesh::getGodFather() const
1864 int MEDCouplingCartesianAMRMesh::getAbsoluteLevel() const
1869 void MEDCouplingCartesianAMRMesh::detachFromFather()
1870 {//not a bug - do nothing
1873 std::vector< std::pair<int,int> > MEDCouplingCartesianAMRMesh::positionRelativeToGodFather(std::vector<int>& st) const
1875 st=_mesh->getCellGridStructure();
1876 return MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(st);
1879 int MEDCouplingCartesianAMRMesh::getAbsoluteLevelRelativeTo(const MEDCouplingCartesianAMRMeshGen *ref) const
1883 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::getAbsoluteLevelRelativeTo : ref is not in the progeny of this !");
1886 std::vector<MEDCouplingCartesianAMRPatchGen *> MEDCouplingCartesianAMRMesh::retrieveGridsAt(int absoluteLev) const
1888 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGen> > rets;
1889 retrieveGridsAtInternal(absoluteLev,rets);
1890 std::vector< MEDCouplingCartesianAMRPatchGen * > ret(rets.size());
1891 for(std::size_t i=0;i<rets.size();i++)
1893 ret[i]=rets[i].retn();
1898 MEDCouplingCartesianAMRMesh *MEDCouplingCartesianAMRMesh::deepCpy(MEDCouplingCartesianAMRMeshGen *father) const
1901 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::deepCpy : specifying a not null father for a God Father object !");
1902 return new MEDCouplingCartesianAMRMesh(*this);
1906 * This method creates a multi level patches split at once.
1907 * This method calls as times as size of \a bso createPatchesFromCriterion. Size of \a bso and size of \a factors must be the same !
1908 * \b WARNING, after the call the number of levels in \a this is equal to bso.size() + 1 !
1911 * \param [in] criterion
1912 * \param [in] factors
1913 * \param [in] eps - See DataArrayDouble::toVectorOfBool for more information about the semantic of eps.
1915 * \sa createPatchesFromCriterion
1917 void MEDCouplingCartesianAMRMesh::createPatchesFromCriterionML(const std::vector<const INTERP_KERNEL::BoxSplittingOptions *>& bso, const DataArrayDouble *criterion, const std::vector< std::vector<int> >& factors, double eps)
1919 std::size_t nbOfLevs(bso.size());
1920 if(nbOfLevs!=factors.size())
1921 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createPatchesFromCriterionML : size of vectors must be the same !");
1925 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createPatchesFromCriterionML : pointers in 1st arg must be not NULL !");
1926 createPatchesFromCriterion(*bso[0],criterion,factors[0],eps);
1927 for(std::size_t i=1;i<nbOfLevs;i++)
1930 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createPatchesFromCriterionML : presence of a NULL BoxSplittingOptions in input vector !");
1932 std::vector<MEDCouplingCartesianAMRPatchGen *> elts(retrieveGridsAt((int)(i)));
1933 std::size_t sz(elts.size());
1934 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGen> > elts2(sz);
1935 std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> > elts3(sz);
1936 for(std::size_t ii=0;ii<sz;ii++)
1939 static const char TMP_STR[]="TMP";
1940 std::vector< std::pair<std::string,int> > fieldNames(1); fieldNames[0].first=TMP_STR; fieldNames[0].second=1;
1941 MEDCouplingAutoRefCountObjectPtr<MEDCouplingAMRAttribute> att(MEDCouplingAMRAttribute::New(this,fieldNames,0));
1943 DataArrayDouble *tmpDa(const_cast<DataArrayDouble *>(att->getFieldOn(this,TMP_STR)));
1944 tmpDa->cpyFrom(*criterion);
1945 att->synchronizeCoarseToFine();
1946 for(std::size_t ii=0;ii<sz;ii++)
1948 const DataArrayDouble *critOnLeaf(att->getFieldOn(const_cast<MEDCouplingCartesianAMRMeshGen *>(elts[ii]->getMesh()),TMP_STR));
1949 elts3[ii]=const_cast<DataArrayDouble *>(critOnLeaf); elts3[ii]->incrRef();
1952 for(std::size_t ii=0;ii<sz;ii++)
1953 const_cast<MEDCouplingCartesianAMRMeshGen *>(elts[ii]->getMesh())->createPatchesFromCriterion(*bso[i],elts3[ii],factors[i],eps);
1957 void MEDCouplingCartesianAMRMesh::getPositionRelativeToInternal(const MEDCouplingCartesianAMRMeshGen *ref, std::vector<int>& ret) const
1962 MEDCouplingCartesianAMRMesh::MEDCouplingCartesianAMRMesh(const MEDCouplingCartesianAMRMesh& other):MEDCouplingCartesianAMRMeshGen(other)
1966 MEDCouplingCartesianAMRMesh::MEDCouplingCartesianAMRMesh(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1967 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop):MEDCouplingCartesianAMRMeshGen(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop)
1971 MEDCouplingCartesianAMRMesh::MEDCouplingCartesianAMRMesh(MEDCouplingIMesh *mesh):MEDCouplingCartesianAMRMeshGen(mesh)
1975 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRMesh::getDirectChildren() const
1977 std::vector<const BigMemoryObject *> ret(MEDCouplingCartesianAMRMeshGen::getDirectChildren());
1981 MEDCouplingCartesianAMRMesh::~MEDCouplingCartesianAMRMesh()