#include "InterpKernelGeo2DBounds.hxx"
#include "InterpKernelGeo2DEdge.txx"
-#include "NormalizedUnstructuredMesh.hxx"
+#include "InterpolationUtils.hxx"
#include <fstream>
#include <sstream>
std::size_t nbOfSeg=std::distance(descBg,descEnd);
for(std::size_t i=0;i<nbOfSeg;i++)
{
- appendEdgeFromCrudeDataArray(i,mapp,isQuad,nodalBg,coords,descBg,descEnd,intersectEdges);
+ appendEdgeFromCrudeDataArray<ALL_FORTRAN_MODE>(descBg[i],i,std::distance(descBg,descEnd),mapp,isQuad,nodalBg,coords,intersectEdges);
}
}
-void QuadraticPolygon::appendEdgeFromCrudeDataArray(std::size_t edgePos, const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad,
- const int *nodalBg, const double *coords,
- const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges)
+template<NumberingPolicy numPol>
+void QuadraticPolygon::appendEdgeFromCrudeDataArray(int edgeIdFortOrC, std::size_t edgePos, std::size_t nbConsituents,
+ const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad,
+ const int *nodalBg, const double *coords,
+ const std::vector<std::vector<int> >& intersectEdges)
{
if(!isQuad)
{
- bool direct=descBg[edgePos]>0;
- int edgeId=abs(descBg[edgePos])-1; // back to C indexing mode
+ bool direct=false;
+ int edgeId=OTT2<int, numPol>::ind2C(edgeIdFortOrC, direct); // back to C indexing mode
const std::vector<int>& subEdge=intersectEdges[edgeId];
std::size_t nbOfSubEdges=subEdge.size()/2;
for(std::size_t j=0;j<nbOfSubEdges;j++)
}
else
{
- std::size_t nbOfSeg=std::distance(descBg,descEnd);
- const double *st=coords+2*(nodalBg[edgePos]);
+ const double *st=coords+2*(nodalBg[edgePos]);
+ const double *endd=coords+2*(nodalBg[(edgePos+1)%nbConsituents]);
+ const double *middle=coords+2*(nodalBg[edgePos+nbConsituents]);
INTERP_KERNEL::Node *st0=new INTERP_KERNEL::Node(st[0],st[1]);
- const double *endd=coords+2*(nodalBg[(edgePos+1)%nbOfSeg]);
INTERP_KERNEL::Node *endd0=new INTERP_KERNEL::Node(endd[0],endd[1]);
- const double *middle=coords+2*(nodalBg[edgePos+nbOfSeg]);
INTERP_KERNEL::Node *middle0=new INTERP_KERNEL::Node(middle[0],middle[1]);
EdgeLin *e1,*e2;
e1=new EdgeLin(st0,middle0);
bool colinearity=inters.areColinears();
delete e1; delete e2;
//
- bool direct=descBg[edgePos]>0;
- int edgeId=abs(descBg[edgePos])-1;
+ bool direct=false;
+ int edgeId=OTT2<int, numPol>::ind2C(edgeIdFortOrC, direct); // back to C indexing mode
const std::vector<int>& subEdge=intersectEdges[edgeId];
std::size_t nbOfSubEdges=subEdge.size()/2;
if(colinearity)
}
}
-void QuadraticPolygon::appendEdgeFromCrudeDataArray2(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad,
- const int *nodalBg, const double *coords,
- const int segId, const std::vector<std::vector<int> >& intersectEdges)
-{
- if(!isQuad)
- {
- bool direct=true;
- int edgeId=segId;
- const std::vector<int>& subEdge=intersectEdges[edgeId];
- std::size_t nbOfSubEdges=subEdge.size()/2;
- for(std::size_t j=0;j<nbOfSubEdges;j++)
- appendSubEdgeFromCrudeDataArray(0,j,direct,edgeId,subEdge,mapp);
- }
- else
- {
- std::size_t nbOfSeg; //=std::distance(descBg,descEnd);
- nbOfSeg = 1;
- const double *st=coords+2*(nodalBg[0]);
- INTERP_KERNEL::Node *st0=new INTERP_KERNEL::Node(st[0],st[1]);
- const double *endd=coords+2*(nodalBg[1]);
- INTERP_KERNEL::Node *endd0=new INTERP_KERNEL::Node(endd[0],endd[1]);
- const double *middle=coords+2*(nodalBg[2]);
- INTERP_KERNEL::Node *middle0=new INTERP_KERNEL::Node(middle[0],middle[1]);
- EdgeLin *e1,*e2;
- e1=new EdgeLin(st0,middle0);
- e2=new EdgeLin(middle0,endd0);
- SegSegIntersector inters(*e1,*e2);
- bool colinearity=inters.areColinears();
- delete e1; delete e2;
- //
- bool direct=true;
- int edgeId=segId;
- const std::vector<int>& subEdge=intersectEdges[edgeId];
- std::size_t nbOfSubEdges=subEdge.size()/2;
- if(colinearity)
- {
- for(std::size_t j=0;j<nbOfSubEdges;j++)
- appendSubEdgeFromCrudeDataArray(0,j,direct,edgeId,subEdge,mapp);
- }
- else
- {
- Edge *e=new EdgeArcCircle(st0,middle0,endd0,true);
- for(std::size_t j=0;j<nbOfSubEdges;j++)
- appendSubEdgeFromCrudeDataArray(e,j,direct,edgeId,subEdge,mapp);
- e->decrRef();
- }
- st0->decrRef(); endd0->decrRef(); middle0->decrRef();
- }
-}
-
-
void QuadraticPolygon::appendSubEdgeFromCrudeDataArray(Edge *baseEdge, std::size_t j, bool direct, int edgeId, const std::vector<int>& subEdge, const std::map<int,INTERP_KERNEL::Node *>& mapp)
{
std::size_t nbOfSubEdges=subEdge.size()/2;
* The information stored in intersectEdges1 and colinear1 is also used to know which intermediate nodes (intersection points with the other polygon) should be considered
* on the edges.
*/
-void QuadraticPolygon::buildFromCrudeDataArray2(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords, const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges,
+void QuadraticPolygon::buildFromCrudeDataArray2(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords,
+ const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges,
const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1,
const std::vector< std::vector<int> >& colinear1,
std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >& alreadyExistingIn2)
std::size_t nbOfSeg=std::distance(descBg,descEnd);
for(std::size_t i=0;i<nbOfSeg;i++)//loop over all edges of pol2
{
- bool direct=descBg[i]>0;
- int edgeId=abs(descBg[i])-1;//current edge id of pol2
- std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >::const_iterator it1=alreadyExistingIn2.find(descBg[i]),it2=alreadyExistingIn2.find(-descBg[i]);
- if(it1!=alreadyExistingIn2.end() || it2!=alreadyExistingIn2.end())
- {
- bool sameDir=(it1!=alreadyExistingIn2.end());
- const std::vector<INTERP_KERNEL::ElementaryEdge *>& edgesAlreadyBuilt=sameDir?(*it1).second:(*it2).second;
- if(sameDir)
- {
- for(std::vector<INTERP_KERNEL::ElementaryEdge *>::const_iterator it3=edgesAlreadyBuilt.begin();it3!=edgesAlreadyBuilt.end();it3++)
- {
- Edge *ee=(*it3)->getPtr(); ee->incrRef();
- pushBack(new ElementaryEdge(ee,(*it3)->getDirection()));
- }
- }
- else
- {
- for(std::vector<INTERP_KERNEL::ElementaryEdge *>::const_reverse_iterator it4=edgesAlreadyBuilt.rbegin();it4!=edgesAlreadyBuilt.rend();it4++)
- {
- Edge *ee=(*it4)->getPtr(); ee->incrRef();
- pushBack(new ElementaryEdge(ee,!(*it4)->getDirection()));
- }
- }
- continue;
- }
- bool directos=colinear1[edgeId].empty();
- std::vector<std::pair<int,std::pair<bool,int> > > idIns1;
- int offset1=0;
- if(!directos)
- {// if the current edge of pol2 has one or more colinear edges part into pol1
- const std::vector<int>& c=colinear1[edgeId];
- std::size_t nbOfEdgesIn1=std::distance(descBg1,descEnd1);
- for(std::size_t j=0;j<nbOfEdgesIn1;j++)
- {
- int edgeId1=abs(descBg1[j])-1;
- if(std::find(c.begin(),c.end(),edgeId1)!=c.end())
- {
- idIns1.push_back(std::pair<int,std::pair<bool,int> >(edgeId1,std::pair<bool,int>(descBg1[j]>0,offset1)));// it exists an edge into pol1 given by tuple (idIn1,direct1) that is colinear at edge 'edgeId' in pol2
- //std::pair<edgeId1); direct1=descBg1[j]>0;
- }
- offset1+=intersectEdges1[edgeId1].size()/2;//offset1 is used to find the INTERP_KERNEL::Edge * instance into pol1 that will be part of edge into pol2
- }
- directos=idIns1.empty();
- }
- if(directos)
- {//no subpart of edge 'edgeId' of pol2 is in pol1 so let's operate the same thing that QuadraticPolygon::buildFromCrudeDataArray method
- std::size_t oldSz=_sub_edges.size();
- appendEdgeFromCrudeDataArray(i,mapp,isQuad,nodalBg,coords,descBg,descEnd,intersectEdges);
- std::size_t newSz=_sub_edges.size();
- std::size_t zeSz=newSz-oldSz;
- alreadyExistingIn2[descBg[i]].resize(zeSz);
- std::list<ElementaryEdge *>::const_reverse_iterator it5=_sub_edges.rbegin();
- for(std::size_t p=0;p<zeSz;p++,it5++)
- alreadyExistingIn2[descBg[i]][zeSz-p-1]=*it5;
- }
- else
- {//there is subpart of edge 'edgeId' of pol2 inside pol1
- const std::vector<int>& subEdge=intersectEdges[edgeId];
- std::size_t nbOfSubEdges=subEdge.size()/2;
- for(std::size_t j=0;j<nbOfSubEdges;j++)
- {
- int idBg=direct?subEdge[2*j]:subEdge[2*nbOfSubEdges-2*j-1];
- int idEnd=direct?subEdge[2*j+1]:subEdge[2*nbOfSubEdges-2*j-2];
- bool direction11,found=false;
- bool direct1;//store if needed the direction in 1
- int offset2;
- std::size_t nbOfSubEdges1;
- for(std::vector<std::pair<int,std::pair<bool,int> > >::const_iterator it=idIns1.begin();it!=idIns1.end() && !found;it++)
- {
- int idIn1=(*it).first;//store if needed the cell id in 1
- direct1=(*it).second.first;
- offset1=(*it).second.second;
- const std::vector<int>& subEdge1PossiblyAlreadyIn1=intersectEdges1[idIn1];
- nbOfSubEdges1=subEdge1PossiblyAlreadyIn1.size()/2;
- offset2=0;
- for(std::size_t k=0;k<nbOfSubEdges1 && !found;k++)
- {//perform a loop on all subedges of pol1 that includes edge 'edgeId' of pol2. For the moment we iterate only on subedges of ['idIn1']... To improve
- if(subEdge1PossiblyAlreadyIn1[2*k]==idBg && subEdge1PossiblyAlreadyIn1[2*k+1]==idEnd)
- { direction11=true; found=true; }
- else if(subEdge1PossiblyAlreadyIn1[2*k]==idEnd && subEdge1PossiblyAlreadyIn1[2*k+1]==idBg)
- { direction11=false; found=true; }
- else
- offset2++;
- }
- }
- if(!found)
- {//the current subedge of edge 'edgeId' of pol2 is not a part of the colinear edge 'idIn1' of pol1 -> build new Edge instance
- //appendEdgeFromCrudeDataArray(j,mapp,isQuad,nodalBg,coords,descBg,descEnd,intersectEdges);
- Node *start=(*mapp.find(idBg)).second;
- Node *end=(*mapp.find(idEnd)).second;
- ElementaryEdge *e=ElementaryEdge::BuildEdgeFromStartEndDir(true,start,end);
- pushBack(e);
- alreadyExistingIn2[descBg[i]].push_back(e);
- }
- else
- {//the current subedge of edge 'edgeId' of pol2 is part of the colinear edge 'idIn1' of pol1 -> reuse Edge instance of pol1
- ElementaryEdge *e=pol1[offset1+(direct1?offset2:nbOfSubEdges1-offset2-1)];
- Edge *ee=e->getPtr();
- ee->incrRef();
- ElementaryEdge *e2=new ElementaryEdge(ee,!(direct1^direction11));
- pushBack(e2);
- alreadyExistingIn2[descBg[i]].push_back(e2);
- }
- }
- }
+ buildFromCrudeDataArrayGen<ALL_FORTRAN_MODE>(descBg[i],i,nbOfSeg,mapp,isQuad,nodalBg,coords,intersectEdges,pol1,descBg1,descEnd1,
+ intersectEdges1,colinear1,alreadyExistingIn2);
}
}
-void QuadraticPolygon::buildFromCrudeDataArray3(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords, const int segId, const std::vector<std::vector<int> >& intersectEdges,
+void QuadraticPolygon::buildFromCrudeDataArrayOneSeg(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords, const int segId, const std::vector<std::vector<int> >& intersectEdges,
const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1,
const std::vector< std::vector<int> >& colinear1,
std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >& alreadyExistingIn2)
{
- bool direct=true;
- int edgeId=segId;//current edge id of pol2
- std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >::const_iterator it1=alreadyExistingIn2.find(segId);
- if(it1!=alreadyExistingIn2.end())
+ buildFromCrudeDataArrayGen<ALL_C_MODE>(segId,0,2,mapp,isQuad,nodalBg,coords,intersectEdges,pol1,descBg1,descEnd1,intersectEdges1,
+ colinear1,alreadyExistingIn2);
+}
+
+template <NumberingPolicy numPol>
+void QuadraticPolygon::buildFromCrudeDataArrayGen(int edgeIdFortOrC, ssize_t edgePos, ssize_t nbConstituents, const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords,
+ const std::vector<std::vector<int> >& intersectEdges,
+ const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1,
+ const std::vector< std::vector<int> >& colinear1,
+ std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >& alreadyExistingIn2)
+{
+ bool direct;
+ int edgeId=OTT2<int, numPol>::ind2C(edgeIdFortOrC, direct); // back to C indexing mode
+ std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >::const_iterator it1=alreadyExistingIn2.find(edgeIdFortOrC),it2=alreadyExistingIn2.find(-edgeIdFortOrC);
+ if(it1!=alreadyExistingIn2.end() || it2!=alreadyExistingIn2.end())
{
- const std::vector<INTERP_KERNEL::ElementaryEdge *>& edgesAlreadyBuilt=(*it1).second;
- for(std::vector<INTERP_KERNEL::ElementaryEdge *>::const_iterator it3=edgesAlreadyBuilt.begin();it3!=edgesAlreadyBuilt.end();it3++)
+ bool sameDir=(it1!=alreadyExistingIn2.end());
+ const std::vector<INTERP_KERNEL::ElementaryEdge *>& edgesAlreadyBuilt=sameDir?(*it1).second:(*it2).second;
+ if(sameDir)
+ {
+ for(std::vector<INTERP_KERNEL::ElementaryEdge *>::const_iterator it3=edgesAlreadyBuilt.begin();it3!=edgesAlreadyBuilt.end();it3++)
+ {
+ Edge *ee=(*it3)->getPtr(); ee->incrRef();
+ pushBack(new ElementaryEdge(ee,(*it3)->getDirection()));
+ }
+ }
+ else
{
- Edge *ee=(*it3)->getPtr(); ee->incrRef();
- pushBack(new ElementaryEdge(ee,(*it3)->getDirection()));
+ for(std::vector<INTERP_KERNEL::ElementaryEdge *>::const_reverse_iterator it4=edgesAlreadyBuilt.rbegin();it4!=edgesAlreadyBuilt.rend();it4++)
+ {
+ Edge *ee=(*it4)->getPtr(); ee->incrRef();
+ pushBack(new ElementaryEdge(ee,!(*it4)->getDirection()));
+ }
}
return;
}
if(directos)
{//no subpart of edge 'edgeId' of pol2 is in pol1 so let's operate the same thing that QuadraticPolygon::buildFromCrudeDataArray method
std::size_t oldSz=_sub_edges.size();
- appendEdgeFromCrudeDataArray2(mapp,isQuad,nodalBg,coords,segId,intersectEdges);
+ appendEdgeFromCrudeDataArray<numPol>(edgeIdFortOrC,edgePos,nbConstituents,mapp,isQuad,nodalBg,coords,intersectEdges);
std::size_t newSz=_sub_edges.size();
std::size_t zeSz=newSz-oldSz;
- alreadyExistingIn2[segId].resize(zeSz);
+ alreadyExistingIn2[edgeIdFortOrC].resize(zeSz);
std::list<ElementaryEdge *>::const_reverse_iterator it5=_sub_edges.rbegin();
for(std::size_t p=0;p<zeSz;p++,it5++)
- alreadyExistingIn2[segId][zeSz-p-1]=*it5;
+ alreadyExistingIn2[edgeIdFortOrC][zeSz-p-1]=*it5;
}
else
{//there is subpart of edge 'edgeId' of pol2 inside pol1
Node *end=(*mapp.find(idEnd)).second;
ElementaryEdge *e=ElementaryEdge::BuildEdgeFromStartEndDir(true,start,end);
pushBack(e);
- alreadyExistingIn2[segId].push_back(e);
+ alreadyExistingIn2[edgeIdFortOrC].push_back(e);
}
else
{//the current subedge of edge 'edgeId' of pol2 is part of the colinear edge 'idIn1' of pol1 -> reuse Edge instance of pol1
ee->incrRef();
ElementaryEdge *e2=new ElementaryEdge(ee,!(direct1^direction11));
pushBack(e2);
- alreadyExistingIn2[segId].push_back(e2);
+ alreadyExistingIn2[edgeIdFortOrC].push_back(e2);
}
}
}
/**
* @param edgeId current edge id of pol2
*/
-void QuadraticPolygon::updateLocOfEdgeFromCrudeDataArray(const int edgeId, const std::vector<std::vector<int> >& intersectEdges,
+void QuadraticPolygon::updateLocOfOneEdgeFromCrudeDataArray(const int edgeId, bool direct, const std::vector<std::vector<int> >& intersectEdges,
const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1,
const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear1) const
{
{
for(std::size_t k=0;k<nbOfSubEdges;k++)
{
- int idBg = subEdge[2*k];
- int idEnd = subEdge[2*k+1];
+ int idBg=direct?subEdge[2*k]:subEdge[2*nbOfSubEdges-2*k-1];
+ int idEnd=direct?subEdge[2*k+1]:subEdge[2*nbOfSubEdges-2*k-2];
int idIn1=edgeId1;
bool direct1=descBg1[j]>0;
const std::vector<int>& subEdge1PossiblyAlreadyIn1=intersectEdges1[idIn1];
* Method expected to be called on pol2. Every params not suffixed by numbered are supposed to refer to pol2 (this).
* Method to find edges that are ON, and mark them correspondingly on pol1.
*/
-void QuadraticPolygon::updateLocOfEdgeFromCrudeDataArray2(const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges,
+void QuadraticPolygon::updateLocOfEdgesFromCrudeDataArray(const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges,
const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1,
const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear1) const
{
- std::size_t nbOfEdgesIn1=std::distance(descBg1,descEnd1);
std::size_t nbOfSeg=std::distance(descBg,descEnd);
for(std::size_t i=0;i<nbOfSeg;i++)//loop over all edges of pol2
{
bool direct=descBg[i]>0;
int edgeId=abs(descBg[i])-1;//current edge id of pol2
- const std::vector<int>& c=colinear1[edgeId];
- if(c.empty())
- continue;
- const std::vector<int>& subEdge=intersectEdges[edgeId];
- std::size_t nbOfSubEdges=subEdge.size()/2;
- //
- int offset1=0;
- for(std::size_t j=0;j<nbOfEdgesIn1;j++)
- {
- int edgeId1=abs(descBg1[j])-1;
- if(std::find(c.begin(),c.end(),edgeId1)!=c.end())
- {
- for(std::size_t k=0;k<nbOfSubEdges;k++)
- {
- int idBg=direct?subEdge[2*k]:subEdge[2*nbOfSubEdges-2*k-1];
- int idEnd=direct?subEdge[2*k+1]:subEdge[2*nbOfSubEdges-2*k-2];
- int idIn1=edgeId1;
- bool direct1=descBg1[j]>0;
- const std::vector<int>& subEdge1PossiblyAlreadyIn1=intersectEdges1[idIn1];
- std::size_t nbOfSubEdges1=subEdge1PossiblyAlreadyIn1.size()/2;
- int offset2=0;
- bool found=false;
- for(std::size_t kk=0;kk<nbOfSubEdges1 && !found;kk++)
- {
- found=(subEdge1PossiblyAlreadyIn1[2*kk]==idBg && subEdge1PossiblyAlreadyIn1[2*kk+1]==idEnd) || \
- (subEdge1PossiblyAlreadyIn1[2*kk]==idEnd && subEdge1PossiblyAlreadyIn1[2*kk+1]==idBg);
- if(!found)
- offset2++;
- }
- if(found)
- {
- ElementaryEdge *e=pol1[offset1+(direct1?offset2:nbOfSubEdges1-offset2-1)];
- e->getPtr()->declareOn();
- }
- }
- }
- offset1+=intersectEdges1[edgeId1].size()/2;//offset1 is used to find the INTERP_KERNEL::Edge * instance into pol1 that will be part of edge into pol2
- }
+ updateLocOfOneEdgeFromCrudeDataArray(edgeId,direct,intersectEdges,pol1,descBg1,descEnd1,intersectEdges1,colinear1);
}
}
for(std::list<ElementaryEdge *>::const_iterator it=_sub_edges.begin();it!=_sub_edges.end();it++,j++,nbOfNodesInPg++)
{
INTERP_KERNEL::Node *node=(*it)->getPtr()->buildRepresentantOfMySelf();
- //node->unApplySimilarity(xBary,yBary,fact);
addCoordsQuadratic.push_back((*node)[0]);
addCoordsQuadratic.push_back((*node)[1]);
conn.push_back(off+j);
}
}
-/** @sa ComposedEdge::initLocationsWithOther()
+/** @sa ComposedEdge::initLocationsWithOther(). Same principle but with several polygons.
+ * This method is not in ComposedEdge because of the trouble to convert a vector<QP> back into a vector<ComposedEdge> ...
*/
-void QuadraticPolygon::initLocationsWithSeveralOthers(const std::vector<QuadraticPolygon> & others) const
+void QuadraticPolygon::InitLocationsWithSeveralOthers(const QuadraticPolygon & first, const std::vector<QuadraticPolygon> & others)
{
std::set<Edge *> s1,s2;
- for(std::list<ElementaryEdge *>::const_iterator it1=_sub_edges.begin();it1!=_sub_edges.end();it1++)
+ for(std::list<ElementaryEdge *>::const_iterator it1=first._sub_edges.begin();it1!=first._sub_edges.end();it1++)
s1.insert((*it1)->getPtr());
std::vector<QuadraticPolygon>::const_iterator it_qp;
for(it_qp=others.begin(); it_qp!=others.end(); it_qp++)
for(std::list<ElementaryEdge *>::const_iterator it2=(*it_qp)._sub_edges.begin();it2!=(*it_qp)._sub_edges.end();it2++)
s2.insert((*it2)->getPtr());
- initLocations();
+ first.initLocations();
for(it_qp=others.begin(); it_qp!=others.end(); it_qp++)
(*it_qp).initLocations();
std::vector<Edge *> s3;
}
}
+void QuadraticPolygon::fullyLocateWithRespectTo(QuadraticPolygon & pol1, const int edgeId, const std::vector<std::vector<int> >& intersectEdges2,
+ const int *descBg1, const int *descEnd1,
+ const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear2)
+{
+ double xBaryBB, yBaryBB;
+ updateLocOfOneEdgeFromCrudeDataArray(edgeId,true,intersectEdges2,pol1, descBg1,descEnd1,intersectEdges1,colinear2);
+ double fact = pol1.normalizeExt(this, xBaryBB, yBaryBB);
+ pol1.performLocatingOperationSlow(*this);
+ pol1.unApplyGlobalSimilarityExt(*this,xBaryBB,yBaryBB,fact);
+}
+
+
/*!
* Given 2 polygons 'pol1' and 'pol2' (localized) the resulting polygons are returned.
*
}
}
}
+
+/**
+ * See ComputeResidual(). The only difference is the way the mapping is handled.
+ * When ComputeResidual() pushes a (-1) into the mapping, this is converted into the (c,cI) format (with cI[k] == cI[k+1]).
+ */
+void QuadraticPolygon::ComputeResidualLineIntersect(const QuadraticPolygon& pol1, const std::set<Edge *>& notUsedInPol1, const std::set<Edge *>& edgesInPol2OnBoundary,
+ const std::map<INTERP_KERNEL::Node *,int>& mapp, int offset, int idThis,
+ std::vector<double>& addCoordsQuadratic, std::vector<int>& conn, std::vector<int>& connI,
+ std::vector<int>& nb1, std::vector<int>& nb2, std::vector<int>& nbI2)
+{
+ size_t oldSz = nb2.size();
+ ComputeResidual(pol1,notUsedInPol1,edgesInPol2OnBoundary,mapp,offset,idThis,addCoordsQuadratic,conn,connI,nb1,nb2);
+ // ComputeResidual has pushed some -1 at the back of cNb2. Replace this with void entries in cNbI2:
+ // TODO: discuss the format with Anthony.
+ size_t newSz = nb2.size();
+ nb2.resize(oldSz);
+ int val = (oldSz == 0) ? 0 : nbI2.back();
+ for(size_t sss=0; sss<newSz-oldSz; sss++, nbI2.push_back(val));
+}
+
#include "InterpKernelGeo2DComposedEdge.hxx"
#include "InterpKernelGeo2DAbstractEdge.hxx"
#include "InterpKernelGeo2DElementaryEdge.hxx"
+#include "NormalizedUnstructuredMesh.hxx"
#include <list>
#include <vector>
const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1,
const std::vector< std::vector<int> >& colinear1,
std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >& alreadyExistingIn2);
- INTERPKERNEL_EXPORT void buildFromCrudeDataArray3(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords, const int segId, const std::vector<std::vector<int> >& intersectEdges,
+ INTERPKERNEL_EXPORT void buildFromCrudeDataArrayOneSeg(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords, const int segId, const std::vector<std::vector<int> >& intersectEdges,
const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1,
const std::vector< std::vector<int> >& colinear1,
std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >& alreadyExistingIn2);
- INTERPKERNEL_EXPORT void updateLocOfEdgeFromCrudeDataArray(const int edgeId, const std::vector<std::vector<int> >& intersectEdges,
- const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1,
- const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear1) const;
- INTERPKERNEL_EXPORT void updateLocOfEdgeFromCrudeDataArray2(const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges, const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear1) const;
- INTERPKERNEL_EXPORT void appendEdgeFromCrudeDataArray(std::size_t edgeId, const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords,
- const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges);
- INTERPKERNEL_EXPORT void appendEdgeFromCrudeDataArray2(const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad,
- const int *nodalBg, const double *coords,
- const int segId, const std::vector<std::vector<int> >& intersectEdges);
+ INTERPKERNEL_EXPORT void updateLocOfOneEdgeFromCrudeDataArray(const int edgeId, bool direct, const std::vector<std::vector<int> >& intersectEdges,
+ const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1,
+ const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear1) const;
+ INTERPKERNEL_EXPORT void updateLocOfEdgesFromCrudeDataArray(const int *descBg, const int *descEnd, const std::vector<std::vector<int> >& intersectEdges, const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear1) const;
INTERPKERNEL_EXPORT void appendSubEdgeFromCrudeDataArray(Edge *baseEdge, std::size_t j, bool direct, int edgeId, const std::vector<int>& subEdge, const std::map<int,INTERP_KERNEL::Node *>& mapp);
INTERPKERNEL_EXPORT void appendCrudeData(const std::map<INTERP_KERNEL::Node *,int>& mapp, int offset, std::vector<double>& addCoordsQuadratic, std::vector<int>& conn, std::vector<int>& connI) const;
INTERPKERNEL_EXPORT void buildPartitionsAbs(QuadraticPolygon& other, std::set<INTERP_KERNEL::Edge *>& edgesThis, std::set<INTERP_KERNEL::Edge *>& edgesBoundaryOther, const std::map<INTERP_KERNEL::Node *,int>& mapp, int idThis, int idOther, int offset,
int offset, std::vector<double>& addCoordsQuadratic,
std::vector<int>& conn,std::vector<int>& connI,
std::vector<int>& nbThis, std::vector<int>& nbOther, std::vector<int>& nbOtherI);
- INTERPKERNEL_EXPORT void initLocationsWithSeveralOthers(const std::vector<QuadraticPolygon> & others) const;
+ INTERPKERNEL_EXPORT static void InitLocationsWithSeveralOthers(const QuadraticPolygon & first, const std::vector<QuadraticPolygon> & others);
//
INTERPKERNEL_EXPORT double intersectWith(const QuadraticPolygon& other) const;
INTERPKERNEL_EXPORT double intersectWith(const QuadraticPolygon& other, double* barycenter) const;
public://Only public for tests reasons
INTERPKERNEL_EXPORT void performLocatingOperation(QuadraticPolygon& pol2) const;
INTERPKERNEL_EXPORT void performLocatingOperationSlow(QuadraticPolygon& pol2) const;
+ INTERPKERNEL_EXPORT void fullyLocateWithRespectTo(QuadraticPolygon & pol1, const int edgeId, const std::vector<std::vector<int> >& intersectEdges2,
+ const int *descBg1, const int *descEnd1,
+ const std::vector<std::vector<int> >& intersectEdges1, const std::vector< std::vector<int> >& colinear2);
INTERPKERNEL_EXPORT static void SplitPolygonsEachOther(QuadraticPolygon& pol1, QuadraticPolygon& pol2, int& nbOfSplits);
INTERPKERNEL_EXPORT static std::vector<QuadraticPolygon *> BuildIntersectionPolygons(const QuadraticPolygon& pol1, const QuadraticPolygon& pol2);
INTERPKERNEL_EXPORT bool haveIAChanceToBeCompletedBy(const QuadraticPolygon& pol1Splitted, const QuadraticPolygon& pol2NotSplitted, bool& direction);
INTERPKERNEL_EXPORT static void ComputeResidual(const QuadraticPolygon& pol1, const std::set<Edge *>& notUsedInPol1, const std::set<Edge *>& edgesInPol2OnBoundary, const std::map<INTERP_KERNEL::Node *,int>& mapp, int offset, int idThis,
std::vector<double>& addCoordsQuadratic, std::vector<int>& conn, std::vector<int>& connI, std::vector<int>& nb1, std::vector<int>& nb2);
+ INTERPKERNEL_EXPORT static void ComputeResidualLineIntersect(const QuadraticPolygon& pol1, const std::set<Edge *>& notUsedInPol1, const std::set<Edge *>& edgesInPol2OnBoundary,
+ const std::map<INTERP_KERNEL::Node *,int>& mapp, int offset, int idThis,
+ std::vector<double>& addCoordsQuadratic, std::vector<int>& conn, std::vector<int>& connI,
+ std::vector<int>& nb1, std::vector<int>& nb2, std::vector<int>& nbI2);
INTERPKERNEL_EXPORT static std::list<QuadraticPolygon *> ZipConsecutiveSegments2(const std::vector<int> & candidates2, const std::vector<QuadraticPolygon> & pol2s,
std::vector<std::vector<int> > & mapZipTo2);
protected:
+ /// @cond INTERNAL
+ template<NumberingPolicy numPol>
+ void appendEdgeFromCrudeDataArray(int edgeIdFortOrC, std::size_t edgePos, std::size_t nbOfSeg,
+ const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad,
+ const int *nodalBg, const double *coords,
+ const std::vector<std::vector<int> >& intersectEdges);
+ template <NumberingPolicy numPol>
+ void buildFromCrudeDataArrayGen(int edgeIdFortOrC, ssize_t edgePos, ssize_t nbConstituents, const std::map<int,INTERP_KERNEL::Node *>& mapp, bool isQuad, const int *nodalBg, const double *coords,
+ const std::vector<std::vector<int> >& intersectEdges,
+ const INTERP_KERNEL::QuadraticPolygon& pol1, const int *descBg1, const int *descEnd1, const std::vector<std::vector<int> >& intersectEdges1,
+ const std::vector< std::vector<int> >& colinear1,
+ std::map<int,std::vector<INTERP_KERNEL::ElementaryEdge *> >& alreadyExistingIn2);
+ /// @endcond
std::list<QuadraticPolygon *> zipConsecutiveInSegments() const;
void dumpInXfigFile(std::ostream& stream, int resolution, const Bounds& box) const;
static void ClosePolygons(std::list<QuadraticPolygon *>& pol2Zip, const QuadraticPolygon& pol1, const QuadraticPolygon& pol2, std::vector<QuadraticPolygon *>& results);
return ret.retn();
}
+/*!
+ * Partitions the first given 2D mesh using the second given 1D mesh as a tool, and
+ * returns a mesh made of polygons.
+ * Thus the final result contains all nodes from m1 plus new nodes. However it doesn't necessarily contains
+ * all nodes from m2.
+ * The meshes should be in 2D space. In
+ * addition, returns three arrays mapping cells of the resulting mesh to cells of the input
+ * meshes. The first array as the same format as in Intersect2DMeshes, the two others are specific since
+ * a resulting cell can be mapped to more than one cell (=a segment) in the tool mesh.
+ * \param [in] m1 - the first input mesh which is to be partioned
+ * \param [in] m2 - the second input mesh which is used as a partition tool.
+ * \param [in] eps - precision used to detect coincident mesh entities.
+ * \param [out] cellNb1 - a new instance of DataArrayInt holding for each result
+ * cell the id of the cell of \a m1 it comes from. The caller is to delete
+ * this array using decrRef() as it is no more needed.
+ * \param [out] cellNb2
+ * \param [out] cellNbI2 - for each 2D cell \a i in the result, gives the list of 1D cells from m2 where
+ * \a i comes from. This is given in surjective format (to be used with \a cellNb2). Both DataArrayInt-s are to
+ * be deleted by the caller.
+ * \return MEDCouplingUMesh * - the result 2D mesh which is a new instance of
+ * MEDCouplingUMesh. The caller is to delete this mesh using decrRef() as it
+ * is no more needed.
+ * \throw If the coordinates array is not set in any of the meshes.
+ * \throw If the nodal connectivity of cells is not defined in any of the meshes.
+ * \throw If any of the meshes is not in 2D space.
+ * \throw If m1 is not a mesh of dimension 2, or m1 is not a mesh of dimension 1
+ * \throw If m2 is not a (piecewise) line (i.e. if a point has more than 2 adjacent segments)
+ */
MEDCouplingUMesh *MEDCouplingUMesh::Intersect2DMeshWith1DLine(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2,
double eps, DataArrayInt *&cellNb1, DataArrayInt *&cellNb2, DataArrayInt *&cellNbI2)
{
ii=0;
for(std::vector<int>::const_iterator it2=candidates2.begin();it2!=candidates2.end();it2++,ii++)
{
- pol1.initLocationsWithOther(pol2s[ii]);
- pol2s[ii].updateLocOfEdgeFromCrudeDataArray2(desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,pol1,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2);
+ INTERP_KERNEL::ComposedEdge::InitLocationsWithOther(pol1,pol2s[ii]);
+ pol2s[ii].updateLocOfEdgesFromCrudeDataArray(desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,pol1,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2);
//MEDCouplingUMeshAssignOnLoc(pol1,pol2,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,desc2+descIndx2[*it2],desc2+descIndx2[*it2+1],intesctEdges2,colinear2);
pol1.buildPartitionsAbs(pol2s[ii],edges1,edgesBoundary2,mapp,i,*it2,offset3,addCoordsQuadratic,cr,crI,cNb1,cNb2);
}
MEDCouplingUMeshBuildQPFromMesh4(coo1,offset1,coo2,offset2,addCoords,*it2,intesctEdges2,
/* output */mapp,mappRev);
// pol2 is the new QP in the final merged result.
- pol2s[ii].buildFromCrudeDataArray3(mappRev,cm2.isQuadratic(),conn2+connI2[*it2]+1,coo2,*it2,intesctEdges2,
+ pol2s[ii].buildFromCrudeDataArrayOneSeg(mappRev,cm2.isQuadratic(),conn2+connI2[*it2]+1,coo2,*it2,intesctEdges2,
pol1,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2, /* output */ edgesIn2ForShare);
}
- pol1.initLocationsWithSeveralOthers(pol2s);
+ INTERP_KERNEL::QuadraticPolygon::InitLocationsWithSeveralOthers(pol1, pol2s);
+ //Locate pol2s[ii] relative to pol1 (this is done the other way around in the 2D/2D intersection)
for(it2 = candidates2.begin(), ii = 0; it2 != candidates2.end(); it2++, ii++)
- {
- pol2s[ii].updateLocOfEdgeFromCrudeDataArray(*it2 ,intesctEdges2,pol1, desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2);
- double xBaryBB, yBaryBB;
- double fact=pol1.normalizeExt(&pol2s[ii], xBaryBB, yBaryBB);
- //Locate pol2s[ii] relative to pol1 (this is done the other way around in the 2D/2D intersection)
- pol1.performLocatingOperationSlow(pol2s[ii]);
- pol1.unApplyGlobalSimilarityExt(pol2s[ii],xBaryBB,yBaryBB,fact);
- }
+ pol2s[ii].fullyLocateWithRespectTo(pol1,*it2,intesctEdges2,desc1+descIndx1[i],desc1+descIndx1[i+1],intesctEdges1,colinear2);
// At this point all edges in (all elements of) pol2s are fully localized.
// Zip consecutive IN segments from the list of candidates, and discard the pieces which are ending nowhere (e.g. an IN polyline
// with an end-point in the middle of the cell).
std::list<INTERP_KERNEL::QuadraticPolygon *> zip_list = INTERP_KERNEL::QuadraticPolygon::ZipConsecutiveSegments2(candidates2, pol2s, mapZipTo2);
// Now the core of the algo - main output is in cr, crI, cNb1, cNb2 and cNbI2.
- INTERP_KERNEL::QuadraticPolygon::BuildPartitionFromZipList(pol1, zip_list, edges1, edgesBoundary2, mapp,
- i, mapZipTo2,
+ INTERP_KERNEL::QuadraticPolygon::BuildPartitionFromZipList(pol1, zip_list, edges1, edgesBoundary2, mapp, i, mapZipTo2,
offset3,addCoordsQuadratic,cr,crI, cNb1,cNb2, cNbI2);
// Deals with remaining (non-consumed) edges from m1: these are the edges that were never touched
// by m2 but that we still want to keep in the final result.
if(!edges1.empty())
{
- size_t oldSz = cNb2.size();
try
{
- INTERP_KERNEL::QuadraticPolygon::ComputeResidual(pol1,edges1,edgesBoundary2,mapp,offset3,i,addCoordsQuadratic,cr,crI,cNb1,cNb2);
+ INTERP_KERNEL::QuadraticPolygon::ComputeResidualLineIntersect(pol1,edges1,edgesBoundary2,mapp,offset3,i,addCoordsQuadratic,cr,crI,cNb1,cNb2,cNbI2);
}
catch(INTERP_KERNEL::Exception& e)
{
std::ostringstream oss; oss << "Error when computing residual of cell #" << i << " in source/m1 mesh ! Maybe the neighbours of this cell in mesh are not well connected !\n" << "The deep reason is the following : " << e.what();
throw INTERP_KERNEL::Exception(oss.str().c_str());
}
- // ComputeResidual has pushed some -1 at the back of cNb2. Replace this with void entries in cNbI2:
- // TODO: discuss the format with Anthony.
- size_t newSz = cNb2.size();
- cNb2.resize(oldSz);
- int val = (oldSz == 0) ? 0 : cNbI2.back();
- for(size_t sss=0; sss<newSz-oldSz; sss++, cNbI2.push_back(val));
}
// Memory clean-up
for(std::map<int,INTERP_KERNEL::Node *>::const_iterator it=mappRev.begin();it!=mappRev.end();it++)
}
void MEDCouplingUMesh::IntersectDescending2DMeshWith1DMesh(const MEDCouplingUMesh *m1, const MEDCouplingUMesh *m2, double eps,
- std::vector< std::vector<int> >& intersectEdge1, std::vector< std::vector<int> >& colinear2, std::vector< std::vector<int> >& subDiv2,
- MEDCouplingUMesh *& m1Desc, DataArrayInt *&desc1, DataArrayInt *&descIndx1, DataArrayInt *&revDesc1, DataArrayInt *&revDescIndx1,
- std::vector<double>& addCoo)
+ std::vector< std::vector<int> >& intersectEdge1, std::vector< std::vector<int> >& colinear2, std::vector< std::vector<int> >& subDiv2,
+ MEDCouplingUMesh *& m1Desc, DataArrayInt *&desc1, DataArrayInt *&descIndx1, DataArrayInt *&revDesc1, DataArrayInt *&revDescIndx1,
+ std::vector<double>& addCoo)
{
static const int SPACEDIM=2;
// Build desc connectivity of first mesh only
