+ if(ghostLev<0)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt : the size of the neighborhood must be >= 0 !");
+ if(!other)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhoodExt : the input patch is NULL !");
+ int lev;
+ const MEDCouplingCartesianAMRMeshGen *com(FindCommonAncestor(this,other,lev));//check that factors are OK
+ if(lev==0)
+ return isInMyNeighborhood(other,ghostLev);
+ std::vector<int> offset(ComputeOffsetFromTwoToOne(com,lev,this,other));
+ const std::vector< std::pair<int,int> >& thisp(getBLTRRange());
+ std::vector< std::pair<int,int> > otherp(other->getBLTRRange());
+ otherp=MEDCouplingStructuredMesh::TranslateCompactFrmt(otherp,offset);
+ return IsInMyNeighborhood(ghostLev,thisp,otherp);
+}
+
+/*!
+ * This method states if \a other patch is in the neighborhood of \a this. The neighborhood zone is defined by \a ghostLev parameter
+ * the must be >= 0. This method works even if \a this and \a other does not share the same father.
+ * \a this is expected to be more refined than \a other. That is to say lev of \a this is greater than level of \a other.
+ *
+ * \param [in] other - The other patch
+ * \param [in] ghostLev - The size of the neighborhood zone.
+ *
+ * \throw if \a this or \a other are invalid (end before start).
+ * \throw if \a ghostLev is \b not >= 0 .
+ * \throw if \a this and \a other have not the same space dimension.
+ * \throw if there is not common ancestor of \a this and \a other.
+ *
+ * \sa isInMyNeighborhoodExt
+ */
+bool MEDCouplingCartesianAMRPatch::isInMyNeighborhoodDiffLev(const MEDCouplingCartesianAMRPatch *other, int ghostLev) const
+{
+ std::vector< std::pair<int,int> > thispp,otherpp;
+ std::vector<int> factors;
+ ComputeZonesOfTwoRelativeToOneDiffLev(ghostLev,this,other,thispp,otherpp,factors);
+ return IsInMyNeighborhood(ghostLev>0?1:0,thispp,otherpp);//1 not ghostLev ! It is not a bug ( I hope :) ) ! Because as \a this is a refinement of \a other ghostLev is supposed to be <= factors
+}
+
+std::vector< std::pair<int,int> > MEDCouplingCartesianAMRPatch::getBLTRRangeRelativeToGF() const
+{
+ std::vector< std::pair<int,int> > ret(_bl_tr);
+ const MEDCouplingCartesianAMRMeshGen *mesh(getMesh());
+ if(!mesh)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::getBLTRRangeRelativeToGF : not valid !");
+ const MEDCouplingCartesianAMRMeshGen *fath(mesh->getFather());
+ if(!fath)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::getBLTRRangeRelativeToGF : not valid 2 !");
+ std::vector<int> factors(fath->getFactors());
+ std::size_t sz(ret.size());
+ for(std::size_t ii=0;ii<sz;ii++)
+ {
+ ret[ii].first*=factors[ii];
+ ret[ii].second*=factors[ii];
+ }
+ const MEDCouplingCartesianAMRMeshGen *oldFather(fath);
+ fath=oldFather->getFather();
+ while(fath)
+ {
+ int pos(fath->getPatchIdFromChildMesh(oldFather));
+ const MEDCouplingCartesianAMRPatch *p(fath->getPatch(pos));
+ const std::vector< std::pair<int,int> >& tmp(p->getBLTRRange());
+ const std::vector<int>& factors2(fath->getFactors());
+ std::transform(factors.begin(),factors.end(),factors2.begin(),factors.begin(),std::multiplies<int>());
+ for(std::size_t ii=0;ii<sz;ii++)
+ {
+ int delta(ret[ii].second-ret[ii].first);
+ ret[ii].first+=tmp[ii].first*factors[ii];
+ ret[ii].second=ret[ii].first+delta;
+ }
+ oldFather=fath;
+ fath=oldFather->getFather();
+ }
+ return ret;
+}
+
+std::vector<int> MEDCouplingCartesianAMRPatch::computeCellGridSt() const
+{
+ const MEDCouplingCartesianAMRMeshGen *m(getMesh());
+ if(!m)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::computeCellGridSt : no mesh held by this !");
+ const MEDCouplingCartesianAMRMeshGen *father(m->getFather());
+ if(!father)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::computeCellGridSt : no father help by underlying mesh !");
+ const std::vector< std::pair<int,int> >& bltr(getBLTRRange());
+ const std::vector<int>& factors(father->getFactors());
+ std::vector<int> ret(MEDCouplingStructuredMesh::GetDimensionsFromCompactFrmt(bltr));
+ std::transform(ret.begin(),ret.end(),factors.begin(),ret.begin(),std::multiplies<int>());
+ return ret;
+}
+
+bool MEDCouplingCartesianAMRPatch::IsInMyNeighborhood(int ghostLev, const std::vector< std::pair<int,int> >& p1, const std::vector< std::pair<int,int> >& p2)
+{
+ std::size_t thispsize(p1.size());
+ if(thispsize!=p2.size())
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : the dimensions must be the same !");
+ for(std::size_t i=0;i<thispsize;i++)
+ {
+ const std::pair<int,int>& thispp(p1[i]);
+ const std::pair<int,int>& otherpp(p2[i]);
+ if(thispp.second<thispp.first)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
+ if(otherpp.second<otherpp.first)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::isInMyNeighborhood : this patch is invalid !");
+ if(otherpp.first==thispp.second+ghostLev-1)
+ continue;
+ if(otherpp.second+ghostLev-1==thispp.first)
+ continue;
+ int start(std::max(thispp.first,otherpp.first)),end(std::min(thispp.second,otherpp.second));
+ if(end<start)
+ return false;
+ }
+ return true;
+}
+
+/*!
+ * \sa FindNeighborsOfSubPatchesOf
+ */
+std::vector< std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > > MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOfSameLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2)
+{
+ if(!p1 || !p2)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOfSameLev : the input pointers must be not NULL !");
+ std::vector< std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > > ret;
+ std::vector< const MEDCouplingCartesianAMRPatch *> p1Work(p1->getMesh()->getPatches()),p2Work(p2->getMesh()->getPatches());
+ while(!p1Work.empty())
+ {
+ std::vector< std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *> > retTmp;
+ std::vector<const MEDCouplingCartesianAMRPatch *> p1Work2,p2Work2;
+ for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it1=p1Work.begin();it1!=p1Work.end();it1++)
+ {
+ for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it2=p2Work.begin();it2!=p2Work.end();it2++)
+ {
+ if((*it1)->isInMyNeighborhoodExt(*it2,ghostLev>0?1:0))//1 not ghostLev ! It is not a bug ( I hope :) ) ! Because as \a this is a refinement of \a other ghostLev is supposed to be <= factors
+ retTmp.push_back(std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *>(*it1,*it2));
+ }
+ std::vector<const MEDCouplingCartesianAMRPatch *> tmp1((*it1)->getMesh()->getPatches());
+ p1Work2.insert(p1Work2.end(),tmp1.begin(),tmp1.end());
+ }
+ for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it2=p2Work.begin();it2!=p2Work.end();it2++)
+ {
+ std::vector<const MEDCouplingCartesianAMRPatch *> tmp2((*it2)->getMesh()->getPatches());
+ p2Work2.insert(p2Work2.end(),tmp2.begin(),tmp2.end());
+ }
+ ret.push_back(retTmp);
+ p1Work=p1Work2;
+ p2Work=p2Work2;
+ }
+ return ret;
+}
+
+/*!
+ * This method returns all pair of patches (pa,pb) so that pb is in the neighborhood of pa (size of neighborhood is \a ghostLev).
+ * 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
+ * FindNeighborsOfSubPatchesOfSameLev.
+ *
+ * \sa FindNeighborsOfSubPatchesOfSameLev
+ */
+void MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOf(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, std::vector< std::pair<const MEDCouplingCartesianAMRPatch *, const MEDCouplingCartesianAMRPatch *> >& ret)
+{
+ if(!p1 || !p2)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::FindNeighborsOfSubPatchesOf : the input pointers must be not NULL !");
+ std::vector< const MEDCouplingCartesianAMRPatch *> p1Work(p1->getMesh()->getPatches());
+ while(!p1Work.empty())
+ {
+ std::vector<const MEDCouplingCartesianAMRPatch *> p1Work2;
+ for(std::vector<const MEDCouplingCartesianAMRPatch *>::const_iterator it0=p1Work.begin();it0!=p1Work.end();it0++)
+ {
+ if((*it0)->isInMyNeighborhoodDiffLev(p2,ghostLev))
+ ret.push_back(std::pair<const MEDCouplingCartesianAMRPatch *,const MEDCouplingCartesianAMRPatch *>(*it0,p2));
+ std::vector<const MEDCouplingCartesianAMRPatch *> tmp2((*it0)->getMesh()->getPatches());
+ p1Work2.insert(p1Work2.end(),tmp2.begin(),tmp2.end());
+ }
+ p1Work=p1Work2;
+ }
+}
+
+/*!
+ * \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.
+ *
+ * \saUpdateNeighborsOfOneWithTwoExt
+ */
+void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwo(int ghostLev, const std::vector<int>& factors, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
+{
+ const std::vector< std::pair<int,int> >& p1BLTR(p1->getBLTRRange());
+ const std::vector< std::pair<int,int> >& p2BLTR(p2->getBLTRRange());
+ UpdateNeighborsOfOneWithTwoInternal(ghostLev,factors,p1BLTR,p2BLTR,dataOnP1,dataOnP2);
+}
+
+/*!
+ * Idem than UpdateNeighborsOfOneWithTwo, except that here \a p1 and \a p2 are not sharing the same direct father.
+ *
+ * \sa UpdateNeighborsOfOneWithTwo
+ */
+void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoExt(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2)
+{
+ const std::vector< std::pair<int,int> >& p1BLTR(p1->getBLTRRange());//p1BLTR=[(10,12),(5,8)]
+ std::vector< std::pair<int,int> > p2BLTR(p2->getBLTRRange());//p2BLTR=[(0,1),(0,5)]
+ int lev(0);
+ const MEDCouplingCartesianAMRMeshGen *ca(FindCommonAncestor(p1,p2,lev));
+ std::vector<int> offset(ComputeOffsetFromTwoToOne(ca,lev,p1,p2));//[12,4]
+ p2BLTR=MEDCouplingStructuredMesh::TranslateCompactFrmt(p2BLTR,offset);//p2BLTR=[(12,13),(4,9)]
+ UpdateNeighborsOfOneWithTwoInternal(ghostLev,p1->getMesh()->getFather()->getFactors(),p1BLTR,p2BLTR,dataOnP1,dataOnP2);
+}
+
+/*!
+ * \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 !
+ */
+void MEDCouplingCartesianAMRPatch::UpdateNeighborsOfOneWithTwoMixedLev(int ghostLev, const MEDCouplingCartesianAMRPatch *p1, const MEDCouplingCartesianAMRPatch *p2, DataArrayDouble *dataOnP1, const DataArrayDouble *dataOnP2, bool isConservative)
+{
+ std::vector< std::pair<int,int> > p1pp,p2pp;
+ std::vector<int> factors;
+ ComputeZonesOfTwoRelativeToOneDiffLev(ghostLev,p1,p2,p1pp,p2pp,factors);
+ //
+ std::vector<int> dimsP2NotRefined(p2->computeCellGridSt());
+ std::vector<int> dimsP2Refined(dimsP2NotRefined);
+ std::transform(dimsP2NotRefined.begin(),dimsP2NotRefined.end(),factors.begin(),dimsP2Refined.begin(),std::multiplies<int>());
+ std::vector< std::pair<int,int> > p2RefinedAbs(MEDCouplingStructuredMesh::GetCompactFrmtFromDimensions(dimsP2NotRefined));
+ std::vector<int> dimsP2RefinedGhost(dimsP2Refined.size());
+ std::transform(dimsP2Refined.begin(),dimsP2Refined.end(),dimsP2RefinedGhost.begin(),std::bind2nd(std::plus<int>(),2*ghostLev));
+ MEDCouplingAutoRefCountObjectPtr<DataArrayDouble> fineP2(DataArrayDouble::New()); fineP2->alloc(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(dimsP2RefinedGhost),dataOnP2->getNumberOfComponents());
+ MEDCouplingIMesh::SpreadCoarseToFineGhost(dataOnP2,dimsP2NotRefined,fineP2,p2RefinedAbs,factors,ghostLev);
+ if(isConservative)
+ {
+ int fact(MEDCouplingStructuredMesh::DeduceNumberOfGivenStructure(factors));
+ std::transform(fineP2->begin(),fineP2->end(),fineP2->getPointer(),std::bind2nd(std::multiplies<double>(),1./((double)fact)));
+ }
+ //
+ UpdateNeighborsOfOneWithTwoInternal(ghostLev,p1->getMesh()->getFather()->getFactors(),p1pp,p2pp,dataOnP1,fineP2);
+}
+
+/*!
+ * \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 !
+ * 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
+ * on the same level as \a p1.
+ */
+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)
+{
+ std::vector<const MEDCouplingCartesianAMRMeshGen *> ancestorsOfThis;
+ const MEDCouplingCartesianAMRMeshGen *work(p1->getMesh()),*work2(0);
+ ancestorsOfThis.push_back(work);
+ while(work)
+ {
+ work=work->getFather();
+ if(work)
+ ancestorsOfThis.push_back(work);
+ }
+ //
+ work=p2->getMesh();
+ bool found(false);
+ std::size_t levThis(0),levOther(0);
+ while(work && !found)
+ {
+ work2=work;
+ work=work->getFather();
+ if(work)
+ {
+ levOther++;
+ std::vector<const MEDCouplingCartesianAMRMeshGen *>::iterator it(std::find(ancestorsOfThis.begin(),ancestorsOfThis.end(),work));
+ if(it!=ancestorsOfThis.end())
+ {
+ levThis=std::distance(ancestorsOfThis.begin(),it);
+ found=true;
+ }
+ }
+ }
+ if(!found)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev : no common ancestor found !");
+ if(levThis<=levOther)
+ throw INTERP_KERNEL::Exception("MEDCouplingCartesianAMRPatch::ComputeZonesOfTwoRelativeToOneDiffLev : this method is not called correctly !");
+ //
+ const MEDCouplingCartesianAMRMeshGen *comAncestor(ancestorsOfThis[levThis]);
+ int idThis(comAncestor->getPatchIdFromChildMesh(ancestorsOfThis[levThis-1])),idOther(comAncestor->getPatchIdFromChildMesh(work2));
+ const MEDCouplingCartesianAMRPatch *thisp(comAncestor->getPatch(idThis)),*otherp(comAncestor->getPatch(idOther));
+ std::vector<int> offset(ComputeOffsetFromTwoToOne(comAncestor,levOther,thisp,otherp));
+ p1Zone=thisp->getBLTRRange(); p2Zone=MEDCouplingStructuredMesh::TranslateCompactFrmt(otherp->getBLTRRange(),offset);
+ factToApplyOn2.resize(p1Zone.size()); std::fill(factToApplyOn2.begin(),factToApplyOn2.end(),1);
+ //
+ std::size_t nbOfTurn(levThis-levOther);
+ for(std::size_t i=0;i<nbOfTurn;i++)
+ {
+ std::vector< std::pair<int,int> > tmp0;
+ MEDCouplingStructuredMesh::ChangeReferenceFromGlobalOfCompactFrmt(p1Zone,p2Zone,tmp0,false);
+ p2Zone=tmp0;
+ const MEDCouplingCartesianAMRMeshGen *curAncestor(ancestorsOfThis[levThis-i]);
+ ApplyFactorsOnCompactFrmt(p2Zone,curAncestor->getFactors());
+ curAncestor=ancestorsOfThis[levThis-1-i];
+ const std::vector<int>& factors(curAncestor->getFactors());
+ std::transform(factToApplyOn2.begin(),factToApplyOn2.end(),factors.begin(),factToApplyOn2.begin(),std::multiplies<int>());
+ int tmpId(curAncestor->getPatchIdFromChildMesh(ancestorsOfThis[levThis-2-i]));
+ p1Zone=curAncestor->getPatch(tmpId)->getBLTRRange();
+ }