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 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatchGen::getMeshSafe() const
60 const MEDCouplingCartesianAMRMeshGen *mesh(_mesh);
62 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatchGen::getMeshSafe const : the mesh is NULL !");
66 MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatchGen::getMeshSafe()
68 MEDCouplingCartesianAMRMeshGen *mesh(_mesh);
70 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatchGen::getMeshSafe : the mesh is NULL !");
74 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRPatchGen::getDirectChildren() const
76 std::vector<const BigMemoryObject *> ret;
77 if((const MEDCouplingCartesianAMRMeshGen *)_mesh)
78 ret.push_back((const MEDCouplingCartesianAMRMeshGen *)_mesh);
83 * \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.
84 * \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,
85 * a the end cell (\b excluded) of the range for the second element of the pair.
87 MEDCouplingCartesianAMRPatch::MEDCouplingCartesianAMRPatch(MEDCouplingCartesianAMRMeshGen *mesh, const std::vector< std::pair<int,int> >& bottomLeftTopRight):MEDCouplingCartesianAMRPatchGen(mesh),_bl_tr(bottomLeftTopRight)
89 int dim((int)bottomLeftTopRight.size()),dimExp(_mesh->getSpaceDimension());
91 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch constructor : space dimension of father and input bottomLeft/topRight size mismatches !");
94 void MEDCouplingCartesianAMRPatch::addPatch(const std::vector< std::pair<int,int> >& bottomLeftTopRight, const std::vector<int>& factors)
96 return getMeshSafe()->addPatch(bottomLeftTopRight,factors);
99 int MEDCouplingCartesianAMRPatch::getNumberOfOverlapedCellsForFather() const
101 return MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(_bl_tr);
105 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
106 * 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 !).
107 * Call isInMyNeighborhoodExt to deal with 2 patches not sharing directly the same father.
109 * \param [in] other - The other patch
110 * \param [in] ghostLev - The size of the neighborhood zone.
112 * \throw if \a this or \a other are invalid (end before start).
113 * \throw if \a ghostLev is \b not >= 0 .
114 * \throw if \a this and \a other have not the same space dimension.
116 * \sa isInMyNeighborhoodExt
118 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhood(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
121 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the size of the neighborhood must be >= 0 !");
123 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the input patch is NULL !");
124 const std::vector< std::pair<int,int> >& thisp(getBLTRRange());
125 const std::vector< std::pair<int,int> >& otherp(other->getBLTRRange());
126 return IsInMyNeighborhood(ghostLev,thisp,otherp);
130 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
131 * 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
132 * ancestor must be the same. If they don't have the same ancestor an exception will be thrown.
134 * \param [in] other - The other patch
135 * \param [in] ghostLev - The size of the neighborhood zone.
137 * \throw if \a this or \a other are invalid (end before start).
138 * \throw if \a ghostLev is \b not >= 0 .
139 * \throw if \a this and \a other have not the same space dimension.
140 * \throw if there is not common ancestor of \a this and \a other.
142 * \sa isInMyNeighborhood, isInMyNeighborhoodDiffLev
144 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
147 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt : the size of the neighborhood must be >= 0 !");
149 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt : the input patch is NULL !");
151 const MEDCouplingCartesianAMRMeshGen *com(FindCommonAncestor(this,other,lev));//check that factors are OK
153 return isInMyNeighborhood(other,ghostLev);
154 std::vector<int> offset(ComputeOffsetFromTwoToOne(com,lev,this,other));
155 const std::vector< std::pair<int,int> >& thisp(getBLTRRange());
156 std::vector< std::pair<int,int> > otherp(other->getBLTRRange());
157 otherp=MEDCouplingStructuredMesh::TranslateCompactFrmt(otherp,offset);
158 return IsInMyNeighborhood(ghostLev,thisp,otherp);
162 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
163 * the must be >= 0. This method works even if \a this and \a other does not share the same father.
164 * 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.
166 * \param [in] other - The other patch
167 * \param [in] ghostLev - The size of the neighborhood zone.
169 * \throw if \a this or \a other are invalid (end before start).
170 * \throw if \a ghostLev is \b not >= 0 .
171 * \throw if \a this and \a other have not the same space dimension.
172 * \throw if there is not common ancestor of \a this and \a other.
174 * \sa isInMyNeighborhoodExt
176 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhoodDiffLev(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
178 std::vector< std::pair<int,int> > thispp,otherpp;
179 std::vector<int> factors;
180 ComputeZonesOfTwoRelativeToOneDiffLev(ghostLev,this,other,thispp,otherpp,factors);
181 return IsInMyNeighborhood(ghostLev,thispp,otherpp);
184 std::vector<int> MEDCouplingCartesianAMRPatch::computeCellGridSt() const
186 const MEDCouplingCartesianAMRMeshGen *m(getMesh());
188 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::computeCellGridSt : no mesh held by this !");
189 const MEDCouplingCartesianAMRMeshGen *father(m->getFather());
191 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::computeCellGridSt : no father help by underlying mesh !");
192 const std::vector< std::pair<int,int> >& bltr(getBLTRRange());
193 const std::vector<int>& factors(father->getFactors());
194 std::vector<int> ret(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(bltr));
195 std::transform(ret.begin(),ret.end(),factors.begin(),ret.begin(),std::multiplies<int>());
199 bool MEDCouplingCartesianAMRPatch::IsInMyNeighborhood(int ghostLev, const std::vector< std::pair<int,int> >& p1, const std::vector< std::pair<int,int> >& p2)
201 std::size_t thispsize(p1.size());
202 if(thispsize!=p2.size())
203 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the dimensions must be the same !");
204 for(std::size_t i=0;i<thispsize;i++)
206 const std::pair<int,int>& thispp(p1[i]);
207 const std::pair<int,int>& otherpp(p2[i]);
208 if(thispp.second<thispp.first)
209 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
210 if(otherpp.second<otherpp.first)
211 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
212 if(otherpp.first==thispp.second+ghostLev-1)
214 if(otherpp.second+ghostLev-1==thispp.first)
216 int start(std::max(thispp.first,otherpp.first)),end(std::min(thispp.second,otherpp.second));
224 * \sa FindNeighborsOfSubPatchesOf
226 std::vector< std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > > MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOfSameLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2)
229 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOfSameLev : the input pointers must be not NULL !");
230 std::vector< std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > > ret;
231 std::vector< const MEDCouplingCartesianAMRPatch *> p1Work(p1->getMesh()->getPatches()),p2Work(p2->getMesh()->getPatches());
232 while(!p1Work.empty())
234 std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > retTmp;
235 std::vector<const MEDCouplingCartesianAMRPatch *> p1Work2,p2Work2;
236 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it1=p1Work.begin();it1!=p1Work.end();it1++)
238 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it2=p2Work.begin();it2!=p2Work.end();it2++)
240 if((*it1)->isInMyNeighborhoodExt(*it2,ghostLev))
241 retTmp.push_back(std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *>(*it1,*it2));
243 std::vector<const MEDCouplingCartesianAMRPatch *> tmp1((*it1)->getMesh()->getPatches());
244 p1Work2.insert(p1Work2.end(),tmp1.begin(),tmp1.end());
246 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it2=p2Work.begin();it2!=p2Work.end();it2++)
248 std::vector<const MEDCouplingCartesianAMRPatch *> tmp2((*it2)->getMesh()->getPatches());
249 p2Work2.insert(p2Work2.end(),tmp2.begin(),tmp2.end());
251 ret.push_back(retTmp);
259 * This method returns all pair of patches (pa,pb) so that pb is in the neighborhood of pa (size of neighborhood is \a ghostLev).
260 * 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
261 * FindNeighborsOfSubPatchesOfSameLev.
263 * \sa FindNeighborsOfSubPatchesOfSameLev
265 void MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOf(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, std::vector< std::pair<const MEDCouplingCartesianAMRPatch *, const MEDCouplingCartesianAMRPatch *> >& ret)
268 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOf : the input pointers must be not NULL !");
269 std::vector< const MEDCouplingCartesianAMRPatch *> p1Work(p1->getMesh()->getPatches());
270 while(!p1Work.empty())
272 std::vector<const MEDCouplingCartesianAMRPatch *> p1Work2;
273 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it0=p1Work.begin();it0!=p1Work.end();it0++)
275 if((*it0)->isInMyNeighborhoodDiffLev(p2,ghostLev))
276 ret.push_back(std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *>(*it0,p2));
277 std::vector<const MEDCouplingCartesianAMRPatch *> tmp2((*it0)->getMesh()->getPatches());
278 p1Work2.insert(p1Work2.end(),tmp2.begin(),tmp2.end());
285 * \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.
287 * \saUpdateNeighborsOfOneWithTwoExt
289 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(int ghostLev, const std::vector<int>& factors, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
291 const std::vector< std::pair<int,int> >& p1BLTR(p1->getBLTRRange());
292 const std::vector< std::pair<int,int> >& p2BLTR(p2->getBLTRRange());
293 UpdateNeighborsOfOneWithTwoInternal(ghostLev,factors,p1BLTR,p2BLTR,dataOnP1,dataOnP2);
297 * Idem than UpdateNeighborsOfOneWithTwo, except that here \a p1 and \a p2 are not sharing the same direct father.
299 * \sa UpdateNeighborsOfOneWithTwo
301 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoExt(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
303 const std::vector< std::pair<int,int> >& p1BLTR(p1->getBLTRRange());//p1BLTR=[(10,12),(5,8)]
304 std::vector< std::pair<int,int> > p2BLTR(p2->getBLTRRange());//p2BLTR=[(0,1),(0,5)]
306 const MEDCouplingCartesianAMRMeshGen *ca(FindCommonAncestor(p1,p2,lev));
307 std::vector<int> offset(ComputeOffsetFromTwoToOne(ca,lev,p1,p2));//[12,4]
308 p2BLTR=MEDCouplingStructuredMesh::TranslateCompactFrmt(p2BLTR,offset);//p2BLTR=[(12,13),(4,9)]
309 UpdateNeighborsOfOneWithTwoInternal(ghostLev,p1->getMesh()->getFather()->getFactors(),p1BLTR,p2BLTR,dataOnP1,dataOnP2);
313 * \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 !
315 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoMixedLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
317 std::vector< std::pair<int,int> > p1pp,p2pp;
318 std::vector<int> factors;
319 ComputeZonesOfTwoRelativeToOneDiffLev(ghostLev,p1,p2,p1pp,p2pp,factors);
321 std::vector<int> dimsP2NotRefined(p2->computeCellGridSt());
322 std::vector<int> dimsP2Refined(dimsP2NotRefined);
323 std::transform(dimsP2NotRefined.begin(),dimsP2NotRefined.end(),factors.begin(),dimsP2Refined.begin(),std::multiplies<int>());
324 std::vector< std::pair<int,int> > p2RefinedAbs(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(dimsP2NotRefined));
325 std::vector<int> dimsP2RefinedGhost(dimsP2Refined.size());
326 std::transform(dimsP2Refined.begin(),dimsP2Refined.end(),dimsP2RefinedGhost.begin(),std::bind2nd(std::plus<int>(),2*ghostLev));
327 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> fineP2(DataArrayDouble::New()); fineP2->alloc(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(dimsP2RefinedGhost),dataOnP2->getNumberOfComponents());
328 MEDCouplingIMesh::SpreadCoarseToFineGhost(dataOnP2,dimsP2NotRefined,fineP2,p2RefinedAbs,factors,ghostLev);
330 UpdateNeighborsOfOneWithTwoInternal(ghostLev,p1->getMesh()->getFather()->getFactors(),p1pp,p2pp,dataOnP1,fineP2);
334 * \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 !
335 * 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
336 * on the same level as \a p1.
338 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)
340 std::vector<const MEDCouplingCartesianAMRMeshGen *> ancestorsOfThis;
341 const MEDCouplingCartesianAMRMeshGen *work(p1->getMesh()),*work2(0);
342 ancestorsOfThis.push_back(work);
345 work=work->getFather();
347 ancestorsOfThis.push_back(work);
352 std::size_t levThis(0),levOther(0);
353 while(work && !found)
356 work=work->getFather();
360 std::vector<const MEDCouplingCartesianAMRMeshGen *>::iterator it(std::find(ancestorsOfThis.begin(),ancestorsOfThis.end(),work));
361 if(it!=ancestorsOfThis.end())
363 levThis=std::distance(ancestorsOfThis.begin(),it);
369 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev : no common ancestor found !");
370 if(levThis<=levOther)
371 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev : this method is not called correctly !");
373 const MEDCouplingCartesianAMRMeshGen *comAncestor(ancestorsOfThis[levThis]);
374 int idThis(comAncestor->getPatchIdFromChildMesh(ancestorsOfThis[levThis-1])),idOther(comAncestor->getPatchIdFromChildMesh(work2));
375 const MEDCouplingCartesianAMRPatch *thisp(comAncestor->getPatch(idThis)),*otherp(comAncestor->getPatch(idOther));
376 std::vector<int> offset(ComputeOffsetFromTwoToOne(comAncestor,levOther,thisp,otherp));
377 p1Zone=thisp->getBLTRRange(); p2Zone=MEDCouplingStructuredMesh::TranslateCompactFrmt(otherp->getBLTRRange(),offset);
378 factToApplyOn2.resize(p1Zone.size()); std::fill(factToApplyOn2.begin(),factToApplyOn2.end(),1);
380 std::size_t nbOfTurn(levThis-levOther);
381 for(std::size_t i=0;i<nbOfTurn;i++)
383 std::vector< std::pair<int,int> > tmp0;
384 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p1Zone,p2Zone,tmp0,false);
386 const MEDCouplingCartesianAMRMeshGen *curAncestor(ancestorsOfThis[levThis-i]);
387 ApplyFactorsOnCompactFrmt(p2Zone,curAncestor->getFactors());
388 curAncestor=ancestorsOfThis[levThis-1-i];
389 const std::vector<int>& factors(curAncestor->getFactors());
390 std::transform(factToApplyOn2.begin(),factToApplyOn2.end(),factors.begin(),factToApplyOn2.begin(),std::multiplies<int>());
391 int tmpId(curAncestor->getPatchIdFromChildMesh(ancestorsOfThis[levThis-2-i]));
392 p1Zone=curAncestor->getPatch(tmpId)->getBLTRRange();
396 std::size_t MEDCouplingCartesianAMRPatch::getHeapMemorySizeWithoutChildren() const
398 std::size_t ret(sizeof(MEDCouplingCartesianAMRPatch));
399 ret+=_bl_tr.capacity()*sizeof(std::pair<int,int>);
403 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatch::FindCommonAncestor(const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, int& lev)
405 const MEDCouplingCartesianAMRMeshGen *f1(p1->_mesh),*f2(p2->_mesh);
407 while(f1!=f2 || f1==0 || f2==0)
409 f1=f1->getFather(); f2=f2->getFather();
410 if(f1->getFactors()!=f2->getFactors())
411 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindCommonAncestor : factors differ !");
415 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindCommonAncestor : no common ancestor between p1 and p2 !");
419 std::vector<int> MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne(const MEDCouplingCartesianAMRMeshGen *comAncestor, int lev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2)
422 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne : this method is useful only for lev > 0 !");
424 int dim(p1->getMesh()->getSpaceDimension());
425 if(p2->getMesh()->getSpaceDimension()!=dim)
426 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne : dimension must be the same !");
427 std::vector< int > ret(dim,0);
428 for(int i=0;i<zeLev;i++)
430 const MEDCouplingCartesianAMRMeshGen *f1(p1->_mesh),*f2(p2->_mesh);
431 const MEDCouplingCartesianAMRPatch *p1h(0),*p2h(0);
432 for(int j=0;j<lev-i;j++)
434 const MEDCouplingCartesianAMRMeshGen *f1tmp(f1->getFather()),*f2tmp(f2->getFather());
435 int pid1(f1tmp->getPatchIdFromChildMesh(f1)),pid2(f2tmp->getPatchIdFromChildMesh(f2));
436 p1h=f1tmp->getPatch(pid1); p2h=f2tmp->getPatch(pid2);
439 std::vector< std::pair<int,int> > p2c(p2h->getBLTRRange());
440 for(int k=0;k<dim;k++)
442 p2c[k].first+=ret[k];
443 p2c[k].second+=ret[k];
445 for(int k=0;k<dim;k++)
447 ret[k]=p2c[k].first-p1h->getBLTRRange()[k].first;
448 ret[k]*=f1->getFactors()[k];
454 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)
455 {//p1=[(1,4),(2,4)] p2=[(4,5),(3,4)]
456 int dim((int)factors.size());
457 std::vector<int> dimsCoarse(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(p1));//[3,2]
458 std::transform(dimsCoarse.begin(),dimsCoarse.end(),factors.begin(),dimsCoarse.begin(),std::multiplies<int>());//[12,8]
459 std::transform(dimsCoarse.begin(),dimsCoarse.end(),dimsCoarse.begin(),std::bind2nd(std::plus<int>(),2*ghostLev));//[14,10]
460 std::vector< std::pair<int,int> > rangeCoarse(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(dimsCoarse));//[(0,14),(0,10)]
461 std::vector<int> fakeFactors(dim,1);
463 std::vector< std::pair<int,int> > tmp0,tmp1,tmp2;
464 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p1,p2,tmp0,false);//tmp0=[(3,4),(1,2)]
465 ApplyFactorsOnCompactFrmt(tmp0,factors);//tmp0=[(12,16),(4,8)]
466 ApplyGhostOnCompactFrmt(tmp0,ghostLev);//tmp0=[(13,17),(5,9)]
467 std::vector< std::pair<int,int> > interstRange(MEDCouplingStructuredMesh::IntersectRanges(tmp0,rangeCoarse));//interstRange=[(13,14),(5,9)]
468 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p2,p1,tmp1,false);//tmp1=[(-3,0),(-1,1)]
469 ApplyFactorsOnCompactFrmt(tmp1,factors);//tmp1=[(-12,-4),(-4,0)]
470 MEDCouplingStructuredMesh::ChangeReferenceToGlobalOfCompactFrmt(tmp1,interstRange,tmp2,false);//tmp2=[(1,2),(1,5)]
472 std::vector< std::pair<int,int> > dimsFine(p2);
473 ApplyFactorsOnCompactFrmt(dimsFine,factors);
474 ApplyAllGhostOnCompactFrmt(dimsFine,ghostLev);
476 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ghostVals(MEDCouplingStructuredMesh::ExtractFieldOfDoubleFrom(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(dimsFine),dataOnP2,tmp2));
477 MEDCouplingIMesh::CondenseFineToCoarse(dimsCoarse,ghostVals,interstRange,fakeFactors,dataOnP1);
481 * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in refined reference.
482 * \param [in] factors - the factors per axis.
484 void MEDCouplingCartesianAMRPatch::ApplyFactorsOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, const std::vector<int>& factors)
486 std::size_t sz(factors.size());
487 if(sz!=partBeforeFact.size())
488 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ApplyFactorsOnCompactFrmt : size of input vectors must be the same !");
489 for(std::size_t i=0;i<sz;i++)
491 partBeforeFact[i].first*=factors[i];
492 partBeforeFact[i].second*=factors[i];
497 * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in ghost reference.
498 * \param [in] ghostSize - the ghost size of zone for all axis.
500 void MEDCouplingCartesianAMRPatch::ApplyGhostOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, int ghostSize)
503 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ApplyGhostOnCompactFrmt : ghost size must be >= 0 !");
504 std::size_t sz(partBeforeFact.size());
505 for(std::size_t i=0;i<sz;i++)
507 partBeforeFact[i].first+=ghostSize;
508 partBeforeFact[i].second+=ghostSize;
513 * This method is different than ApplyGhostOnCompactFrmt. The \a partBeforeFact parameter is enlarger contrary to ApplyGhostOnCompactFrmt.
515 * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in ghost reference.
516 * \param [in] ghostSize - the ghost size of zone for all axis.
518 void MEDCouplingCartesianAMRPatch::ApplyAllGhostOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, int ghostSize)
521 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ApplyAllGhostOnCompactFrmt : ghost size must be >= 0 !");
522 std::size_t sz(partBeforeFact.size());
523 for(std::size_t i=0;i<sz;i++)
525 partBeforeFact[i].first-=ghostSize;
526 partBeforeFact[i].second+=ghostSize;
530 MEDCouplingCartesianAMRPatchGF::MEDCouplingCartesianAMRPatchGF(MEDCouplingCartesianAMRMesh *mesh):MEDCouplingCartesianAMRPatchGen(mesh)
534 std::size_t MEDCouplingCartesianAMRPatchGF::getHeapMemorySizeWithoutChildren() const
536 return sizeof(MEDCouplingCartesianAMRPatchGF);
539 MEDCouplingDataForGodFather::MEDCouplingDataForGodFather(MEDCouplingCartesianAMRMeshGen *gf):_gf(gf),_tlc(gf)
542 throw INTERP_KERNEL::Exception("MEDCouplingDataForGodFather constructor : A data has to be attached to a AMR Mesh instance !");
546 void MEDCouplingDataForGodFather::checkGodFatherFrozen() const
551 bool MEDCouplingDataForGodFather::changeGodFather(MEDCouplingCartesianAMRMeshGen *gf)
553 bool ret(_tlc.keepTrackOfNewTL(gf));
565 int MEDCouplingCartesianAMRMeshGen::getSpaceDimension() const
567 return _mesh->getSpaceDimension();
570 void MEDCouplingCartesianAMRMeshGen::setFactors(const std::vector<int>& newFactors)
572 if(getSpaceDimension()!=(int)newFactors.size())
573 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::setFactors : size of input factors is not equal to the space dimension !");
579 if(_factors==newFactors)
581 if(!_patches.empty())
582 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::setFactors : modification of factors is not allowed when presence of patches !");
587 int MEDCouplingCartesianAMRMeshGen::getMaxNumberOfLevelsRelativeToThis() const
590 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
591 ret=std::max(ret,(*it)->getMaxNumberOfLevelsRelativeToThis()+1);
596 * This method returns the number of cells of \a this with the help of the MEDCouplingIMesh instance representing \a this.
597 * The patches in \a this are ignored here.
598 * \sa getNumberOfCellsAtCurrentLevelGhost, getNumberOfCellsRecursiveWithOverlap
600 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsAtCurrentLevel() const
602 return _mesh->getNumberOfCells();
606 * 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
607 * to take into account of the ghost cells for future computation.
608 * The patches in \a this are ignored here.
610 * \sa getNumberOfCellsAtCurrentLevel
612 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsAtCurrentLevelGhost(int ghostLev) const
614 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> tmp(_mesh->buildWithGhost(ghostLev));
615 return tmp->getNumberOfCells();
619 * This method returns the number of cells including the current level but \b also \b including recursively all cells of other levels
620 * starting from this. The set of cells which size is returned here are generally overlapping each other.
622 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsRecursiveWithOverlap() const
624 int ret(_mesh->getNumberOfCells());
625 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
627 ret+=(*it)->getNumberOfCellsRecursiveWithOverlap();
633 * This method returns the max number of cells covering all the space without overlapping.
634 * It returns the number of cells of the mesh with the highest resolution.
635 * The returned value is equal to the number of cells of mesh returned by buildUnstructured.
637 * \sa buildUnstructured
639 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsRecursiveWithoutOverlap() const
641 int ret(_mesh->getNumberOfCells());
642 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
644 ret-=(*it)->getNumberOfOverlapedCellsForFather();
645 ret+=(*it)->getNumberOfCellsRecursiveWithoutOverlap();
650 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMeshGen::getFather() const
655 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMeshGen::getGodFather() const
660 return _father->getGodFather();
664 * This method returns the level of \a this. 0 for god father. -1 for children of god father ...
666 int MEDCouplingCartesianAMRMeshGen::getAbsoluteLevel() const
671 return _father->getAbsoluteLevel()-1;
675 * This method returns grids relative to god father to specified level \a absoluteLev.
677 * \return std::vector<MEDCouplingCartesianAMRPatchGen *> - objects in vector are to be managed (decrRef) by the caller.
679 std::vector<MEDCouplingCartesianAMRPatchGen *> MEDCouplingCartesianAMRMeshGen::retrieveGridsAt(int absoluteLev) const
682 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::retrieveGridsAt : absolute level must be >=0 !");
684 return getGodFather()->retrieveGridsAt(absoluteLev);
686 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGen> > rets;
687 retrieveGridsAtInternal(absoluteLev,rets);
688 std::vector< MEDCouplingCartesianAMRPatchGen * > ret(rets.size());
689 for(std::size_t i=0;i<rets.size();i++)
691 ret[i]=rets[i].retn();
696 void MEDCouplingCartesianAMRMeshGen::detachFromFather()
703 * \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,
704 * a the end cell (\b excluded) of the range for the second element of the pair.
705 * \param [in] factors The factor of refinement per axis (different from 0).
707 void MEDCouplingCartesianAMRMeshGen::addPatch(const std::vector< std::pair<int,int> >& bottomLeftTopRight, const std::vector<int>& factors)
709 checkFactorsAndIfNotSetAssign(factors);
710 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> mesh(static_cast<MEDCouplingIMesh *>(_mesh->buildStructuredSubPart(bottomLeftTopRight)));
711 mesh->refineWithFactor(factors);
712 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRMeshSub> zeMesh(new MEDCouplingCartesianAMRMeshSub(this,mesh));
713 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> elt(new MEDCouplingCartesianAMRPatch(zeMesh,bottomLeftTopRight));
714 _patches.push_back(elt);
720 class InternalPatch : public RefCountObjectOnly
723 InternalPatch():_nb_of_true(0) { }
724 int getDimension() const { return (int)_part.size(); }
725 double getEfficiency() const { return (double)_nb_of_true/(double)_crit.size(); }
726 int getNumberOfCells() const { return (int)_crit.size(); }
727 void setNumberOfTrue(int nboft) { _nb_of_true=nboft; }
728 std::vector<bool>& getCriterion() { return _crit; }
729 const std::vector<bool>& getConstCriterion() const { return _crit; }
730 void setPart(const std::vector< std::pair<int,int> >& part) { _part=part; }
731 std::vector< std::pair<int,int> >& getPart() { return _part; }
732 const std::vector< std::pair<int,int> >& getConstPart() const { return _part; }
733 bool presenceOfTrue() const { return _nb_of_true>0; }
734 std::vector<int> computeCGS() const { return MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(_part); }
735 std::vector< std::vector<int> > computeSignature() const { return MEDCouplingStructuredMesh::ComputeSignaturePerAxisOf(computeCGS(),getConstCriterion()); }
736 double getEfficiencyPerAxis(int axisId) const { return (double)_nb_of_true/((double)(_part[axisId].second-_part[axisId].first)); }
737 void zipToFitOnCriterion();
738 void updateNumberOfTrue() const;
739 MEDCouplingAutoRefCountObjectPtr<InternalPatch> extractPart(const std::vector< std::pair<int,int> >&partInGlobal) const;
740 MEDCouplingAutoRefCountObjectPtr<InternalPatch> deepCpy() const;
744 mutable int _nb_of_true;
745 std::vector<bool> _crit;
747 std::vector< std::pair<int,int> > _part;
750 void InternalPatch::zipToFitOnCriterion()
752 std::vector<int> cgs(computeCGS());
753 std::vector<bool> newCrit;
754 std::vector< std::pair<int,int> > newPart,newPart2;
755 int newNbOfTrue(MEDCouplingStructuredMesh::FindMinimalPartOf(cgs,_crit,newCrit,newPart));
756 MEDCouplingStructuredMesh::ChangeReferenceToGlobalOfCompactFrmt(_part,newPart,newPart2);
757 if(newNbOfTrue!=_nb_of_true)
758 throw INTERP_KERNEL::Exception("InternalPatch::zipToFitOnCrit : internal error !");
759 _crit=newCrit; _part=newPart2;
762 void InternalPatch::updateNumberOfTrue() const
764 _nb_of_true=(int)std::count(_crit.begin(),_crit.end(),true);
767 MEDCouplingAutoRefCountObjectPtr<InternalPatch> InternalPatch::extractPart(const std::vector< std::pair<int,int> >&partInGlobal) const
769 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(new InternalPatch);
770 std::vector<int> cgs(computeCGS());
771 std::vector< std::pair<int,int> > newPart;
772 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(_part,partInGlobal,newPart);
773 MEDCouplingStructuredMesh::ExtractFieldOfBoolFrom(cgs,_crit,newPart,ret->getCriterion());
774 ret->setPart(partInGlobal);
775 ret->updateNumberOfTrue();
779 MEDCouplingAutoRefCountObjectPtr<InternalPatch> InternalPatch::deepCpy() const
781 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(new InternalPatch);
786 void DissectBigPatch(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int rangeOfAxisId, bool& cutFound, int& cutPlace)
788 cutFound=false; cutPlace=-1;
789 std::vector<double> ratio(rangeOfAxisId-1);
790 for(int id=0;id<rangeOfAxisId-1;id++)
792 double efficiency[2];
795 std::vector< std::pair<int,int> > rectH(patchToBeSplit->getConstPart());
797 rectH[axisId].second=patchToBeSplit->getConstPart()[axisId].first+id;
799 rectH[axisId].first=patchToBeSplit->getConstPart()[axisId].first+id;
800 MEDCouplingAutoRefCountObjectPtr<InternalPatch> p(patchToBeSplit->deepCpy());
801 p->zipToFitOnCriterion();
803 efficiency[h]=p->getEfficiencyPerAxis(axisId);
805 ratio[id]=std::max(efficiency[0],efficiency[1])/std::min(efficiency[0],efficiency[1]);
807 int minCellDirection(bso.getMinCellDirection()),indexMin(-1);
808 int dimRatioInner(rangeOfAxisId-1-2*(minCellDirection-1));
809 std::vector<double> ratio_inner(dimRatioInner);
810 double minRatio(1.e10);
811 for(int i=0; i<dimRatioInner; i++)
813 if(ratio[minCellDirection-1+i]<minRatio)
815 minRatio=ratio[minCellDirection-1+i];
816 indexMin=i+minCellDirection;
819 cutFound=true; cutPlace=indexMin+patchToBeSplit->getConstPart()[axisId].first-1;
822 void FindHole(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int& axisId, bool& cutFound, int& cutPlace)
824 cutPlace=-1; cutFound=false;
825 int minCellDirection(bso.getMinCellDirection());
826 const int dim(patchToBeSplit->getDimension());
827 std::vector< std::vector<int> > signatures(patchToBeSplit->computeSignature());
828 for(int id=0;id<dim;id++)
830 const std::vector<int>& signature(signatures[id]);
831 std::vector<int> hole;
832 std::vector<double> distance;
833 int len((int)signature.size());
834 for(int i=0;i<len;i++)
836 if(len>= 2*minCellDirection && i >= minCellDirection-1 && i <= len-minCellDirection-1)
840 double center(((double)len/2.));
841 for(std::size_t i=0;i<hole.size();i++)
842 distance.push_back(fabs(hole[i]+1.-center));
844 std::size_t posDistanceMin(std::distance(distance.begin(),std::min_element(distance.begin(),distance.end())));
847 cutPlace=hole[posDistanceMin]+patchToBeSplit->getConstPart()[axisId].first+1;
853 void FindInflection(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, bool& cutFound, int& cutPlace, int& axisId)
855 cutFound=false; cutPlace=-1;// do not set axisId before to be sure that cutFound was set to true
857 const std::vector< std::pair<int,int> >& part(patchToBeSplit->getConstPart());
858 int sign,minCellDirection(bso.getMinCellDirection());
859 const int dim(patchToBeSplit->getDimension());
861 std::vector<int> zeroCrossDims(dim,-1);
862 std::vector<int> zeroCrossVals(dim,-1);
863 std::vector< std::vector<int> > signatures(patchToBeSplit->computeSignature());
864 for (int id=0;id<dim;id++)
866 const std::vector<int>& signature(signatures[id]);
868 std::vector<int> derivate_second_order,gradient_absolute,signe_change,zero_cross,edge,max_cross_list ;
869 std::vector<double> distance ;
871 for (unsigned int i=1;i<signature.size()-1;i++)
872 derivate_second_order.push_back(signature[i-1]-2*signature[i]+signature[i+1]) ;
874 // Gradient absolute value
875 for ( unsigned int i=1;i<derivate_second_order.size();i++)
876 gradient_absolute.push_back(fabs(derivate_second_order[i]-derivate_second_order[i-1])) ;
877 if(derivate_second_order.empty())
879 for (unsigned int i=0;i<derivate_second_order.size()-1;i++)
881 if (derivate_second_order[i]*derivate_second_order[i+1] < 0 )
883 if (derivate_second_order[i]*derivate_second_order[i+1] > 0 )
885 if (derivate_second_order[i]*derivate_second_order[i+1] == 0 )
887 if ( sign==0 || sign==-1 )
888 if ( i >= (unsigned int)minCellDirection-2 && i <= signature.size()-minCellDirection-2 )
890 zero_cross.push_back(i) ;
891 edge.push_back(gradient_absolute[i]) ;
893 signe_change.push_back(sign) ;
895 if ( zero_cross.size() > 0 )
897 int max_cross=*max_element(edge.begin(),edge.end()) ;
898 for (unsigned int i=0;i<edge.size();i++)
899 if (edge[i]==max_cross)
900 max_cross_list.push_back(zero_cross[i]+1) ;
902 double center((signature.size()/2.0));
903 for (unsigned int i=0;i<max_cross_list.size();i++)
904 distance.push_back(fabs(max_cross_list[i]+1-center));
906 float distance_min=*min_element(distance.begin(),distance.end()) ;
907 int pos_distance_min=find(distance.begin(),distance.end(),distance_min)-distance.begin();
908 int best_place = max_cross_list[pos_distance_min] + part[id].first ;
911 zeroCrossDims[id] = best_place ;
912 zeroCrossVals[id] = max_cross ;
915 derivate_second_order.clear() ;
916 gradient_absolute.clear() ;
917 signe_change.clear() ;
920 max_cross_list.clear() ;
924 if ( zeroCrossDims[0]!=-1 || zeroCrossDims[1]!=-1 )
926 int max_cross_dims = *max_element(zeroCrossVals.begin(),zeroCrossVals.end()) ;
928 if (zeroCrossVals[0]==max_cross_dims && zeroCrossVals[1]==max_cross_dims )
930 int nl_left(part[0].second-part[0].first);
931 int nc_left(part[1].second-part[1].first);
932 if ( nl_left >= nc_left )
938 max_cross_dims=std::find(zeroCrossVals.begin(),zeroCrossVals.end(),max_cross_dims)-zeroCrossVals.begin();
940 cutPlace=zeroCrossDims[max_cross_dims];
941 axisId=max_cross_dims ;
945 void TryAction4(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int rangeOfAxisId, bool& cutFound, int& cutPlace)
948 if(patchToBeSplit->getEfficiency()<=bso.getEffeciencySnd())
950 if(rangeOfAxisId>=2*bso.getMinCellDirection())
953 cutPlace=rangeOfAxisId/2+patchToBeSplit->getConstPart()[axisId].first-1;
958 if(patchToBeSplit->getNumberOfCells()>bso.getMaxCells())
960 DissectBigPatch(bso,patchToBeSplit,axisId,rangeOfAxisId,cutFound,cutPlace);
965 MEDCouplingAutoRefCountObjectPtr<InternalPatch> DealWithNoCut(const InternalPatch *patch)
967 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(const_cast<InternalPatch *>(patch));
972 void DealWithCut(const InternalPatch *patchToBeSplit, int axisId, int cutPlace, std::vector<MEDCouplingAutoRefCountObjectPtr<InternalPatch> >& listOfPatches)
974 MEDCouplingAutoRefCountObjectPtr<InternalPatch> leftPart,rightPart;
975 std::vector< std::pair<int,int> > rect(patchToBeSplit->getConstPart());
976 std::vector< std::pair<int,int> > leftRect(rect),rightRect(rect);
977 leftRect[axisId].second=cutPlace+1;
978 rightRect[axisId].first=cutPlace+1;
979 leftPart=patchToBeSplit->extractPart(leftRect);
980 rightPart=patchToBeSplit->extractPart(rightRect);
981 leftPart->zipToFitOnCriterion(); rightPart->zipToFitOnCriterion();
982 listOfPatches.push_back(leftPart);
983 listOfPatches.push_back(rightPart);
989 * 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.
990 * This method only create patches at level 0 relative to \a this.
992 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const std::vector<bool>& criterion, const std::vector<int>& factors)
994 int nbCells(getNumberOfCellsAtCurrentLevel());
995 if(nbCells!=(int)criterion.size())
996 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 !");
998 std::vector<int> cgs(_mesh->getCellGridStructure());
999 std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> > listOfPatches,listOfPatchesOK;
1001 MEDCouplingAutoRefCountObjectPtr<InternalPatch> p(new InternalPatch);
1002 p->setNumberOfTrue(MEDCouplingStructuredMesh::FindMinimalPartOf(cgs,criterion,p->getCriterion(),p->getPart()));
1003 if(p->presenceOfTrue())
1004 listOfPatches.push_back(p);
1005 while(!listOfPatches.empty())
1007 std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> > listOfPatchesTmp;
1008 for(std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> >::iterator it=listOfPatches.begin();it!=listOfPatches.end();it++)
1011 int axisId,rangeOfAxisId,cutPlace;
1013 MEDCouplingStructuredMesh::FindTheWidestAxisOfGivenRangeInCompactFrmt((*it)->getConstPart(),axisId,rangeOfAxisId);
1014 if((*it)->getEfficiency()>=bso.getEffeciency() && (*it)->getNumberOfCells()<bso.getMaxCells())
1015 { listOfPatchesOK.push_back(DealWithNoCut(*it)); continue; }//action 1
1016 FindHole(bso,*it,axisId,cutFound,cutPlace);//axisId overwritten here if FindHole equal to true !
1018 { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 2
1019 FindInflection(bso,*it,cutFound,cutPlace,axisId);//axisId overwritten here if cutFound equal to true !
1021 { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 3
1022 TryAction4(bso,*it,axisId,rangeOfAxisId,cutFound,cutPlace);
1024 { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 4
1025 listOfPatchesOK.push_back(DealWithNoCut(*it));
1027 listOfPatches=listOfPatchesTmp;
1029 for(std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> >::const_iterator it=listOfPatchesOK.begin();it!=listOfPatchesOK.end();it++)
1030 addPatch((*it)->getConstPart(),factors);
1035 * 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.
1036 * This method only create patches at level 0 relative to \a this.
1038 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const DataArrayByte *criterion, const std::vector<int>& factors)
1040 if(!criterion || !criterion->isAllocated())
1041 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createPatchesFromCriterion : the criterion DataArrayByte instance must be allocated and not NULL !");
1042 std::vector<bool> crit(criterion->toVectorOfBool());//check that criterion has one component.
1043 createPatchesFromCriterion(bso,crit,factors);
1047 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const DataArrayDouble *criterion, const std::vector<int>& factors, double eps)
1050 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion : null criterion pointer !");
1051 std::vector<bool> inp(criterion->toVectorOfBool(eps));
1052 createPatchesFromCriterion(bso,inp,factors);
1056 * This method creates a multi level patches split at once.
1057 * This method calls as times as size of \a bso createPatchesFromCriterion. Size of \a bso and size of \a factors must be the same !
1060 * \param [in] criterion
1061 * \param [in] factors
1062 * \param [in] eps - See DataArrayDouble::toVectorOfBool for more information about the semantic of eps.
1064 * \sa createPatchesFromCriterion
1066 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterionML(const std::vector<const INTERP_KERNEL::BoxSplittingOptions *>& bso, const DataArrayDouble *criterion, const std::vector< std::vector<int> >& factors, double eps)
1068 std::size_t nbOfLevs(bso.size());
1069 if(nbOfLevs!=factors.size())
1070 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterionML : size of vectors must be the same !");
1074 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterionML : pointers in 1st arg must be not NULL !");
1075 createPatchesFromCriterion(*bso[0],criterion,factors[0],eps);
1076 for(std::size_t i=1;i<nbOfLevs;i++)
1079 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterionML : presence of a NULL BoxSplittingOptions in input vector !");
1081 std::vector<MEDCouplingCartesianAMRPatchGen *> elts(retrieveGridsAt((int)(nbOfLevs-1)));
1082 std::size_t sz(elts.size());
1083 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGen> > elts2(sz);
1084 std::vector< MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> > elts3(sz);
1085 for(std::size_t ii=0;ii<sz;ii++)
1088 static const char TMP_STR[]="TMP";
1089 std::vector< std::pair<std::string,int> > fieldNames(1); fieldNames[0].first=TMP_STR; fieldNames[0].second=1;
1090 MEDCouplingAutoRefCountObjectPtr<MEDCouplingAMRAttribute> att(MEDCouplingAMRAttribute::New(this,fieldNames,0));
1092 DataArrayDouble *tmpDa(const_cast<DataArrayDouble *>(att->getFieldOn(this,TMP_STR)));
1093 tmpDa->cpyFrom(*criterion);
1094 att->synchronizeCoarseToFine();
1095 for(std::size_t ii=0;ii<sz;ii++)
1097 const DataArrayDouble *critOnLeaf(att->getFieldOn(const_cast<MEDCouplingCartesianAMRMeshGen *>(elts[ii]->getMesh()),TMP_STR));
1098 elts3[ii]=const_cast<DataArrayDouble *>(critOnLeaf); elts3[ii]->incrRef();
1101 for(std::size_t ii=0;ii<sz;ii++)
1102 const_cast<MEDCouplingCartesianAMRMeshGen *>(elts[ii]->getMesh())->createPatchesFromCriterion(*bso[i],elts3[ii],factors[i],eps);
1106 void MEDCouplingCartesianAMRMeshGen::removeAllPatches()
1112 void MEDCouplingCartesianAMRMeshGen::removePatch(int patchId)
1114 checkPatchId(patchId);
1115 int sz((int)_patches.size()),j(0);
1116 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> > patches(sz-1);
1117 for(int i=0;i<sz;i++)
1119 patches[j++]=_patches[i];
1120 (const_cast<MEDCouplingCartesianAMRMeshGen *>(_patches[patchId]->getMesh()))->detachFromFather();
1125 int MEDCouplingCartesianAMRMeshGen::getNumberOfPatches() const
1127 return (int)_patches.size();
1130 int MEDCouplingCartesianAMRMeshGen::getPatchIdFromChildMesh(const MEDCouplingCartesianAMRMeshGen *mesh) const
1133 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ret++)
1135 if((*it)->getMesh()==mesh)
1138 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::getPatchIdFromChildMesh : no such a mesh in my direct progeny !");
1141 std::vector< const MEDCouplingCartesianAMRPatch *> MEDCouplingCartesianAMRMeshGen::getPatches() const
1143 std::size_t sz(_patches.size());
1144 std::vector< const MEDCouplingCartesianAMRPatch *> ret(sz);
1145 for(std::size_t i=0;i<sz;i++)
1150 const MEDCouplingCartesianAMRPatch *MEDCouplingCartesianAMRMeshGen::getPatch(int patchId) const
1152 checkPatchId(patchId);
1153 return _patches[patchId];
1157 * This method states if patch2 (with id \a patchId2) is in the neighborhood of patch1 (with id \a patchId1).
1158 * The neighborhood size is defined by \a ghostLev in the reference of \a this ( \b not in the reference of patches !).
1160 bool MEDCouplingCartesianAMRMeshGen::isPatchInNeighborhoodOf(int patchId1, int patchId2, int ghostLev) const
1162 const MEDCouplingCartesianAMRPatch *p1(getPatch(patchId1)),*p2(getPatch(patchId2));
1163 return p1->isInMyNeighborhood(p2,GetGhostLevelInFineRef(ghostLev,_factors));
1167 * 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.
1168 * 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
1169 * defined by the patch with id \a patchId.
1171 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1172 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1173 * \return DataArrayDouble * - The array of the cell field on the requested patch
1175 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1176 * \throw if \a cellFieldOnThis is NULL or not allocated
1177 * \sa fillCellFieldOnPatch, MEDCouplingIMesh::SpreadCoarseToFine
1179 DataArrayDouble *MEDCouplingCartesianAMRMeshGen::createCellFieldOnPatch(int patchId, const DataArrayDouble *cellFieldOnThis) const
1181 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1182 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatch : the input cell field array is NULL or not allocated !");
1183 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1184 const MEDCouplingIMesh *fine(patch->getMesh()->getImageMesh());
1185 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(fine->getNumberOfCells(),cellFieldOnThis->getNumberOfComponents());
1186 ret->copyStringInfoFrom(*cellFieldOnThis);
1187 MEDCouplingIMesh::SpreadCoarseToFine(cellFieldOnThis,_mesh->getCellGridStructure(),ret,patch->getBLTRRange(),getFactors());
1192 * This method is equivalent to MEDCouplingCartesianAMRMesh::createCellFieldOnPatch except that here instead of creating a new instance
1193 * it fills the value into the \a cellFieldOnPatch data.
1195 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1196 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1197 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1199 * \sa createCellFieldOnPatch, fillCellFieldComingFromPatch
1201 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatch(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch) const
1203 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1204 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatch : the input cell field array is NULL or not allocated !");
1205 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1206 MEDCouplingIMesh::SpreadCoarseToFine(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors());
1210 * This method is the generalization of fillCellFieldOnPatch method. This method only projects coarse to fine without considering the
1211 * potential neighbor patches covered by the ghost cells of patch with id \a patchId.
1213 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1214 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1215 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1216 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1218 * \sa fillCellFieldOnPatch, fillCellFieldOnPatchGhostAdv
1220 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchGhost(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch, int ghostLev) const
1222 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1223 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatchGhost : the input cell field array is NULL or not allocated !");
1224 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1225 MEDCouplingIMesh::SpreadCoarseToFineGhost(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),ghostLev);
1229 * This method is equivalent to fillCellFieldOnPatchGhost except that here \b ONLY \b the \b ghost \b zone will be updated
1230 * in \a cellFieldOnPatch.
1232 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1233 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1234 * \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.
1235 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1237 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyOnGhostZone(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch, int ghostLev) const
1239 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1240 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchOnlyOnGhostZone : the input cell field array is NULL or not allocated !");
1241 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1242 MEDCouplingIMesh::SpreadCoarseToFineGhostZone(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),ghostLev);
1246 * This method is a refinement of fillCellFieldOnPatchGhost. fillCellFieldOnPatchGhost is first called.
1247 * Then for all other patches than those pointed by \a patchId that overlap the ghost zone of the patch impact the ghost zone adequately.
1249 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1250 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1251 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1252 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1253 * \param [in] arrsOnPatches - \b WARNING arrsOnPatches[patchId] is \b NOT \b const. All others are const.
1255 * \sa fillCellFieldOnPatchOnlyGhostAdv
1257 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchGhostAdv(int patchId, const DataArrayDouble *cellFieldOnThis, int ghostLev, const std::vector<const DataArrayDouble *>& arrsOnPatches) const
1259 int nbp(getNumberOfPatches());
1260 if(nbp!=(int)arrsOnPatches.size())
1262 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchGhostAdv : there are " << nbp << " patches in this and " << arrsOnPatches.size() << " arrays in the last parameter !";
1263 throw INTERP_KERNEL::Exception(oss.str().c_str());
1265 DataArrayDouble *theFieldToFill(const_cast<DataArrayDouble *>(arrsOnPatches[patchId]));
1266 // first, do as usual
1267 fillCellFieldOnPatchGhost(patchId,cellFieldOnThis,theFieldToFill,ghostLev);
1268 fillCellFieldOnPatchOnlyGhostAdv(patchId,ghostLev,arrsOnPatches);
1272 * This method updates the patch with id \a patchId considering the only the all the patches in \a this to fill ghost zone.
1273 * So \b warning, the DataArrayDouble instance \a arrsOnPatches[patchId] is non const.
1275 * \sa getPatchIdsInTheNeighborhoodOf
1277 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyGhostAdv(int patchId, int ghostLev, const std::vector<const DataArrayDouble *>& arrsOnPatches) const
1279 int nbp(getNumberOfPatches());
1280 if(nbp!=(int)arrsOnPatches.size())
1282 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchOnlyGhostAdv : there are " << nbp << " patches in this and " << arrsOnPatches.size() << " arrays in the last parameter !";
1283 throw INTERP_KERNEL::Exception(oss.str().c_str());
1285 const MEDCouplingCartesianAMRPatch *refP(getPatch(patchId));
1286 DataArrayDouble *theFieldToFill(const_cast<DataArrayDouble *>(arrsOnPatches[patchId]));
1287 std::vector<int> ids(getPatchIdsInTheNeighborhoodOf(patchId,ghostLev));
1288 for(std::vector<int>::const_iterator it=ids.begin();it!=ids.end();it++)
1290 const MEDCouplingCartesianAMRPatch *otherP(getPatch(*it));
1291 MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(ghostLev,_factors,refP,otherP,theFieldToFill,arrsOnPatches[*it]);
1295 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyOnGhostZoneWith(int ghostLev, const MEDCouplingCartesianAMRPatch *patchToBeModified, const MEDCouplingCartesianAMRPatch *neighborPatch, DataArrayDouble *cellFieldOnPatch, const DataArrayDouble *cellFieldNeighbor) const
1297 MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(ghostLev,_factors,patchToBeModified,neighborPatch,cellFieldOnPatch,cellFieldNeighbor);
1301 * This method updates \a cellFieldOnThis part of values coming from the cell field \a cellFieldOnPatch lying on patch having id \a patchId.
1303 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1304 * \param [in] cellFieldOnPatch - The array of the cell field on patch with id \a patchId.
1305 * \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.
1307 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1308 * \throw if \a cellFieldOnPatch is NULL or not allocated
1309 * \sa createCellFieldOnPatch, MEDCouplingIMesh::CondenseFineToCoarse,fillCellFieldComingFromPatchGhost
1311 void MEDCouplingCartesianAMRMeshGen::fillCellFieldComingFromPatch(int patchId, const DataArrayDouble *cellFieldOnPatch, DataArrayDouble *cellFieldOnThis) const
1313 if(!cellFieldOnPatch || !cellFieldOnPatch->isAllocated())
1314 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch : the input cell field array is NULL or not allocated !");
1315 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1316 MEDCouplingIMesh::CondenseFineToCoarse(_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),cellFieldOnThis);
1320 * This method is the extension of MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch managing the ghost cells. If this
1321 * method is called with \a ghostLev equal to 0 it behaves exactly as MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch.
1323 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1324 * \param [in] cellFieldOnPatch - The array of the cell field on patch with id \a patchId.
1325 * \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.
1326 * \param [in] ghostLev The size of ghost zone (must be >= 0 !)
1328 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1329 * \throw if \a cellFieldOnPatch is NULL or not allocated
1330 * \sa fillCellFieldComingFromPatch
1332 void MEDCouplingCartesianAMRMeshGen::fillCellFieldComingFromPatchGhost(int patchId, const DataArrayDouble *cellFieldOnPatch, DataArrayDouble *cellFieldOnThis, int ghostLev) const
1334 if(!cellFieldOnPatch || !cellFieldOnPatch->isAllocated())
1335 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatchGhost : the input cell field array is NULL or not allocated !");
1336 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1337 MEDCouplingIMesh::CondenseFineToCoarseGhost(_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),cellFieldOnThis,ghostLev);
1341 * This method finds all patches (located by their ids) that are in the neighborhood of patch with id \a patchId. The neighborhood size is
1342 * defined by ghostLev.
1344 * \param [in] patchId - the id of the considered patch.
1345 * \param [in] ghostLev - the size of the neighborhood.
1346 * \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.
1348 DataArrayInt *MEDCouplingCartesianAMRMeshGen::findPatchesInTheNeighborhoodOf(int patchId, int ghostLev) const
1350 int nbp(getNumberOfPatches());
1351 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(0,1);
1352 for(int i=0;i<nbp;i++)
1355 if(isPatchInNeighborhoodOf(i,patchId,ghostLev))
1356 ret->pushBackSilent(i);
1361 MEDCouplingUMesh *MEDCouplingCartesianAMRMeshGen::buildUnstructured() const
1363 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> part(_mesh->buildUnstructured());
1364 std::vector<bool> bs(_mesh->getNumberOfCells(),false);
1365 std::vector<int> cgs(_mesh->getCellGridStructure());
1366 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> > msSafe(_patches.size()+1);
1368 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ii++)
1370 MEDCouplingStructuredMesh::SwitchOnIdsFrom(cgs,(*it)->getBLTRRange(),bs);
1371 msSafe[ii+1]=(*it)->getMesh()->buildUnstructured();
1373 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsOff(DataArrayInt::BuildListOfSwitchedOff(bs));
1374 msSafe[0]=static_cast<MEDCouplingUMesh *>(part->buildPartOfMySelf(eltsOff->begin(),eltsOff->end(),false));
1375 std::vector< const MEDCouplingUMesh * > ms(msSafe.size());
1376 for(std::size_t i=0;i<msSafe.size();i++)
1378 return MEDCouplingUMesh::MergeUMeshes(ms);
1382 * This method returns a mesh containing as cells that there is patches at the current level.
1383 * The patches are seen like 'boxes' that is too say the refinement will not appear here.
1385 * \return MEDCoupling1SGTUMesh * - A new object to be managed by the caller containing as cells as there are patches in \a this.
1387 MEDCoupling1SGTUMesh *MEDCouplingCartesianAMRMeshGen::buildMeshFromPatchEnvelop() const
1389 std::vector<const MEDCoupling1SGTUMesh *> cells;
1390 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > cellsSafe;
1391 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1393 const MEDCouplingCartesianAMRPatch *patch(*it);
1396 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> cell(patch->getMesh()->getImageMesh()->asSingleCell());
1397 MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> cell1SGT(cell->build1SGTUnstructured());
1398 cellsSafe.push_back(cell1SGT); cells.push_back(cell1SGT);
1401 return MEDCoupling1SGTUMesh::Merge1SGTUMeshes(cells);
1404 MEDCoupling1SGTUMesh *MEDCouplingCartesianAMRMeshGen::buildMeshOfDirectChildrenOnly() const
1406 std::vector<const MEDCoupling1SGTUMesh *> patches;
1407 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > patchesSafe;
1408 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1410 const MEDCouplingCartesianAMRPatch *patch(*it);
1413 MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> patchMesh(patch->getMesh()->getImageMesh()->build1SGTUnstructured());
1414 patchesSafe.push_back(patchMesh); patches.push_back(patchMesh);
1417 return MEDCoupling1SGTUMesh::Merge1SGTUMeshes(patches);
1421 * This method works same as buildUnstructured except that arrays are given in input to build a field on cell in output.
1422 * \return MEDCouplingFieldDouble * - a newly created instance the caller has reponsability to deal with.
1423 * \sa buildUnstructured
1425 MEDCouplingFieldDouble *MEDCouplingCartesianAMRMeshGen::buildCellFieldOnRecurseWithoutOverlapWithoutGhost(int ghostSz, const std::vector<const DataArrayDouble *>& recurseArrs) const
1427 if(recurseArrs.empty())
1428 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::buildCellFieldOnRecurseWithoutOverlapWithoutGhost : array is empty ! Should never happen !");
1430 std::vector<bool> bs(_mesh->getNumberOfCells(),false);
1431 std::vector<int> cgs(_mesh->getCellGridStructure());
1432 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> > msSafe(_patches.size()+1);
1434 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ii++)
1436 MEDCouplingStructuredMesh::SwitchOnIdsFrom(cgs,(*it)->getBLTRRange(),bs);
1437 std::vector<const DataArrayDouble *> tmpArrs(extractSubTreeFromGlobalFlatten((*it)->getMesh(),recurseArrs));
1438 msSafe[ii+1]=(*it)->getMesh()->buildCellFieldOnRecurseWithoutOverlapWithoutGhost(ghostSz,tmpArrs);
1440 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsOff(DataArrayInt::BuildListOfSwitchedOff(bs));
1442 MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret(MEDCouplingFieldDouble::New(ON_CELLS));
1443 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr2(extractGhostFrom(ghostSz,recurseArrs[0]));
1444 arr2=arr2->selectByTupleIdSafe(eltsOff->begin(),eltsOff->end());
1445 ret->setArray(arr2);
1446 ret->setName(arr2->getName());
1447 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> part(_mesh->buildUnstructured());
1448 MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> mesh(part->buildPartOfMySelf(eltsOff->begin(),eltsOff->end(),false));
1452 std::vector< const MEDCouplingFieldDouble * > ms(msSafe.size());
1453 for(std::size_t i=0;i<msSafe.size();i++)
1456 return MEDCouplingFieldDouble::MergeFields(ms);
1460 * This method extracts from \arr arr the part inside \a arr by cutting the \a ghostSz external part.
1461 * \arr is expected to be an array having a number of tuples equal to \c getImageMesh()->buildWithGhost(ghostSz).
1463 DataArrayDouble *MEDCouplingCartesianAMRMeshGen::extractGhostFrom(int ghostSz, const DataArrayDouble *arr) const
1465 std::vector<int> st(_mesh->getCellGridStructure());
1466 std::vector< std::pair<int,int> > p(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(st));
1467 std::transform(st.begin(),st.end(),st.begin(),std::bind2nd(std::plus<int>(),2*ghostSz));
1468 MEDCouplingCartesianAMRPatch::ApplyGhostOnCompactFrmt(p,ghostSz);
1469 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(MEDCouplingStructuredMesh::ExtractFieldOfDoubleFrom(st,arr,p));
1474 * This method returns all the patches in \a this not equal to \a patchId that are in neighborhood of patch with id \a patchId.
1476 * \sa fillCellFieldOnPatchOnlyGhostAdv
1478 std::vector<int> MEDCouplingCartesianAMRMeshGen::getPatchIdsInTheNeighborhoodOf(int patchId, int ghostLev) const
1480 std::vector<int> ret;
1481 int nbp(getNumberOfPatches());
1483 for(int i=0;i<nbp;i++)
1486 if(isPatchInNeighborhoodOf(i,patchId,ghostLev))
1492 MEDCouplingCartesianAMRMeshGen::MEDCouplingCartesianAMRMeshGen(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1493 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop):_father(0)
1495 _mesh=MEDCouplingIMesh::New(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop);
1498 MEDCouplingCartesianAMRMeshGen::MEDCouplingCartesianAMRMeshGen(MEDCouplingCartesianAMRMeshGen *father, MEDCouplingIMesh *mesh):_father(father)
1501 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen(MEDCouplingIMesh *mesh) constructor : empty father !");
1503 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen(MEDCouplingIMesh *mesh) constructor : The input mesh is null !");
1504 mesh->checkCoherency();
1505 _mesh=mesh; _mesh->incrRef();
1508 void MEDCouplingCartesianAMRMeshGen::checkPatchId(int patchId) const
1510 int sz(getNumberOfPatches());
1511 if(patchId<0 || patchId>=sz)
1513 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMeshGen::checkPatchId : invalid patchId (" << patchId << ") ! Must be in [0," << sz << ") !";
1514 throw INTERP_KERNEL::Exception(oss.str().c_str());
1518 void MEDCouplingCartesianAMRMeshGen::checkFactorsAndIfNotSetAssign(const std::vector<int>& factors)
1520 if(getSpaceDimension()!=(int)factors.size())
1521 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::checkFactorsAndIfNotSetAssign : invalid size of factors ! size must be equal to the spaceDimension !");
1522 if(_factors.empty())
1528 if(_factors!=factors)
1529 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::checkFactorsAndIfNotSetAssign : the factors ");
1533 void MEDCouplingCartesianAMRMeshGen::retrieveGridsAtInternal(int lev, std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGen> >& grids) const
1537 const MEDCouplingCartesianAMRMesh *thisc(dynamic_cast<const MEDCouplingCartesianAMRMesh *>(this));//tony
1538 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGF> elt(new MEDCouplingCartesianAMRPatchGF(const_cast<MEDCouplingCartesianAMRMesh *>(thisc)));
1539 grids.push_back(DynamicCastSafe<MEDCouplingCartesianAMRPatchGF,MEDCouplingCartesianAMRPatchGen>(elt));
1543 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1545 const MEDCouplingCartesianAMRPatch *pt(*it);
1548 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> tmp1(*it);
1549 grids.push_back(DynamicCastSafe<MEDCouplingCartesianAMRPatch,MEDCouplingCartesianAMRPatchGen>(tmp1));
1555 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1557 const MEDCouplingCartesianAMRPatch *pt(*it);
1559 pt->getMesh()->retrieveGridsAtInternal(lev-1,grids);
1564 int MEDCouplingCartesianAMRMeshGen::GetGhostLevelInFineRef(int ghostLev, const std::vector<int>& factors)
1567 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::GetGhostLevelInFineRef : the ghost size must be >=0 !");
1569 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::GetGhostLevelInFineRef : no factors defined !");
1570 int ghostLevInPatchRef;
1572 ghostLevInPatchRef=0;
1575 ghostLevInPatchRef=(ghostLev-1)/factors[0]+1;
1576 for(std::size_t i=0;i<factors.size();i++)
1577 ghostLevInPatchRef=std::max(ghostLevInPatchRef,(ghostLev-1)/factors[i]+1);
1579 return ghostLevInPatchRef;
1583 * This method returns a sub set of \a all. The subset is defined by the \a head in the tree defined by \a this.
1584 * Elements in \a all are expected to be sorted from god father to most refined structure.
1586 std::vector<const DataArrayDouble *> MEDCouplingCartesianAMRMeshGen::extractSubTreeFromGlobalFlatten(const MEDCouplingCartesianAMRMeshGen *head, const std::vector<const DataArrayDouble *>& all) const
1588 int maxLev(getMaxNumberOfLevelsRelativeToThis());
1589 std::vector<const DataArrayDouble *> ret;
1590 std::vector<const MEDCouplingCartesianAMRMeshGen *> meshes(1,this);
1592 for(int i=0;i<maxLev;i++)
1594 std::vector<const MEDCouplingCartesianAMRMeshGen *> meshesTmp;
1595 for(std::vector<const MEDCouplingCartesianAMRMeshGen *>::const_iterator it=meshes.begin();it!=meshes.end();it++)
1597 if((*it)==head || head->isObjectInTheProgeny(*it))
1598 ret.push_back(all[kk]);
1600 std::vector< const MEDCouplingCartesianAMRPatch *> ps((*it)->getPatches());
1601 for(std::vector< const MEDCouplingCartesianAMRPatch *>::const_iterator it0=ps.begin();it0!=ps.end();it0++)
1603 const MEDCouplingCartesianAMRMeshGen *mesh((*it0)->getMesh());
1604 meshesTmp.push_back(mesh);
1610 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::extractSubTreeFromGlobalFlatten : the size of input vector is not compatible with number of leaves in this !");
1614 std::size_t MEDCouplingCartesianAMRMeshGen::getHeapMemorySizeWithoutChildren() const
1616 return sizeof(MEDCouplingCartesianAMRMeshGen);
1619 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRMeshGen::getDirectChildren() const
1621 std::vector<const BigMemoryObject *> ret;
1622 if((const MEDCouplingIMesh *)_mesh)
1623 ret.push_back((const MEDCouplingIMesh *)_mesh);
1624 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1626 if((const MEDCouplingCartesianAMRPatch*)*it)
1627 ret.push_back((const MEDCouplingCartesianAMRPatch*)*it);
1632 void MEDCouplingCartesianAMRMeshGen::updateTime() const
1634 if((const MEDCouplingIMesh *)_mesh)
1635 updateTimeWith(*_mesh);
1636 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1638 const MEDCouplingCartesianAMRPatch *elt(*it);
1641 const MEDCouplingCartesianAMRMeshGen *mesh(elt->getMesh());
1643 updateTimeWith(*mesh);
1647 MEDCouplingCartesianAMRMeshSub::MEDCouplingCartesianAMRMeshSub(MEDCouplingCartesianAMRMeshGen *father, MEDCouplingIMesh *mesh):MEDCouplingCartesianAMRMeshGen(father,mesh)
1651 MEDCouplingCartesianAMRMesh *MEDCouplingCartesianAMRMesh::New(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1652 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop)
1654 return new MEDCouplingCartesianAMRMesh(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop);
1657 MEDCouplingCartesianAMRMesh::MEDCouplingCartesianAMRMesh(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1658 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop):MEDCouplingCartesianAMRMeshGen(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop)
1662 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRMesh::getDirectChildren() const
1664 std::vector<const BigMemoryObject *> ret(MEDCouplingCartesianAMRMeshGen::getDirectChildren());
1668 MEDCouplingCartesianAMRMesh::~MEDCouplingCartesianAMRMesh()