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 "MEDCouplingFieldDouble.hxx"
23 #include "MEDCoupling1GTUMesh.hxx"
24 #include "MEDCouplingIMesh.hxx"
25 #include "MEDCouplingUMesh.hxx"
31 using namespace ParaMEDMEM;
35 int MEDCouplingCartesianAMRPatchGen::getNumberOfCellsRecursiveWithOverlap() const
37 return _mesh->getNumberOfCellsRecursiveWithOverlap();
40 int MEDCouplingCartesianAMRPatchGen::getNumberOfCellsRecursiveWithoutOverlap() const
42 return _mesh->getNumberOfCellsRecursiveWithoutOverlap();
45 int MEDCouplingCartesianAMRPatchGen::getMaxNumberOfLevelsRelativeToThis() const
47 return _mesh->getMaxNumberOfLevelsRelativeToThis();
50 MEDCouplingCartesianAMRPatchGen::MEDCouplingCartesianAMRPatchGen(MEDCouplingCartesianAMRMeshGen *mesh)
53 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatchGen constructor : input mesh is NULL !");
54 _mesh=mesh; _mesh->incrRef();
57 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatchGen::getMeshSafe() const
59 const MEDCouplingCartesianAMRMeshGen *mesh(_mesh);
61 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatchGen::getMeshSafe const : the mesh is NULL !");
65 MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatchGen::getMeshSafe()
67 MEDCouplingCartesianAMRMeshGen *mesh(_mesh);
69 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatchGen::getMeshSafe : the mesh is NULL !");
73 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRPatchGen::getDirectChildren() const
75 std::vector<const BigMemoryObject *> ret;
76 if((const MEDCouplingCartesianAMRMeshGen *)_mesh)
77 ret.push_back((const MEDCouplingCartesianAMRMeshGen *)_mesh);
82 * \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.
83 * \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,
84 * a the end cell (\b excluded) of the range for the second element of the pair.
86 MEDCouplingCartesianAMRPatch::MEDCouplingCartesianAMRPatch(MEDCouplingCartesianAMRMeshGen *mesh, const std::vector< std::pair<int,int> >& bottomLeftTopRight):MEDCouplingCartesianAMRPatchGen(mesh),_bl_tr(bottomLeftTopRight)
88 int dim((int)bottomLeftTopRight.size()),dimExp(_mesh->getSpaceDimension());
90 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch constructor : space dimension of father and input bottomLeft/topRight size mismatches !");
93 void MEDCouplingCartesianAMRPatch::addPatch(const std::vector< std::pair<int,int> >& bottomLeftTopRight, const std::vector<int>& factors)
95 return getMeshSafe()->addPatch(bottomLeftTopRight,factors);
98 int MEDCouplingCartesianAMRPatch::getNumberOfOverlapedCellsForFather() const
100 return MEDCouplingStructuredMesh::DeduceNumberOfGivenRangeInCompactFrmt(_bl_tr);
104 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
105 * 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 !).
106 * Call isInMyNeighborhoodExt to deal with 2 patches not sharing directly the same father.
108 * \param [in] other - The other patch
109 * \param [in] ghostLev - The size of the neighborhood zone.
111 * \throw if \a this or \a other are invalid (end before start).
112 * \throw if \a ghostLev is \b not >= 0 .
113 * \throw if \a this and \a other have not the same space dimension.
115 * \sa isInMyNeighborhoodExt
117 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhood(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
120 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the size of the neighborhood must be >= 0 !");
122 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the input patch is NULL !");
123 const std::vector< std::pair<int,int> >& thisp(getBLTRRange());
124 const std::vector< std::pair<int,int> >& otherp(other->getBLTRRange());
125 return IsInMyNeighborhood(ghostLev,thisp,otherp);
129 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
130 * 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
131 * ancestor must be the same. If they don't have the same ancestor an exception will be thrown.
133 * \param [in] other - The other patch
134 * \param [in] ghostLev - The size of the neighborhood zone.
136 * \throw if \a this or \a other are invalid (end before start).
137 * \throw if \a ghostLev is \b not >= 0 .
138 * \throw if \a this and \a other have not the same space dimension.
139 * \throw if there is not common ancestor of \a this and \a other.
141 * \sa isInMyNeighborhood, isInMyNeighborhoodDiffLev
143 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
146 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt : the size of the neighborhood must be >= 0 !");
148 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt : the input patch is NULL !");
150 const MEDCouplingCartesianAMRMeshGen *com(FindCommonAncestor(this,other,lev));//check that factors are OK
152 return isInMyNeighborhood(other,ghostLev);
153 std::vector<int> offset(ComputeOffsetFromTwoToOne(com,lev,this,other));
154 const std::vector< std::pair<int,int> >& thisp(getBLTRRange());
155 std::vector< std::pair<int,int> > otherp(other->getBLTRRange());
156 otherp=MEDCouplingStructuredMesh::TranslateCompactFrmt(otherp,offset);
157 return IsInMyNeighborhood(ghostLev,thisp,otherp);
161 * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
162 * the must be >= 0. This method works even if \a this and \a other does not share the same father.
163 * 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.
165 * \param [in] other - The other patch
166 * \param [in] ghostLev - The size of the neighborhood zone.
168 * \throw if \a this or \a other are invalid (end before start).
169 * \throw if \a ghostLev is \b not >= 0 .
170 * \throw if \a this and \a other have not the same space dimension.
171 * \throw if there is not common ancestor of \a this and \a other.
173 * \sa isInMyNeighborhoodExt
175 bool MEDCouplingCartesianAMRPatch::isInMyNeighborhoodDiffLev(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
177 std::vector< std::pair<int,int> > thispp,otherpp;
178 std::vector<int> factors;
179 ComputeZonesOfTwoRelativeToOneDiffLev(ghostLev,this,other,thispp,otherpp,factors);
180 return IsInMyNeighborhood(ghostLev,thispp,otherpp);
183 std::vector<int> MEDCouplingCartesianAMRPatch::computeCellGridSt() const
185 const MEDCouplingCartesianAMRMeshGen *m(getMesh());
187 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::computeCellGridSt : no mesh held by this !");
188 const MEDCouplingCartesianAMRMeshGen *father(m->getFather());
190 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::computeCellGridSt : no father help by underlying mesh !");
191 const std::vector< std::pair<int,int> >& bltr(getBLTRRange());
192 const std::vector<int>& factors(father->getFactors());
193 std::vector<int> ret(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(bltr));
194 std::transform(ret.begin(),ret.end(),factors.begin(),ret.begin(),std::multiplies<int>());
198 bool MEDCouplingCartesianAMRPatch::IsInMyNeighborhood(int ghostLev, const std::vector< std::pair<int,int> >& p1, const std::vector< std::pair<int,int> >& p2)
200 std::size_t thispsize(p1.size());
201 if(thispsize!=p2.size())
202 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the dimensions must be the same !");
203 for(std::size_t i=0;i<thispsize;i++)
205 const std::pair<int,int>& thispp(p1[i]);
206 const std::pair<int,int>& otherpp(p2[i]);
207 if(thispp.second<thispp.first)
208 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
209 if(otherpp.second<otherpp.first)
210 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
211 if(otherpp.first==thispp.second+ghostLev-1)
213 if(otherpp.second+ghostLev-1==thispp.first)
215 int start(std::max(thispp.first,otherpp.first)),end(std::min(thispp.second,otherpp.second));
223 * \sa FindNeighborsOfSubPatchesOf
225 std::vector< std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > > MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOfSameLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2)
228 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOfSameLev : the input pointers must be not NULL !");
229 std::vector< std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > > ret;
230 std::vector< const MEDCouplingCartesianAMRPatch *> p1Work(p1->getMesh()->getPatches()),p2Work(p2->getMesh()->getPatches());
231 while(!p1Work.empty())
233 std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > retTmp;
234 std::vector<const MEDCouplingCartesianAMRPatch *> p1Work2,p2Work2;
235 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it1=p1Work.begin();it1!=p1Work.end();it1++)
237 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it2=p2Work.begin();it2!=p2Work.end();it2++)
239 if((*it1)->isInMyNeighborhoodExt(*it2,ghostLev))
240 retTmp.push_back(std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *>(*it1,*it2));
242 std::vector<const MEDCouplingCartesianAMRPatch *> tmp1((*it1)->getMesh()->getPatches());
243 p1Work2.insert(p1Work2.end(),tmp1.begin(),tmp1.end());
245 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it2=p2Work.begin();it2!=p2Work.end();it2++)
247 std::vector<const MEDCouplingCartesianAMRPatch *> tmp2((*it2)->getMesh()->getPatches());
248 p2Work2.insert(p2Work2.end(),tmp2.begin(),tmp2.end());
250 ret.push_back(retTmp);
258 * This method returns all pair of patches (pa,pb) so that pb is in the neighborhood of pa (size of neighborhood is \a ghostLev).
259 * 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
260 * FindNeighborsOfSubPatchesOfSameLev.
262 * \sa FindNeighborsOfSubPatchesOfSameLev
264 void MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOf(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, std::vector< std::pair<const MEDCouplingCartesianAMRPatch *, const MEDCouplingCartesianAMRPatch *> >& ret)
267 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOf : the input pointers must be not NULL !");
268 std::vector< const MEDCouplingCartesianAMRPatch *> p1Work(p1->getMesh()->getPatches());
269 while(!p1Work.empty())
271 std::vector<const MEDCouplingCartesianAMRPatch *> p1Work2;
272 for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it0=p1Work.begin();it0!=p1Work.end();it0++)
274 if((*it0)->isInMyNeighborhoodDiffLev(p2,ghostLev))
275 ret.push_back(std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *>(*it0,p2));
276 std::vector<const MEDCouplingCartesianAMRPatch *> tmp2((*it0)->getMesh()->getPatches());
277 p1Work2.insert(p1Work2.end(),tmp2.begin(),tmp2.end());
284 * \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.
286 * \saUpdateNeighborsOfOneWithTwoExt
288 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(int ghostLev, const std::vector<int>& factors, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
290 const std::vector< std::pair<int,int> >& p1BLTR(p1->getBLTRRange());
291 const std::vector< std::pair<int,int> >& p2BLTR(p2->getBLTRRange());
292 UpdateNeighborsOfOneWithTwoInternal(ghostLev,factors,p1BLTR,p2BLTR,dataOnP1,dataOnP2);
296 * Idem than UpdateNeighborsOfOneWithTwo, except that here \a p1 and \a p2 are not sharing the same direct father.
298 * \sa UpdateNeighborsOfOneWithTwo
300 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoExt(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
302 const std::vector< std::pair<int,int> >& p1BLTR(p1->getBLTRRange());//p1BLTR=[(10,12),(5,8)]
303 std::vector< std::pair<int,int> > p2BLTR(p2->getBLTRRange());//p2BLTR=[(0,1),(0,5)]
305 const MEDCouplingCartesianAMRMeshGen *ca(FindCommonAncestor(p1,p2,lev));
306 std::vector<int> offset(ComputeOffsetFromTwoToOne(ca,lev,p1,p2));//[12,4]
307 p2BLTR=MEDCouplingStructuredMesh::TranslateCompactFrmt(p2BLTR,offset);//p2BLTR=[(12,13),(4,9)]
308 UpdateNeighborsOfOneWithTwoInternal(ghostLev,p1->getMesh()->getFather()->getFactors(),p1BLTR,p2BLTR,dataOnP1,dataOnP2);
312 * \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 !
314 void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoMixedLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
316 std::vector< std::pair<int,int> > p1pp,p2pp;
317 std::vector<int> factors;
318 ComputeZonesOfTwoRelativeToOneDiffLev(ghostLev,p1,p2,p1pp,p2pp,factors);
320 std::vector<int> dimsP2NotRefined(p2->computeCellGridSt());
321 std::vector<int> dimsP2Refined(dimsP2NotRefined);
322 std::transform(dimsP2NotRefined.begin(),dimsP2NotRefined.end(),factors.begin(),dimsP2Refined.begin(),std::multiplies<int>());
323 std::vector< std::pair<int,int> > p2RefinedAbs(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(dimsP2NotRefined));
324 std::vector<int> dimsP2RefinedGhost(dimsP2Refined.size());
325 std::transform(dimsP2Refined.begin(),dimsP2Refined.end(),dimsP2RefinedGhost.begin(),std::bind2nd(std::plus<int>(),2*ghostLev));
326 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> fineP2(DataArrayDouble::New()); fineP2->alloc(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(dimsP2RefinedGhost),dataOnP2->getNumberOfComponents());
327 MEDCouplingIMesh::SpreadCoarseToFineGhost(dataOnP2,dimsP2NotRefined,fineP2,p2RefinedAbs,factors,ghostLev);
329 UpdateNeighborsOfOneWithTwoInternal(ghostLev,p1->getMesh()->getFather()->getFactors(),p1pp,p2pp,dataOnP1,fineP2);
333 * \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 !
334 * 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
335 * on the same level as \a p1.
337 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)
339 std::vector<const MEDCouplingCartesianAMRMeshGen *> ancestorsOfThis;
340 const MEDCouplingCartesianAMRMeshGen *work(p1->getMesh()),*work2(0);
341 ancestorsOfThis.push_back(work);
344 work=work->getFather();
346 ancestorsOfThis.push_back(work);
351 std::size_t levThis(0),levOther(0);
352 while(work && !found)
355 work=work->getFather();
359 std::vector<const MEDCouplingCartesianAMRMeshGen *>::iterator it(std::find(ancestorsOfThis.begin(),ancestorsOfThis.end(),work));
360 if(it!=ancestorsOfThis.end())
362 levThis=std::distance(ancestorsOfThis.begin(),it);
368 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev : no common ancestor found !");
369 if(levThis<=levOther)
370 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev : this method is not called correctly !");
372 const MEDCouplingCartesianAMRMeshGen *comAncestor(ancestorsOfThis[levThis]);
373 int idThis(comAncestor->getPatchIdFromChildMesh(ancestorsOfThis[levThis-1])),idOther(comAncestor->getPatchIdFromChildMesh(work2));
374 const MEDCouplingCartesianAMRPatch *thisp(comAncestor->getPatch(idThis)),*otherp(comAncestor->getPatch(idOther));
375 std::vector<int> offset(ComputeOffsetFromTwoToOne(comAncestor,levOther,thisp,otherp));
376 p1Zone=thisp->getBLTRRange(); p2Zone=MEDCouplingStructuredMesh::TranslateCompactFrmt(otherp->getBLTRRange(),offset);
377 factToApplyOn2.resize(p1Zone.size()); std::fill(factToApplyOn2.begin(),factToApplyOn2.end(),1);
379 std::size_t nbOfTurn(levThis-levOther);
380 for(std::size_t i=0;i<nbOfTurn;i++)
382 std::vector< std::pair<int,int> > tmp0;
383 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p1Zone,p2Zone,tmp0,false);
385 const MEDCouplingCartesianAMRMeshGen *curAncestor(ancestorsOfThis[levThis-i]);
386 ApplyFactorsOnCompactFrmt(p2Zone,curAncestor->getFactors());
387 curAncestor=ancestorsOfThis[levThis-1-i];
388 const std::vector<int>& factors(curAncestor->getFactors());
389 std::transform(factToApplyOn2.begin(),factToApplyOn2.end(),factors.begin(),factToApplyOn2.begin(),std::multiplies<int>());
390 int tmpId(curAncestor->getPatchIdFromChildMesh(ancestorsOfThis[levThis-2-i]));
391 p1Zone=curAncestor->getPatch(tmpId)->getBLTRRange();
395 std::size_t MEDCouplingCartesianAMRPatch::getHeapMemorySizeWithoutChildren() const
397 std::size_t ret(sizeof(MEDCouplingCartesianAMRPatch));
398 ret+=_bl_tr.capacity()*sizeof(std::pair<int,int>);
402 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRPatch::FindCommonAncestor(const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, int& lev)
404 const MEDCouplingCartesianAMRMeshGen *f1(p1->_mesh),*f2(p2->_mesh);
406 while(f1!=f2 || f1==0 || f2==0)
408 f1=f1->getFather(); f2=f2->getFather();
409 if(f1->getFactors()!=f2->getFactors())
410 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindCommonAncestor : factors differ !");
414 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindCommonAncestor : no common ancestor between p1 and p2 !");
418 std::vector<int> MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne(const MEDCouplingCartesianAMRMeshGen *comAncestor, int lev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2)
421 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne : this method is useful only for lev > 0 !");
423 int dim(p1->getMesh()->getSpaceDimension());
424 if(p2->getMesh()->getSpaceDimension()!=dim)
425 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeOffsetFromTwoToOne : dimension must be the same !");
426 std::vector< int > ret(dim,0);
427 for(int i=0;i<zeLev;i++)
429 const MEDCouplingCartesianAMRMeshGen *f1(p1->_mesh),*f2(p2->_mesh);
430 const MEDCouplingCartesianAMRPatch *p1h(0),*p2h(0);
431 for(int j=0;j<lev-i;j++)
433 const MEDCouplingCartesianAMRMeshGen *f1tmp(f1->getFather()),*f2tmp(f2->getFather());
434 int pid1(f1tmp->getPatchIdFromChildMesh(f1)),pid2(f2tmp->getPatchIdFromChildMesh(f2));
435 p1h=f1tmp->getPatch(pid1); p2h=f2tmp->getPatch(pid2);
438 std::vector< std::pair<int,int> > p2c(p2h->getBLTRRange());
439 for(int k=0;k<dim;k++)
441 p2c[k].first+=ret[k];
442 p2c[k].second+=ret[k];
444 for(int k=0;k<dim;k++)
446 ret[k]=p2c[k].first-p1h->getBLTRRange()[k].first;
447 ret[k]*=f1->getFactors()[k];
453 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)
454 {//p1=[(1,4),(2,4)] p2=[(4,5),(3,4)]
455 int dim((int)factors.size());
456 std::vector<int> dimsCoarse(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(p1));//[3,2]
457 std::transform(dimsCoarse.begin(),dimsCoarse.end(),factors.begin(),dimsCoarse.begin(),std::multiplies<int>());//[12,8]
458 std::transform(dimsCoarse.begin(),dimsCoarse.end(),dimsCoarse.begin(),std::bind2nd(std::plus<int>(),2*ghostLev));//[14,10]
459 std::vector< std::pair<int,int> > rangeCoarse(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(dimsCoarse));//[(0,14),(0,10)]
460 std::vector<int> fakeFactors(dim,1);
462 std::vector< std::pair<int,int> > tmp0,tmp1,tmp2;
463 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p1,p2,tmp0,false);//tmp0=[(3,4),(1,2)]
464 ApplyFactorsOnCompactFrmt(tmp0,factors);//tmp0=[(12,16),(4,8)]
465 ApplyGhostOnCompactFrmt(tmp0,ghostLev);//tmp0=[(13,17),(5,9)]
466 std::vector< std::pair<int,int> > interstRange(MEDCouplingStructuredMesh::IntersectRanges(tmp0,rangeCoarse));//interstRange=[(13,14),(5,9)]
467 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p2,p1,tmp1,false);//tmp1=[(-3,0),(-1,1)]
468 ApplyFactorsOnCompactFrmt(tmp1,factors);//tmp1=[(-12,-4),(-4,0)]
469 MEDCouplingStructuredMesh::ChangeReferenceToGlobalOfCompactFrmt(tmp1,interstRange,tmp2,false);//tmp2=[(1,2),(1,5)]
471 std::vector< std::pair<int,int> > dimsFine(p2);
472 ApplyFactorsOnCompactFrmt(dimsFine,factors);
473 ApplyAllGhostOnCompactFrmt(dimsFine,ghostLev);
475 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ghostVals(MEDCouplingStructuredMesh::ExtractFieldOfDoubleFrom(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(dimsFine),dataOnP2,tmp2));
476 MEDCouplingIMesh::CondenseFineToCoarse(dimsCoarse,ghostVals,interstRange,fakeFactors,dataOnP1);
480 * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in refined reference.
481 * \param [in] factors - the factors per axis.
483 void MEDCouplingCartesianAMRPatch::ApplyFactorsOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, const std::vector<int>& factors)
485 std::size_t sz(factors.size());
486 if(sz!=partBeforeFact.size())
487 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ApplyFactorsOnCompactFrmt : size of input vectors must be the same !");
488 for(std::size_t i=0;i<sz;i++)
490 partBeforeFact[i].first*=factors[i];
491 partBeforeFact[i].second*=factors[i];
496 * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in ghost reference.
497 * \param [in] ghostSize - the ghost size of zone for all axis.
499 void MEDCouplingCartesianAMRPatch::ApplyGhostOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, int ghostSize)
502 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ApplyGhostOnCompactFrmt : ghost size must be >= 0 !");
503 std::size_t sz(partBeforeFact.size());
504 for(std::size_t i=0;i<sz;i++)
506 partBeforeFact[i].first+=ghostSize;
507 partBeforeFact[i].second+=ghostSize;
512 * This method is different than ApplyGhostOnCompactFrmt. The \a partBeforeFact parameter is enlarger contrary to ApplyGhostOnCompactFrmt.
514 * \param [in,out] partBeforeFact - the part of a image mesh in compact format that will be put in ghost reference.
515 * \param [in] ghostSize - the ghost size of zone for all axis.
517 void MEDCouplingCartesianAMRPatch::ApplyAllGhostOnCompactFrmt(std::vector< std::pair<int,int> >& partBeforeFact, int ghostSize)
520 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ApplyAllGhostOnCompactFrmt : ghost size must be >= 0 !");
521 std::size_t sz(partBeforeFact.size());
522 for(std::size_t i=0;i<sz;i++)
524 partBeforeFact[i].first-=ghostSize;
525 partBeforeFact[i].second+=ghostSize;
529 MEDCouplingCartesianAMRPatchGF::MEDCouplingCartesianAMRPatchGF(MEDCouplingCartesianAMRMesh *mesh):MEDCouplingCartesianAMRPatchGen(mesh)
533 std::size_t MEDCouplingCartesianAMRPatchGF::getHeapMemorySizeWithoutChildren() const
535 return sizeof(MEDCouplingCartesianAMRPatchGF);
538 MEDCouplingDataForGodFather::MEDCouplingDataForGodFather(MEDCouplingCartesianAMRMesh *gf):_gf(gf),_tlc(gf)
541 throw INTERP_KERNEL::Exception("MEDCouplingDataForGodFather constructor : A data has to be attached to a AMR Mesh instance !");
545 void MEDCouplingDataForGodFather::checkGodFatherFrozen() const
550 bool MEDCouplingDataForGodFather::changeGodFather(MEDCouplingCartesianAMRMesh *gf)
552 bool ret(_tlc.keepTrackOfNewTL(gf));
562 int MEDCouplingCartesianAMRMeshGen::getSpaceDimension() const
564 return _mesh->getSpaceDimension();
567 void MEDCouplingCartesianAMRMeshGen::setFactors(const std::vector<int>& newFactors)
569 if(getSpaceDimension()!=(int)newFactors.size())
570 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::setFactors : size of input factors is not equal to the space dimension !");
576 if(_factors==newFactors)
578 if(!_patches.empty())
579 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::setFactors : modification of factors is not allowed when presence of patches !");
584 int MEDCouplingCartesianAMRMeshGen::getMaxNumberOfLevelsRelativeToThis() const
587 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
588 ret=std::max(ret,(*it)->getMaxNumberOfLevelsRelativeToThis()+1);
593 * This method returns the number of cells of \a this with the help of the MEDCouplingIMesh instance representing \a this.
594 * The patches in \a this are ignored here.
595 * \sa getNumberOfCellsAtCurrentLevelGhost, getNumberOfCellsRecursiveWithOverlap
597 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsAtCurrentLevel() const
599 return _mesh->getNumberOfCells();
603 * 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
604 * to take into account of the ghost cells for future computation.
605 * The patches in \a this are ignored here.
607 * \sa getNumberOfCellsAtCurrentLevel
609 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsAtCurrentLevelGhost(int ghostLev) const
611 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> tmp(_mesh->buildWithGhost(ghostLev));
612 return tmp->getNumberOfCells();
616 * This method returns the number of cells including the current level but \b also \b including recursively all cells of other levels
617 * starting from this. The set of cells which size is returned here are generally overlapping each other.
619 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsRecursiveWithOverlap() const
621 int ret(_mesh->getNumberOfCells());
622 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
624 ret+=(*it)->getNumberOfCellsRecursiveWithOverlap();
630 * This method returns the max number of cells covering all the space without overlapping.
631 * It returns the number of cells of the mesh with the highest resolution.
632 * The returned value is equal to the number of cells of mesh returned by buildUnstructured.
634 * \sa buildUnstructured
636 int MEDCouplingCartesianAMRMeshGen::getNumberOfCellsRecursiveWithoutOverlap() const
638 int ret(_mesh->getNumberOfCells());
639 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
641 ret-=(*it)->getNumberOfOverlapedCellsForFather();
642 ret+=(*it)->getNumberOfCellsRecursiveWithoutOverlap();
647 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMeshGen::getFather() const
652 const MEDCouplingCartesianAMRMeshGen *MEDCouplingCartesianAMRMeshGen::getGodFather() const
657 return _father->getGodFather();
661 * This method returns the level of \a this. 0 for god father. -1 for children of god father ...
663 int MEDCouplingCartesianAMRMeshGen::getAbsoluteLevel() const
668 return _father->getAbsoluteLevel()-1;
672 * This method returns grids relative to god father to specified level \a absoluteLev.
674 * \return std::vector<MEDCouplingCartesianAMRPatchGen *> - objects in vector are to be managed (decrRef) by the caller.
676 std::vector<MEDCouplingCartesianAMRPatchGen *> MEDCouplingCartesianAMRMeshGen::retrieveGridsAt(int absoluteLev) const
679 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::retrieveGridsAt : absolute level must be >=0 !");
681 return getGodFather()->retrieveGridsAt(absoluteLev);
683 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGen> > rets;
684 retrieveGridsAtInternal(absoluteLev,rets);
685 std::vector< MEDCouplingCartesianAMRPatchGen * > ret(rets.size());
686 for(std::size_t i=0;i<rets.size();i++)
688 ret[i]=rets[i].retn();
693 void MEDCouplingCartesianAMRMeshGen::detachFromFather()
700 * \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,
701 * a the end cell (\b excluded) of the range for the second element of the pair.
702 * \param [in] factors The factor of refinement per axis (different from 0).
704 void MEDCouplingCartesianAMRMeshGen::addPatch(const std::vector< std::pair<int,int> >& bottomLeftTopRight, const std::vector<int>& factors)
706 checkFactorsAndIfNotSetAssign(factors);
707 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> mesh(static_cast<MEDCouplingIMesh *>(_mesh->buildStructuredSubPart(bottomLeftTopRight)));
708 mesh->refineWithFactor(factors);
709 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRMeshSub> zeMesh(new MEDCouplingCartesianAMRMeshSub(this,mesh));
710 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> elt(new MEDCouplingCartesianAMRPatch(zeMesh,bottomLeftTopRight));
711 _patches.push_back(elt);
717 class InternalPatch : public RefCountObjectOnly
720 InternalPatch():_nb_of_true(0) { }
721 int getDimension() const { return (int)_part.size(); }
722 double getEfficiency() const { return (double)_nb_of_true/(double)_crit.size(); }
723 int getNumberOfCells() const { return (int)_crit.size(); }
724 void setNumberOfTrue(int nboft) { _nb_of_true=nboft; }
725 std::vector<bool>& getCriterion() { return _crit; }
726 const std::vector<bool>& getConstCriterion() const { return _crit; }
727 void setPart(const std::vector< std::pair<int,int> >& part) { _part=part; }
728 std::vector< std::pair<int,int> >& getPart() { return _part; }
729 const std::vector< std::pair<int,int> >& getConstPart() const { return _part; }
730 bool presenceOfTrue() const { return _nb_of_true>0; }
731 std::vector<int> computeCGS() const { return MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(_part); }
732 std::vector< std::vector<int> > computeSignature() const { return MEDCouplingStructuredMesh::ComputeSignaturePerAxisOf(computeCGS(),getConstCriterion()); }
733 double getEfficiencyPerAxis(int axisId) const { return (double)_nb_of_true/((double)(_part[axisId].second-_part[axisId].first)); }
734 void zipToFitOnCriterion();
735 void updateNumberOfTrue() const;
736 MEDCouplingAutoRefCountObjectPtr<InternalPatch> extractPart(const std::vector< std::pair<int,int> >&partInGlobal) const;
737 MEDCouplingAutoRefCountObjectPtr<InternalPatch> deepCpy() const;
741 mutable int _nb_of_true;
742 std::vector<bool> _crit;
744 std::vector< std::pair<int,int> > _part;
747 void InternalPatch::zipToFitOnCriterion()
749 std::vector<int> cgs(computeCGS());
750 std::vector<bool> newCrit;
751 std::vector< std::pair<int,int> > newPart,newPart2;
752 int newNbOfTrue(MEDCouplingStructuredMesh::FindMinimalPartOf(cgs,_crit,newCrit,newPart));
753 MEDCouplingStructuredMesh::ChangeReferenceToGlobalOfCompactFrmt(_part,newPart,newPart2);
754 if(newNbOfTrue!=_nb_of_true)
755 throw INTERP_KERNEL::Exception("InternalPatch::zipToFitOnCrit : internal error !");
756 _crit=newCrit; _part=newPart2;
759 void InternalPatch::updateNumberOfTrue() const
761 _nb_of_true=(int)std::count(_crit.begin(),_crit.end(),true);
764 MEDCouplingAutoRefCountObjectPtr<InternalPatch> InternalPatch::extractPart(const std::vector< std::pair<int,int> >&partInGlobal) const
766 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(new InternalPatch);
767 std::vector<int> cgs(computeCGS());
768 std::vector< std::pair<int,int> > newPart;
769 MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(_part,partInGlobal,newPart);
770 MEDCouplingStructuredMesh::ExtractFieldOfBoolFrom(cgs,_crit,newPart,ret->getCriterion());
771 ret->setPart(partInGlobal);
772 ret->updateNumberOfTrue();
776 MEDCouplingAutoRefCountObjectPtr<InternalPatch> InternalPatch::deepCpy() const
778 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(new InternalPatch);
783 void DissectBigPatch(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int rangeOfAxisId, bool& cutFound, int& cutPlace)
785 cutFound=false; cutPlace=-1;
786 std::vector<double> ratio(rangeOfAxisId-1);
787 for(int id=0;id<rangeOfAxisId-1;id++)
789 double efficiency[2];
792 std::vector< std::pair<int,int> > rectH(patchToBeSplit->getConstPart());
794 rectH[axisId].second=patchToBeSplit->getConstPart()[axisId].first+id;
796 rectH[axisId].first=patchToBeSplit->getConstPart()[axisId].first+id;
797 MEDCouplingAutoRefCountObjectPtr<InternalPatch> p(patchToBeSplit->deepCpy());
798 p->zipToFitOnCriterion();
800 efficiency[h]=p->getEfficiencyPerAxis(axisId);
802 ratio[id]=std::max(efficiency[0],efficiency[1])/std::min(efficiency[0],efficiency[1]);
804 int minCellDirection(bso.getMinCellDirection()),indexMin(-1);
805 int dimRatioInner(rangeOfAxisId-1-2*(minCellDirection-1));
806 std::vector<double> ratio_inner(dimRatioInner);
807 double minRatio(1.e10);
808 for(int i=0; i<dimRatioInner; i++)
810 if(ratio[minCellDirection-1+i]<minRatio)
812 minRatio=ratio[minCellDirection-1+i];
813 indexMin=i+minCellDirection;
816 cutFound=true; cutPlace=indexMin+patchToBeSplit->getConstPart()[axisId].first-1;
819 void FindHole(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int& axisId, bool& cutFound, int& cutPlace)
821 cutPlace=-1; cutFound=false;
822 int minCellDirection(bso.getMinCellDirection());
823 const int dim(patchToBeSplit->getDimension());
824 std::vector< std::vector<int> > signatures(patchToBeSplit->computeSignature());
825 for(int id=0;id<dim;id++)
827 const std::vector<int>& signature(signatures[id]);
828 std::vector<int> hole;
829 std::vector<double> distance;
830 int len((int)signature.size());
831 for(int i=0;i<len;i++)
833 if(len>= 2*minCellDirection && i >= minCellDirection-1 && i <= len-minCellDirection-1)
837 double center(((double)len/2.));
838 for(std::size_t i=0;i<hole.size();i++)
839 distance.push_back(fabs(hole[i]+1.-center));
841 std::size_t posDistanceMin(std::distance(distance.begin(),std::min_element(distance.begin(),distance.end())));
844 cutPlace=hole[posDistanceMin]+patchToBeSplit->getConstPart()[axisId].first+1;
850 void FindInflection(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, bool& cutFound, int& cutPlace, int& axisId)
852 cutFound=false; cutPlace=-1;// do not set axisId before to be sure that cutFound was set to true
854 const std::vector< std::pair<int,int> >& part(patchToBeSplit->getConstPart());
855 int sign,minCellDirection(bso.getMinCellDirection());
856 const int dim(patchToBeSplit->getDimension());
858 std::vector<int> zeroCrossDims(dim,-1);
859 std::vector<int> zeroCrossVals(dim,-1);
860 std::vector< std::vector<int> > signatures(patchToBeSplit->computeSignature());
861 for (int id=0;id<dim;id++)
863 const std::vector<int>& signature(signatures[id]);
865 std::vector<int> derivate_second_order,gradient_absolute,signe_change,zero_cross,edge,max_cross_list ;
866 std::vector<double> distance ;
868 for (unsigned int i=1;i<signature.size()-1;i++)
869 derivate_second_order.push_back(signature[i-1]-2*signature[i]+signature[i+1]) ;
871 // Gradient absolute value
872 for ( unsigned int i=1;i<derivate_second_order.size();i++)
873 gradient_absolute.push_back(fabs(derivate_second_order[i]-derivate_second_order[i-1])) ;
874 if(derivate_second_order.empty())
876 for (unsigned int i=0;i<derivate_second_order.size()-1;i++)
878 if (derivate_second_order[i]*derivate_second_order[i+1] < 0 )
880 if (derivate_second_order[i]*derivate_second_order[i+1] > 0 )
882 if (derivate_second_order[i]*derivate_second_order[i+1] == 0 )
884 if ( sign==0 || sign==-1 )
885 if ( i >= (unsigned int)minCellDirection-2 && i <= signature.size()-minCellDirection-2 )
887 zero_cross.push_back(i) ;
888 edge.push_back(gradient_absolute[i]) ;
890 signe_change.push_back(sign) ;
892 if ( zero_cross.size() > 0 )
894 int max_cross=*max_element(edge.begin(),edge.end()) ;
895 for (unsigned int i=0;i<edge.size();i++)
896 if (edge[i]==max_cross)
897 max_cross_list.push_back(zero_cross[i]+1) ;
899 double center((signature.size()/2.0));
900 for (unsigned int i=0;i<max_cross_list.size();i++)
901 distance.push_back(fabs(max_cross_list[i]+1-center));
903 float distance_min=*min_element(distance.begin(),distance.end()) ;
904 int pos_distance_min=find(distance.begin(),distance.end(),distance_min)-distance.begin();
905 int best_place = max_cross_list[pos_distance_min] + part[id].first ;
908 zeroCrossDims[id] = best_place ;
909 zeroCrossVals[id] = max_cross ;
912 derivate_second_order.clear() ;
913 gradient_absolute.clear() ;
914 signe_change.clear() ;
917 max_cross_list.clear() ;
921 if ( zeroCrossDims[0]!=-1 || zeroCrossDims[1]!=-1 )
923 int max_cross_dims = *max_element(zeroCrossVals.begin(),zeroCrossVals.end()) ;
925 if (zeroCrossVals[0]==max_cross_dims && zeroCrossVals[1]==max_cross_dims )
927 int nl_left(part[0].second-part[0].first);
928 int nc_left(part[1].second-part[1].first);
929 if ( nl_left >= nc_left )
935 max_cross_dims=std::find(zeroCrossVals.begin(),zeroCrossVals.end(),max_cross_dims)-zeroCrossVals.begin();
937 cutPlace=zeroCrossDims[max_cross_dims];
938 axisId=max_cross_dims ;
942 void TryAction4(const INTERP_KERNEL::BoxSplittingOptions& bso, const InternalPatch *patchToBeSplit, int axisId, int rangeOfAxisId, bool& cutFound, int& cutPlace)
945 if(patchToBeSplit->getEfficiency()<=bso.getEffeciencySnd())
947 if(rangeOfAxisId>=2*bso.getMinCellDirection())
950 cutPlace=rangeOfAxisId/2+patchToBeSplit->getConstPart()[axisId].first-1;
955 if(patchToBeSplit->getNumberOfCells()>bso.getMaxCells())
957 DissectBigPatch(bso,patchToBeSplit,axisId,rangeOfAxisId,cutFound,cutPlace);
962 MEDCouplingAutoRefCountObjectPtr<InternalPatch> DealWithNoCut(const InternalPatch *patch)
964 MEDCouplingAutoRefCountObjectPtr<InternalPatch> ret(const_cast<InternalPatch *>(patch));
969 void DealWithCut(const InternalPatch *patchToBeSplit, int axisId, int cutPlace, std::vector<MEDCouplingAutoRefCountObjectPtr<InternalPatch> >& listOfPatches)
971 MEDCouplingAutoRefCountObjectPtr<InternalPatch> leftPart,rightPart;
972 std::vector< std::pair<int,int> > rect(patchToBeSplit->getConstPart());
973 std::vector< std::pair<int,int> > leftRect(rect),rightRect(rect);
974 leftRect[axisId].second=cutPlace+1;
975 rightRect[axisId].first=cutPlace+1;
976 leftPart=patchToBeSplit->extractPart(leftRect);
977 rightPart=patchToBeSplit->extractPart(rightRect);
978 leftPart->zipToFitOnCriterion(); rightPart->zipToFitOnCriterion();
979 listOfPatches.push_back(leftPart);
980 listOfPatches.push_back(rightPart);
986 * 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.
988 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const std::vector<bool>& criterion, const std::vector<int>& factors)
990 int nbCells(getNumberOfCellsAtCurrentLevel());
991 if(nbCells!=(int)criterion.size())
992 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 !");
994 std::vector<int> cgs(_mesh->getCellGridStructure());
995 std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> > listOfPatches,listOfPatchesOK;
997 MEDCouplingAutoRefCountObjectPtr<InternalPatch> p(new InternalPatch);
998 p->setNumberOfTrue(MEDCouplingStructuredMesh::FindMinimalPartOf(cgs,criterion,p->getCriterion(),p->getPart()));
999 if(p->presenceOfTrue())
1000 listOfPatches.push_back(p);
1001 while(!listOfPatches.empty())
1003 std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> > listOfPatchesTmp;
1004 for(std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> >::iterator it=listOfPatches.begin();it!=listOfPatches.end();it++)
1007 int axisId,rangeOfAxisId,cutPlace;
1009 MEDCouplingStructuredMesh::FindTheWidestAxisOfGivenRangeInCompactFrmt((*it)->getConstPart(),axisId,rangeOfAxisId);
1010 if((*it)->getEfficiency()>=bso.getEffeciency() && (*it)->getNumberOfCells()<bso.getMaxCells())
1011 { listOfPatchesOK.push_back(DealWithNoCut(*it)); continue; }//action 1
1012 FindHole(bso,*it,axisId,cutFound,cutPlace);//axisId overwritten here if FindHole equal to true !
1014 { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 2
1015 FindInflection(bso,*it,cutFound,cutPlace,axisId);//axisId overwritten here if cutFound equal to true !
1017 { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 3
1018 TryAction4(bso,*it,axisId,rangeOfAxisId,cutFound,cutPlace);
1020 { DealWithCut(*it,axisId,cutPlace,listOfPatchesTmp); continue; }//action 4
1021 listOfPatchesOK.push_back(DealWithNoCut(*it));
1023 listOfPatches=listOfPatchesTmp;
1025 for(std::vector< MEDCouplingAutoRefCountObjectPtr<InternalPatch> >::const_iterator it=listOfPatchesOK.begin();it!=listOfPatchesOK.end();it++)
1026 addPatch((*it)->getConstPart(),factors);
1031 * 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.
1033 void MEDCouplingCartesianAMRMeshGen::createPatchesFromCriterion(const INTERP_KERNEL::BoxSplittingOptions& bso, const DataArrayByte *criterion, const std::vector<int>& factors)
1035 if(!criterion || !criterion->isAllocated())
1036 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createPatchesFromCriterion : the criterion DataArrayByte instance must be allocated and not NULL !");
1037 std::vector<bool> crit(criterion->toVectorOfBool());//check that criterion has one component.
1038 createPatchesFromCriterion(bso,crit,factors);
1042 void MEDCouplingCartesianAMRMeshGen::removeAllPatches()
1048 void MEDCouplingCartesianAMRMeshGen::removePatch(int patchId)
1050 checkPatchId(patchId);
1051 int sz((int)_patches.size()),j(0);
1052 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> > patches(sz-1);
1053 for(int i=0;i<sz;i++)
1055 patches[j++]=_patches[i];
1056 (const_cast<MEDCouplingCartesianAMRMeshGen *>(_patches[patchId]->getMesh()))->detachFromFather();
1061 int MEDCouplingCartesianAMRMeshGen::getNumberOfPatches() const
1063 return (int)_patches.size();
1066 int MEDCouplingCartesianAMRMeshGen::getPatchIdFromChildMesh(const MEDCouplingCartesianAMRMeshGen *mesh) const
1069 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ret++)
1071 if((*it)->getMesh()==mesh)
1074 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::getPatchIdFromChildMesh : no such a mesh in my direct progeny !");
1077 std::vector< const MEDCouplingCartesianAMRPatch *> MEDCouplingCartesianAMRMeshGen::getPatches() const
1079 std::size_t sz(_patches.size());
1080 std::vector< const MEDCouplingCartesianAMRPatch *> ret(sz);
1081 for(std::size_t i=0;i<sz;i++)
1086 const MEDCouplingCartesianAMRPatch *MEDCouplingCartesianAMRMeshGen::getPatch(int patchId) const
1088 checkPatchId(patchId);
1089 return _patches[patchId];
1093 * This method states if patch2 (with id \a patchId2) is in the neighborhood of patch1 (with id \a patchId1).
1094 * The neighborhood size is defined by \a ghostLev in the reference of \a this ( \b not in the reference of patches !).
1096 bool MEDCouplingCartesianAMRMeshGen::isPatchInNeighborhoodOf(int patchId1, int patchId2, int ghostLev) const
1098 const MEDCouplingCartesianAMRPatch *p1(getPatch(patchId1)),*p2(getPatch(patchId2));
1099 return p1->isInMyNeighborhood(p2,GetGhostLevelInFineRef(ghostLev,_factors));
1103 * 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.
1104 * 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
1105 * defined by the patch with id \a patchId.
1107 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1108 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1109 * \return DataArrayDouble * - The array of the cell field on the requested patch
1111 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1112 * \throw if \a cellFieldOnThis is NULL or not allocated
1113 * \sa fillCellFieldOnPatch, MEDCouplingIMesh::SpreadCoarseToFine
1115 DataArrayDouble *MEDCouplingCartesianAMRMeshGen::createCellFieldOnPatch(int patchId, const DataArrayDouble *cellFieldOnThis) const
1117 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1118 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatch : the input cell field array is NULL or not allocated !");
1119 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1120 const MEDCouplingIMesh *fine(patch->getMesh()->getImageMesh());
1121 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(DataArrayDouble::New()); ret->alloc(fine->getNumberOfCells(),cellFieldOnThis->getNumberOfComponents());
1122 ret->copyStringInfoFrom(*cellFieldOnThis);
1123 MEDCouplingIMesh::SpreadCoarseToFine(cellFieldOnThis,_mesh->getCellGridStructure(),ret,patch->getBLTRRange(),getFactors());
1128 * This method is equivalent to MEDCouplingCartesianAMRMesh::createCellFieldOnPatch except that here instead of creating a new instance
1129 * it fills the value into the \a cellFieldOnPatch data.
1131 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1132 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1133 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1135 * \sa createCellFieldOnPatch, fillCellFieldComingFromPatch
1137 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatch(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch) const
1139 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1140 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatch : the input cell field array is NULL or not allocated !");
1141 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1142 MEDCouplingIMesh::SpreadCoarseToFine(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors());
1146 * This method is the generalization of fillCellFieldOnPatch method. This method only projects coarse to fine without considering the
1147 * potential neighbor patches covered by the ghost cells of patch with id \a patchId.
1149 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1150 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1151 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1152 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1154 * \sa fillCellFieldOnPatch, fillCellFieldOnPatchGhostAdv
1156 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchGhost(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch, int ghostLev) const
1158 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1159 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::createCellFieldOnPatchGhost : the input cell field array is NULL or not allocated !");
1160 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1161 MEDCouplingIMesh::SpreadCoarseToFineGhost(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),ghostLev);
1165 * This method is equivalent to fillCellFieldOnPatchGhost except that here \b ONLY \b the \b ghost \b zone will be updated
1166 * in \a cellFieldOnPatch.
1168 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1169 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1170 * \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.
1171 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1173 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyOnGhostZone(int patchId, const DataArrayDouble *cellFieldOnThis, DataArrayDouble *cellFieldOnPatch, int ghostLev) const
1175 if(!cellFieldOnThis || !cellFieldOnThis->isAllocated())
1176 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchOnlyOnGhostZone : the input cell field array is NULL or not allocated !");
1177 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1178 MEDCouplingIMesh::SpreadCoarseToFineGhostZone(cellFieldOnThis,_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),ghostLev);
1182 * This method is a refinement of fillCellFieldOnPatchGhost. fillCellFieldOnPatchGhost is first called.
1183 * Then for all other patches than those pointed by \a patchId that overlap the ghost zone of the patch impact the ghost zone adequately.
1185 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1186 * \param [in] cellFieldOnThis - The array of the cell field on \c this->getImageMesh() to be projected to patch having id \a patchId.
1187 * \param [in,out] cellFieldOnPatch - The array of the cell field on the requested patch to be filled.
1188 * \param [in] ghostLev - The size of the ghost zone (must be >=0 !)
1189 * \param [in] arrsOnPatches - \b WARNING arrsOnPatches[patchId] is \b NOT \b const. All others are const.
1191 * \sa fillCellFieldOnPatchOnlyGhostAdv
1193 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchGhostAdv(int patchId, const DataArrayDouble *cellFieldOnThis, int ghostLev, const std::vector<const DataArrayDouble *>& arrsOnPatches) const
1195 int nbp(getNumberOfPatches());
1196 if(nbp!=(int)arrsOnPatches.size())
1198 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchGhostAdv : there are " << nbp << " patches in this and " << arrsOnPatches.size() << " arrays in the last parameter !";
1199 throw INTERP_KERNEL::Exception(oss.str().c_str());
1201 DataArrayDouble *theFieldToFill(const_cast<DataArrayDouble *>(arrsOnPatches[patchId]));
1202 // first, do as usual
1203 fillCellFieldOnPatchGhost(patchId,cellFieldOnThis,theFieldToFill,ghostLev);
1204 fillCellFieldOnPatchOnlyGhostAdv(patchId,ghostLev,arrsOnPatches);
1208 * This method updates the patch with id \a patchId considering the only the all the patches in \a this to fill ghost zone.
1209 * So \b warning, the DataArrayDouble instance \a arrsOnPatches[patchId] is non const.
1211 * \sa getPatchIdsInTheNeighborhoodOf
1213 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyGhostAdv(int patchId, int ghostLev, const std::vector<const DataArrayDouble *>& arrsOnPatches) const
1215 int nbp(getNumberOfPatches());
1216 if(nbp!=(int)arrsOnPatches.size())
1218 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMesh::fillCellFieldOnPatchOnlyGhostAdv : there are " << nbp << " patches in this and " << arrsOnPatches.size() << " arrays in the last parameter !";
1219 throw INTERP_KERNEL::Exception(oss.str().c_str());
1221 const MEDCouplingCartesianAMRPatch *refP(getPatch(patchId));
1222 DataArrayDouble *theFieldToFill(const_cast<DataArrayDouble *>(arrsOnPatches[patchId]));
1223 std::vector<int> ids(getPatchIdsInTheNeighborhoodOf(patchId,ghostLev));
1224 for(std::vector<int>::const_iterator it=ids.begin();it!=ids.end();it++)
1226 const MEDCouplingCartesianAMRPatch *otherP(getPatch(*it));
1227 MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(ghostLev,_factors,refP,otherP,theFieldToFill,arrsOnPatches[*it]);
1231 void MEDCouplingCartesianAMRMeshGen::fillCellFieldOnPatchOnlyOnGhostZoneWith(int ghostLev, const MEDCouplingCartesianAMRPatch *patchToBeModified, const MEDCouplingCartesianAMRPatch *neighborPatch, DataArrayDouble *cellFieldOnPatch, const DataArrayDouble *cellFieldNeighbor) const
1233 MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(ghostLev,_factors,patchToBeModified,neighborPatch,cellFieldOnPatch,cellFieldNeighbor);
1237 * This method updates \a cellFieldOnThis part of values coming from the cell field \a cellFieldOnPatch lying on patch having id \a patchId.
1239 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1240 * \param [in] cellFieldOnPatch - The array of the cell field on patch with id \a patchId.
1241 * \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.
1243 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1244 * \throw if \a cellFieldOnPatch is NULL or not allocated
1245 * \sa createCellFieldOnPatch, MEDCouplingIMesh::CondenseFineToCoarse,fillCellFieldComingFromPatchGhost
1247 void MEDCouplingCartesianAMRMeshGen::fillCellFieldComingFromPatch(int patchId, const DataArrayDouble *cellFieldOnPatch, DataArrayDouble *cellFieldOnThis) const
1249 if(!cellFieldOnPatch || !cellFieldOnPatch->isAllocated())
1250 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch : the input cell field array is NULL or not allocated !");
1251 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1252 MEDCouplingIMesh::CondenseFineToCoarse(_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),cellFieldOnThis);
1256 * This method is the extension of MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch managing the ghost cells. If this
1257 * method is called with \a ghostLev equal to 0 it behaves exactly as MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatch.
1259 * \param [in] patchId - The id of the patch \a cellFieldOnThis has to be put on.
1260 * \param [in] cellFieldOnPatch - The array of the cell field on patch with id \a patchId.
1261 * \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.
1262 * \param [in] ghostLev The size of ghost zone (must be >= 0 !)
1264 * \throw if \a patchId is not in [ 0 , \c this->getNumberOfPatches() )
1265 * \throw if \a cellFieldOnPatch is NULL or not allocated
1266 * \sa fillCellFieldComingFromPatch
1268 void MEDCouplingCartesianAMRMeshGen::fillCellFieldComingFromPatchGhost(int patchId, const DataArrayDouble *cellFieldOnPatch, DataArrayDouble *cellFieldOnThis, int ghostLev) const
1270 if(!cellFieldOnPatch || !cellFieldOnPatch->isAllocated())
1271 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::fillCellFieldComingFromPatchGhost : the input cell field array is NULL or not allocated !");
1272 const MEDCouplingCartesianAMRPatch *patch(getPatch(patchId));
1273 MEDCouplingIMesh::CondenseFineToCoarseGhost(_mesh->getCellGridStructure(),cellFieldOnPatch,patch->getBLTRRange(),getFactors(),cellFieldOnThis,ghostLev);
1277 * This method finds all patches (located by their ids) that are in the neighborhood of patch with id \a patchId. The neighborhood size is
1278 * defined by ghostLev.
1280 * \param [in] patchId - the id of the considered patch.
1281 * \param [in] ghostLev - the size of the neighborhood.
1282 * \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.
1284 DataArrayInt *MEDCouplingCartesianAMRMeshGen::findPatchesInTheNeighborhoodOf(int patchId, int ghostLev) const
1286 int nbp(getNumberOfPatches());
1287 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> ret(DataArrayInt::New()); ret->alloc(0,1);
1288 for(int i=0;i<nbp;i++)
1291 if(isPatchInNeighborhoodOf(i,patchId,ghostLev))
1292 ret->pushBackSilent(i);
1297 MEDCouplingUMesh *MEDCouplingCartesianAMRMeshGen::buildUnstructured() const
1299 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> part(_mesh->buildUnstructured());
1300 std::vector<bool> bs(_mesh->getNumberOfCells(),false);
1301 std::vector<int> cgs(_mesh->getCellGridStructure());
1302 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> > msSafe(_patches.size()+1);
1304 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ii++)
1306 MEDCouplingStructuredMesh::SwitchOnIdsFrom(cgs,(*it)->getBLTRRange(),bs);
1307 msSafe[ii+1]=(*it)->getMesh()->buildUnstructured();
1309 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsOff(DataArrayInt::BuildListOfSwitchedOff(bs));
1310 msSafe[0]=static_cast<MEDCouplingUMesh *>(part->buildPartOfMySelf(eltsOff->begin(),eltsOff->end(),false));
1311 std::vector< const MEDCouplingUMesh * > ms(msSafe.size());
1312 for(std::size_t i=0;i<msSafe.size();i++)
1314 return MEDCouplingUMesh::MergeUMeshes(ms);
1318 * This method returns a mesh containing as cells that there is patches at the current level.
1319 * The patches are seen like 'boxes' that is too say the refinement will not appear here.
1321 * \return MEDCoupling1SGTUMesh * - A new object to be managed by the caller containing as cells as there are patches in \a this.
1323 MEDCoupling1SGTUMesh *MEDCouplingCartesianAMRMeshGen::buildMeshFromPatchEnvelop() const
1325 std::vector<const MEDCoupling1SGTUMesh *> cells;
1326 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > cellsSafe;
1327 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1329 const MEDCouplingCartesianAMRPatch *patch(*it);
1332 MEDCouplingAutoRefCountObjectPtr<MEDCouplingIMesh> cell(patch->getMesh()->getImageMesh()->asSingleCell());
1333 MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> cell1SGT(cell->build1SGTUnstructured());
1334 cellsSafe.push_back(cell1SGT); cells.push_back(cell1SGT);
1337 return MEDCoupling1SGTUMesh::Merge1SGTUMeshes(cells);
1340 MEDCoupling1SGTUMesh *MEDCouplingCartesianAMRMeshGen::buildMeshOfDirectChildrenOnly() const
1342 std::vector<const MEDCoupling1SGTUMesh *> patches;
1343 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> > patchesSafe;
1344 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1346 const MEDCouplingCartesianAMRPatch *patch(*it);
1349 MEDCouplingAutoRefCountObjectPtr<MEDCoupling1SGTUMesh> patchMesh(patch->getMesh()->getImageMesh()->build1SGTUnstructured());
1350 patchesSafe.push_back(patchMesh); patches.push_back(patchMesh);
1353 return MEDCoupling1SGTUMesh::Merge1SGTUMeshes(patches);
1357 * This method works same as buildUnstructured except that arrays are given in input to build a field on cell in output.
1358 * \return MEDCouplingFieldDouble * - a newly created instance the caller has reponsability to deal with.
1359 * \sa buildUnstructured
1361 MEDCouplingFieldDouble *MEDCouplingCartesianAMRMeshGen::buildCellFieldOnRecurseWithoutOverlapWithoutGhost(int ghostSz, const std::vector<const DataArrayDouble *>& recurseArrs) const
1363 if(recurseArrs.empty())
1364 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::buildCellFieldOnRecurseWithoutOverlapWithoutGhost : array is empty ! Should never happen !");
1366 std::vector<bool> bs(_mesh->getNumberOfCells(),false);
1367 std::vector<int> cgs(_mesh->getCellGridStructure());
1368 std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> > msSafe(_patches.size()+1);
1370 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++,ii++)
1372 MEDCouplingStructuredMesh::SwitchOnIdsFrom(cgs,(*it)->getBLTRRange(),bs);
1373 std::vector<const DataArrayDouble *> tmpArrs(extractSubTreeFromGlobalFlatten((*it)->getMesh(),recurseArrs));
1374 msSafe[ii+1]=(*it)->getMesh()->buildCellFieldOnRecurseWithoutOverlapWithoutGhost(ghostSz,tmpArrs);
1376 MEDCouplingAutoRefCountObjectPtr<DataArrayInt> eltsOff(DataArrayInt::BuildListOfSwitchedOff(bs));
1378 MEDCouplingAutoRefCountObjectPtr<MEDCouplingFieldDouble> ret(MEDCouplingFieldDouble::New(ON_CELLS));
1379 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> arr2(extractGhostFrom(ghostSz,recurseArrs[0]));
1380 arr2=arr2->selectByTupleIdSafe(eltsOff->begin(),eltsOff->end());
1381 ret->setArray(arr2);
1382 ret->setName(arr2->getName());
1383 MEDCouplingAutoRefCountObjectPtr<MEDCouplingUMesh> part(_mesh->buildUnstructured());
1384 MEDCouplingAutoRefCountObjectPtr<MEDCouplingMesh> mesh(part->buildPartOfMySelf(eltsOff->begin(),eltsOff->end(),false));
1388 std::vector< const MEDCouplingFieldDouble * > ms(msSafe.size());
1389 for(std::size_t i=0;i<msSafe.size();i++)
1392 return MEDCouplingFieldDouble::MergeFields(ms);
1396 * This method extracts from \arr arr the part inside \a arr by cutting the \a ghostSz external part.
1397 * \arr is expected to be an array having a number of tuples equal to \c getImageMesh()->buildWithGhost(ghostSz).
1399 DataArrayDouble *MEDCouplingCartesianAMRMeshGen::extractGhostFrom(int ghostSz, const DataArrayDouble *arr) const
1401 std::vector<int> st(_mesh->getCellGridStructure());
1402 std::vector< std::pair<int,int> > p(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(st));
1403 std::transform(st.begin(),st.end(),st.begin(),std::bind2nd(std::plus<int>(),2*ghostSz));
1404 MEDCouplingCartesianAMRPatch::ApplyGhostOnCompactFrmt(p,ghostSz);
1405 MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> ret(MEDCouplingStructuredMesh::ExtractFieldOfDoubleFrom(st,arr,p));
1410 * This method returns all the patches in \a this not equal to \a patchId that are in neighborhood of patch with id \a patchId.
1412 * \sa fillCellFieldOnPatchOnlyGhostAdv
1414 std::vector<int> MEDCouplingCartesianAMRMeshGen::getPatchIdsInTheNeighborhoodOf(int patchId, int ghostLev) const
1416 std::vector<int> ret;
1417 int nbp(getNumberOfPatches());
1419 for(int i=0;i<nbp;i++)
1422 if(isPatchInNeighborhoodOf(i,patchId,ghostLev))
1428 MEDCouplingCartesianAMRMeshGen::MEDCouplingCartesianAMRMeshGen(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1429 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop):_father(0)
1431 _mesh=MEDCouplingIMesh::New(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop);
1434 MEDCouplingCartesianAMRMeshGen::MEDCouplingCartesianAMRMeshGen(MEDCouplingCartesianAMRMeshGen *father, MEDCouplingIMesh *mesh):_father(father)
1437 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen(MEDCouplingIMesh *mesh) constructor : empty father !");
1439 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen(MEDCouplingIMesh *mesh) constructor : The input mesh is null !");
1440 mesh->checkCoherency();
1441 _mesh=mesh; _mesh->incrRef();
1444 void MEDCouplingCartesianAMRMeshGen::checkPatchId(int patchId) const
1446 int sz(getNumberOfPatches());
1447 if(patchId<0 || patchId>=sz)
1449 std::ostringstream oss; oss << "MEDCouplingCartesianAMRMeshGen::checkPatchId : invalid patchId (" << patchId << ") ! Must be in [0," << sz << ") !";
1450 throw INTERP_KERNEL::Exception(oss.str().c_str());
1454 void MEDCouplingCartesianAMRMeshGen::checkFactorsAndIfNotSetAssign(const std::vector<int>& factors)
1456 if(getSpaceDimension()!=(int)factors.size())
1457 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::checkFactorsAndIfNotSetAssign : invalid size of factors ! size must be equal to the spaceDimension !");
1458 if(_factors.empty())
1464 if(_factors!=factors)
1465 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::checkFactorsAndIfNotSetAssign : the factors ");
1469 void MEDCouplingCartesianAMRMeshGen::retrieveGridsAtInternal(int lev, std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGen> >& grids) const
1473 const MEDCouplingCartesianAMRMesh *thisc(dynamic_cast<const MEDCouplingCartesianAMRMesh *>(this));//tony
1474 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatchGF> elt(new MEDCouplingCartesianAMRPatchGF(const_cast<MEDCouplingCartesianAMRMesh *>(thisc)));
1475 grids.push_back(DynamicCastSafe<MEDCouplingCartesianAMRPatchGF,MEDCouplingCartesianAMRPatchGen>(elt));
1479 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1481 const MEDCouplingCartesianAMRPatch *pt(*it);
1484 MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> tmp1(*it);
1485 grids.push_back(DynamicCastSafe<MEDCouplingCartesianAMRPatch,MEDCouplingCartesianAMRPatchGen>(tmp1));
1491 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1493 const MEDCouplingCartesianAMRPatch *pt(*it);
1495 pt->getMesh()->retrieveGridsAtInternal(lev-1,grids);
1500 int MEDCouplingCartesianAMRMeshGen::GetGhostLevelInFineRef(int ghostLev, const std::vector<int>& factors)
1503 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::GetGhostLevelInFineRef : the ghost size must be >=0 !");
1505 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::GetGhostLevelInFineRef : no factors defined !");
1506 int ghostLevInPatchRef;
1508 ghostLevInPatchRef=0;
1511 ghostLevInPatchRef=(ghostLev-1)/factors[0]+1;
1512 for(std::size_t i=0;i<factors.size();i++)
1513 ghostLevInPatchRef=std::max(ghostLevInPatchRef,(ghostLev-1)/factors[i]+1);
1515 return ghostLevInPatchRef;
1519 * This method returns a sub set of \a all. The subset is defined by the \a head in the tree defined by \a this.
1520 * Elements in \a all are expected to be sorted from god father to most refined structure.
1522 std::vector<const DataArrayDouble *> MEDCouplingCartesianAMRMeshGen::extractSubTreeFromGlobalFlatten(const MEDCouplingCartesianAMRMeshGen *head, const std::vector<const DataArrayDouble *>& all) const
1524 int maxLev(getMaxNumberOfLevelsRelativeToThis());
1525 std::vector<const DataArrayDouble *> ret;
1526 std::vector<const MEDCouplingCartesianAMRMeshGen *> meshes(1,this);
1528 for(int i=0;i<maxLev;i++)
1530 std::vector<const MEDCouplingCartesianAMRMeshGen *> meshesTmp;
1531 for(std::vector<const MEDCouplingCartesianAMRMeshGen *>::const_iterator it=meshes.begin();it!=meshes.end();it++)
1533 if((*it)==head || head->isObjectInTheProgeny(*it))
1534 ret.push_back(all[kk]);
1536 std::vector< const MEDCouplingCartesianAMRPatch *> ps((*it)->getPatches());
1537 for(std::vector< const MEDCouplingCartesianAMRPatch *>::const_iterator it0=ps.begin();it0!=ps.end();it0++)
1539 const MEDCouplingCartesianAMRMeshGen *mesh((*it0)->getMesh());
1540 meshesTmp.push_back(mesh);
1546 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMeshGen::extractSubTreeFromGlobalFlatten : the size of input vector is not compatible with number of leaves in this !");
1550 std::size_t MEDCouplingCartesianAMRMeshGen::getHeapMemorySizeWithoutChildren() const
1552 return sizeof(MEDCouplingCartesianAMRMeshGen);
1555 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRMeshGen::getDirectChildren() const
1557 std::vector<const BigMemoryObject *> ret;
1558 if((const MEDCouplingIMesh *)_mesh)
1559 ret.push_back((const MEDCouplingIMesh *)_mesh);
1560 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1562 if((const MEDCouplingCartesianAMRPatch*)*it)
1563 ret.push_back((const MEDCouplingCartesianAMRPatch*)*it);
1568 void MEDCouplingCartesianAMRMeshGen::updateTime() const
1570 if((const MEDCouplingIMesh *)_mesh)
1571 updateTimeWith(*_mesh);
1572 for(std::vector< MEDCouplingAutoRefCountObjectPtr<MEDCouplingCartesianAMRPatch> >::const_iterator it=_patches.begin();it!=_patches.end();it++)
1574 const MEDCouplingCartesianAMRPatch *elt(*it);
1577 const MEDCouplingCartesianAMRMeshGen *mesh(elt->getMesh());
1579 updateTimeWith(*mesh);
1583 MEDCouplingCartesianAMRMeshSub::MEDCouplingCartesianAMRMeshSub(MEDCouplingCartesianAMRMeshGen *father, MEDCouplingIMesh *mesh):MEDCouplingCartesianAMRMeshGen(father,mesh)
1587 MEDCouplingCartesianAMRMesh *MEDCouplingCartesianAMRMesh::New(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1588 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop)
1590 return new MEDCouplingCartesianAMRMesh(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop);
1593 void MEDCouplingCartesianAMRMesh::setData(MEDCouplingDataForGodFather *data)
1595 MEDCouplingDataForGodFather *myData(_data);
1599 myData->changeGodFather(0);
1605 void MEDCouplingCartesianAMRMesh::allocData() const
1611 void MEDCouplingCartesianAMRMesh::deallocData() const
1617 MEDCouplingCartesianAMRMesh::MEDCouplingCartesianAMRMesh(const std::string& meshName, int spaceDim, const int *nodeStrctStart, const int *nodeStrctStop,
1618 const double *originStart, const double *originStop, const double *dxyzStart, const double *dxyzStop):MEDCouplingCartesianAMRMeshGen(meshName,spaceDim,nodeStrctStart,nodeStrctStop,originStart,originStop,dxyzStart,dxyzStop)
1622 std::vector<const BigMemoryObject *> MEDCouplingCartesianAMRMesh::getDirectChildren() const
1624 std::vector<const BigMemoryObject *> ret(MEDCouplingCartesianAMRMeshGen::getDirectChildren());
1625 const MEDCouplingDataForGodFather *pt(_data);
1631 void MEDCouplingCartesianAMRMesh::checkData() const
1633 const MEDCouplingDataForGodFather *data(_data);
1635 throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRMesh::checkData : no data set !");
1638 MEDCouplingCartesianAMRMesh::~MEDCouplingCartesianAMRMesh()
1640 MEDCouplingDataForGodFather *data(_data);
1642 data->changeGodFather(0);