src/INTERP_KERNEL/Geometric2D/InterpKernelGeo2DComposedEdge.cxx

   1 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D
   2 //
   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.
   7 //
   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.
  12 //
  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
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19 // Author : Anthony Geay (CEA/DEN)
  20
  21 #include "InterpKernelGeo2DComposedEdge.hxx"
  22 #include "InterpKernelGeo2DElementaryEdge.hxx"
  23 #include "InterpKernelGeo2DEdgeArcCircle.hxx"
  24 #include "InterpKernelGeo2DEdgeInfLin.hxx"
  25 #include "InterpKernelException.hxx"
  26
  27 #include <algorithm>
  28 #include <memory>
  29 #include <iterator>
  30 #include <set>
  31
  32 using namespace INTERP_KERNEL;
  33
  34 ComposedEdge::ComposedEdge(const ComposedEdge& other)
  35 {
  36   for(std::list<ElementaryEdge *>::const_iterator iter=other._sub_edges.begin();iter!=other._sub_edges.end();iter++)
  37     _sub_edges.push_back((*iter)->clone());
  38 }
  39
  40 ComposedEdge::~ComposedEdge()
  41 {
  42   clearAll(_sub_edges.begin());
  43 }
  44
  45 void ComposedEdge::setValueAt(int i, Edge *e, bool direction)
  46 {
  47   std::list<ElementaryEdge*>::iterator it=_sub_edges.begin();
  48   for(int j=0;j<i;j++)
  49     it++;
  50   delete *it;
  51   *it=new ElementaryEdge(e,direction);
  52 }
  53
  54 struct AbsEdgeCmp
  55 {
  56   AbsEdgeCmp(ElementaryEdge *b):_b1(b) { }
  57   bool operator()(ElementaryEdge *a) { return a->getPtr()==_b1->getPtr();}
  58
  59   ElementaryEdge *_b1;
  60 };
  61
  62 double ComposedEdge::getCommonLengthWith(const ComposedEdge& other) const
  63 {
  64   double ret=0.;
  65   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
  66     {
  67       if(find_if(other._sub_edges.begin(),other._sub_edges.end(),AbsEdgeCmp(*iter))!=other._sub_edges.end())
  68         {
  69           const ElementaryEdge *tmp=static_cast<const ElementaryEdge *>(*iter);
  70           ret+=tmp->getCurveLength();
  71         }
  72     }
  73   return ret;
  74 }
  75
  76 void ComposedEdge::clear()
  77 {
  78   clearAll(_sub_edges.begin());
  79   _sub_edges.clear();
  80 }
  81
  82 void ComposedEdge::pushBack(Edge *edge, bool direction)
  83 {
  84   _sub_edges.push_back(new ElementaryEdge(edge,direction));
  85 }
  86
  87 void ComposedEdge::pushBack(ElementaryEdge *elem)
  88 {
  89   _sub_edges.push_back(elem);
  90 }
  91
  92 void ComposedEdge::pushBack(ComposedEdge *elem)
  93 {
  94   std::list<ElementaryEdge *> *elemsOfElem=elem->getListBehind();
  95   _sub_edges.insert(_sub_edges.end(),elemsOfElem->begin(),elemsOfElem->end());
  96 }
  97
  98 ElementaryEdge *ComposedEdge::operator[](int i) const
  99 {
 100   std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();
 101   for(int ii=0;ii<i;ii++)
 102     iter++;
 103   return *iter;
 104 }
 105
 106 void ComposedEdge::reverse()
 107 {
 108   _sub_edges.reverse();
 109   for(std::list<ElementaryEdge *>::iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 110     (*iter)->reverse();
 111 }
 112
 113 bool ComposedEdge::presenceOfOn() const
 114 {
 115   bool ret=false;
 116   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end() && !ret;iter++)
 117     ret=((*iter)->getLoc()==FULL_ON_1);
 118   return ret;
 119 }
 120
 121 bool ComposedEdge::presenceOfQuadraticEdge() const
 122 {
 123   bool ret=false;
 124   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end() && !ret;iter++)
 125     {
 126       Edge *e=(*iter)->getPtr();
 127       if(e)
 128         ret=dynamic_cast<EdgeArcCircle*>(e)!=0;
 129     }
 130   return ret;
 131 }
 132
 133 void ComposedEdge::initLocations() const
 134 {
 135   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 136     (*iter)->initLocations();
 137 }
 138
 139 ComposedEdge *ComposedEdge::clone() const
 140 {
 141   return new ComposedEdge(*this);
 142 }
 143
 144 bool ComposedEdge::isNodeIn(Node *n) const
 145 {
 146   bool ret=false;
 147   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end() && !ret;iter++)
 148     ret=(*iter)->isNodeIn(n);
 149   return ret;
 150 }
 151
 152 /*!
 153  * This method computes the area of 'this'.
 154  * By definition :
 155  * \f[
 156  * Area=\int_{Polygon} dS
 157  * \f]
 158  * Thanks to Green's theorem we have.
 159  * \f[
 160  * \int_{Polygon} x \cdot dS=\sum_{0 \leq i < nb of edges} -\int_{Edge_{i}}ydx=\sum_{0 \leq i < nb of edges} AreaOfZone_{Edge_{i}}
 161  * \f]
 162  * Where \f$ AreaOfZone_{i} \f$ is computed virtually by INTERP_KERNEL::Edge::getAreaOfZone with following formula :
 163  * \f[
 164  * AreaOfZone_{i}=\int_{Edge_{i}} -ydx
 165  * \f]
 166  */
 167 double ComposedEdge::getArea() const
 168 {
 169   double ret=0.;
 170   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 171     ret+=(*iter)->getAreaOfZone();
 172   return ret;
 173 }
 174
 175 double ComposedEdge::getPerimeter() const
 176 {
 177   double ret=0.;
 178   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 179     ret+=(*iter)->getCurveLength();
 180   return ret;
 181 }
 182
 183 double ComposedEdge::getHydraulicDiameter() const
 184 {
 185   return 4*fabs(getArea())/getPerimeter();
 186 }
 187
 188 /*!
 189  * This method computes barycenter of 'this' by returning xG in bary[0] and yG in bary[1].
 190  * By definition :
 191  * \f[
 192  * Area \cdot x_{G}=\int_{Polygon} x \cdot dS
 193  * \f]
 194  * \f[
 195  * Area \cdot y_{G}=\int_{Polygon} y \cdot dS
 196  * \f]
 197  * Thanks to Green's theorem we have.
 198  * \f[
 199  * \int_{Polygon} x \cdot dS=\sum_{0 \leq i < nb of edges} -\int_{Edge_{i}}yxdx
 200  * \f]
 201  * \f[
 202  * \int_{Polygon} y \cdot dS=\sum_{0 \leq i < nb of edges} -\int_{Edge_{i}}\frac{y^{2}}{2}dx
 203  * \f]
 204  * Area is computed using the same principle than described in INTERP_KERNEL::ComposedEdge::getArea method.
 205  * \f$ -\int_{Edge_{i}}yxdx \f$ and \f$ -\int_{Edge_{i}}\frac{y^{2}}{2}dx \f$ are computed virtually with INTERP_KERNEL::Edge::getBarycenterOfZone.
 206  */
 207 void ComposedEdge::getBarycenter(double *bary) const
 208 {
 209   bary[0]=0.;
 210   bary[1]=0.;
 211   double area=0.;
 212   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 213     {
 214       (*iter)->getBarycenterOfZone(bary);
 215       area+=(*iter)->getAreaOfZone();
 216     }
 217   bary[0]/=area;
 218   bary[1]/=area;
 219 }
 220
 221 /*!
 222  * Idem ComposedEdge::getBarycenter except that the special case where _sub_edges==1 is dealt here.
 223  */
 224 void ComposedEdge::getBarycenterGeneral(double *bary) const
 225 {
 226   if(_sub_edges.empty())
 227     throw INTERP_KERNEL::Exception("ComposedEdge::getBarycenterGeneral called on an empty polygon !");
 228   if(_sub_edges.size()>2)
 229     return getBarycenter(bary);
 230   double w;
 231   _sub_edges.back()->getBarycenter(bary,w);
 232 }
 233
 234 double ComposedEdge::normalize(ComposedEdge *other, double& xBary, double& yBary)
 235 {
 236   Bounds b;
 237   b.prepareForAggregation();
 238   fillBounds(b);
 239   other->fillBounds(b);
 240   double dimChar=b.getCaracteristicDim();
 241   b.getBarycenter(xBary,yBary);
 242   applyGlobalSimilarity(xBary,yBary,dimChar);
 243   other->applyGlobalSimilarity(xBary,yBary,dimChar);
 244   return dimChar;
 245 }
 246
 247 /*!
 248  * This method operates the opposite operation than ComposedEdge::applyGlobalSimilarity.
 249  */
 250 void ComposedEdge::unApplyGlobalSimilarityExt(ComposedEdge& other, double xBary, double yBary, double fact)
 251 {
 252   std::set<Node *> allNodes;
 253   getAllNodes(allNodes);
 254   other.getAllNodes(allNodes);
 255   for(std::set<Node *>::iterator iter=allNodes.begin();iter!=allNodes.end();iter++)
 256     (*iter)->unApplySimilarity(xBary,yBary,fact);
 257   for(std::list<ElementaryEdge *>::iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 258     (*iter)->unApplySimilarity(xBary,yBary,fact);
 259 }
 260
 261 double ComposedEdge::normalizeExt(ComposedEdge *other, double& xBary, double& yBary)
 262 {
 263   Bounds b;
 264   b.prepareForAggregation();
 265   fillBounds(b);
 266   other->fillBounds(b);
 267   double dimChar=b.getCaracteristicDim();
 268   b.getBarycenter(xBary,yBary);
 269   applyGlobalSimilarity2(other,xBary,yBary,dimChar);
 270   return dimChar;
 271 }
 272
 273 void ComposedEdge::dumpInXfigFile(std::ostream& stream, int resolution, const Bounds& box) const
 274 {
 275   stream.precision(10);
 276   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 277     (*iter)->dumpInXfigFile(stream,resolution,box);
 278 }
 279
 280 Node *ComposedEdge::getEndNode() const
 281 {
 282   return _sub_edges.back()->getEndNode();
 283 }
 284
 285 Node *ComposedEdge::getStartNode() const
 286 {
 287   return _sub_edges.front()->getStartNode();
 288 }
 289
 290 bool ComposedEdge::changeEndNodeWith(Node *node) const
 291 {
 292   return _sub_edges.back()->changeEndNodeWith(node);
 293 }
 294
 295 bool ComposedEdge::changeStartNodeWith(Node *node) const
 296 {
 297   return _sub_edges.front()->changeStartNodeWith(node);
 298 }
 299
 300 void ComposedEdge::fillBounds(Bounds& output) const
 301 {
 302   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 303     (*iter)->fillBounds(output);
 304 }
 305
 306 /*!
 307  * \b WARNING : applies similarity \b ONLY on edges without any change on Nodes. To perform a global similarity call applyGlobalSimilarity.
 308  */
 309 void ComposedEdge::applySimilarity(double xBary, double yBary, double dimChar)
 310 {
 311   for(std::list<ElementaryEdge *>::iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 312     (*iter)->applySimilarity(xBary,yBary,dimChar);
 313 }
 314
 315 /*!
 316  * Perform Similarity transformation on all elements of this Nodes and Edges.
 317  */
 318 void ComposedEdge::applyGlobalSimilarity(double xBary, double yBary, double dimChar)
 319 {
 320   std::set<Node *> allNodes;
 321   getAllNodes(allNodes);
 322   for(std::set<Node *>::iterator iter=allNodes.begin();iter!=allNodes.end();iter++)
 323     (*iter)->applySimilarity(xBary,yBary,dimChar);
 324   for(std::list<ElementaryEdge *>::iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 325     (*iter)->applySimilarity(xBary,yBary,dimChar);
 326 }
 327
 328 /*!
 329  * Perform Similarity transformation on all elements of this Nodes and Edges on 'this' and 'other'.
 330  * Nodes can be shared between 'this' and 'other'.
 331  */
 332 void ComposedEdge::applyGlobalSimilarity2(ComposedEdge *other, double xBary, double yBary, double dimChar)
 333 {
 334   std::set<Node *> allNodes;
 335   getAllNodes(allNodes);
 336   std::set<Node *> allNodes2;
 337   other->getAllNodes(allNodes2);
 338   for(std::set<Node *>::const_iterator it=allNodes2.begin();it!=allNodes2.end();it++)
 339     if(allNodes.find(*it)!=allNodes.end())
 340       (*it)->declareOn();
 341   allNodes.insert(allNodes2.begin(),allNodes2.end());
 342   for(std::set<Node *>::iterator iter=allNodes.begin();iter!=allNodes.end();iter++)
 343     (*iter)->applySimilarity(xBary,yBary,dimChar);
 344   for(std::list<ElementaryEdge *>::iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 345     (*iter)->applySimilarity(xBary,yBary,dimChar);
 346   for(std::list<ElementaryEdge *>::iterator iter=other->_sub_edges.begin();iter!=other->_sub_edges.end();iter++)
 347     (*iter)->applySimilarity(xBary,yBary,dimChar);
 348 }
 349
 350 /*!
 351  * This method append to param 'partConsidered' the part of length of subedges IN or ON.
 352  * @param partConsidered INOUT param.
 353  */
 354 void ComposedEdge::dispatchPerimeter(double& partConsidered) const
 355 {
 356   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 357     {
 358       TypeOfEdgeLocInPolygon loc=(*iter)->getLoc();
 359       if(loc==FULL_IN_1 || loc==FULL_ON_1)
 360         partConsidered+=(*iter)->getCurveLength();
 361     }
 362 }
 363
 364 /*!
 365  * Idem dispatchPerimeterExcl except that when a subedge is declared as ON this subedge is counted in commonPart.
 366  */
 367 void ComposedEdge::dispatchPerimeterExcl(double& partConsidered, double& commonPart) const
 368 {
 369   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 370     {
 371       TypeOfEdgeLocInPolygon loc=(*iter)->getLoc();
 372       if(loc==FULL_IN_1)
 373         partConsidered+=(*iter)->getCurveLength();
 374       if(loc==FULL_ON_1)
 375         commonPart+=(*iter)->getCurveLength();
 376     }
 377 }
 378
 379 void ComposedEdge::getAllNodes(std::set<Node *>& output) const
 380 {
 381   std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();
 382   for(;iter!=_sub_edges.end();iter++)
 383     (*iter)->getAllNodes(output);
 384 }
 385
 386 void ComposedEdge::getBarycenter(double *bary, double& weigh) const
 387 {
 388   weigh=0.; bary[0]=0.; bary[1]=0.;
 389   double tmp1,tmp2[2];
 390   for(std::list<ElementaryEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 391     {
 392       (*iter)->getBarycenter(tmp2,tmp1);
 393       weigh+=tmp1;
 394       bary[0]+=tmp1*tmp2[0];
 395       bary[1]+=tmp1*tmp2[1];
 396     }
 397   bary[0]/=weigh;
 398   bary[1]/=weigh;
 399 }
 400
 401 bool ComposedEdge::isInOrOut(Node *nodeToTest) const
 402 {
 403   Bounds b; b.prepareForAggregation();
 404   fillBounds(b);
 405   if(b.nearlyWhere((*nodeToTest)[0],(*nodeToTest)[1])==OUT)
 406     return false;
 407   // searching for e1
 408   std::set<Node *> nodes;
 409   getAllNodes(nodes);
 410   std::set<double> radialDistributionOfNodes;
 411   std::set<Node *>::const_iterator iter;
 412   for(iter=nodes.begin();iter!=nodes.end();iter++)
 413     radialDistributionOfNodes.insert(nodeToTest->getSlope(*(*iter)));
 414   std::vector<double> radialDistrib(radialDistributionOfNodes.begin(),radialDistributionOfNodes.end());
 415   radialDistributionOfNodes.clear();
 416   std::vector<double> radialDistrib2(radialDistrib.size());
 417   copy(radialDistrib.begin()+1,radialDistrib.end(),radialDistrib2.begin());
 418   radialDistrib2.back()=M_PI+radialDistrib.front();
 419   std::vector<double> radialDistrib3(radialDistrib.size());
 420   std::transform(radialDistrib2.begin(),radialDistrib2.end(),radialDistrib.begin(),radialDistrib3.begin(),std::minus<double>());
 421   std::vector<double>::iterator iter3=max_element(radialDistrib3.begin(),radialDistrib3.end());
 422   int i=iter3-radialDistrib3.begin();
 423   // ok for e1 - Let's go.
 424   EdgeInfLin *e1=new EdgeInfLin(nodeToTest,radialDistrib[i]+radialDistrib3[i]/2.);
 425   double ref=e1->getCharactValue(*nodeToTest);
 426   std::set< IntersectElement > inOutSwitch;
 427   for(std::list<ElementaryEdge *>::const_iterator iter4=_sub_edges.begin();iter4!=_sub_edges.end();iter4++)
 428     {
 429       ElementaryEdge *val=(*iter4);
 430       if(val)
 431         {
 432           Edge *e=val->getPtr();
 433           std::auto_ptr<EdgeIntersector> intersc(Edge::BuildIntersectorWith(e1,e));
 434           bool obviousNoIntersection,areOverlapped;
 435           intersc->areOverlappedOrOnlyColinears(0,obviousNoIntersection,areOverlapped);
 436           if(obviousNoIntersection)
 437             {
 438               continue;
 439             }
 440           if(!areOverlapped)
 441             {
 442               std::list< IntersectElement > listOfIntesc=intersc->getIntersectionsCharacteristicVal();
 443               for(std::list< IntersectElement >::iterator iter2=listOfIntesc.begin();iter2!=listOfIntesc.end();iter2++)
 444                 if((*iter2).isIncludedByBoth())
 445                   inOutSwitch.insert(*iter2);
 446               }
 447           //if overlapped we can forget
 448         }
 449       else
 450         throw Exception("Invalid use of ComposedEdge::isInOrOut : only one level supported !");
 451     }
 452   e1->decrRef();
 453   bool ret=false;
 454   for(std::set< IntersectElement >::iterator iter4=inOutSwitch.begin();iter4!=inOutSwitch.end();iter4++)
 455     {
 456       if((*iter4).getVal1()<ref)
 457         {
 458           if((*iter4).getNodeOnly()->getLoc()==ON_1)
 459             ret=!ret;
 460         }
 461       else
 462         break;
 463     }
 464   return ret;
 465 }
 466
 467 /*bool ComposedEdge::isInOrOut(Node *aNodeOn, Node *nodeToTest) const
 468 {
 469
 470   EdgeInfLin *e1=new EdgeInfLin(aNodeOn,nodeToTest);
 471   double ref=e1->getCharactValue(*nodeToTest);
 472   set< IntersectElement > inOutSwitch;
 473   for(vector<AbstractEdge *>::const_iterator iter=_sub_edges.begin();iter!=_sub_edges.end();iter++)
 474     {
 475       ElementaryEdge *val=dynamic_cast<ElementaryEdge *>(*iter);
 476       if(val)
 477         {
 478           Edge *e=val->getPtr();
 479           auto_ptr<Intersector> intersc(Edge::buildIntersectorWith(e1,e));
 480           bool obviousNoIntersection,areOverlapped;
 481           intersc->areOverlappedOrOnlyColinears(0,obviousNoIntersection,areOverlapped);
 482           if(obviousNoIntersection)
 483             {
 484               continue;
 485             }
 486           if(!areOverlapped)
 487             {
 488               list< IntersectElement > listOfIntesc=intersc->getIntersectionsCharacteristicVal();
 489               for(list< IntersectElement >::iterator iter2=listOfIntesc.begin();iter2!=listOfIntesc.end();iter2++)
 490                 if((*iter2).isIncludedByBoth())
 491                   inOutSwitch.insert(*iter2);
 492               }
 493           //if overlapped we can forget
 494         }
 495       else
 496         throw Exception("Invalid use of ComposedEdge::isInOrOut : only one level supported !");
 497     }
 498   e1->decrRef();
 499   bool ret=false;
 500   for(set< IntersectElement >::iterator iter=inOutSwitch.begin();iter!=inOutSwitch.end();iter++)
 501     {
 502       if((*iter).getVal1()<ref)
 503         {
 504           if((*iter).getNodeOnly()->getLoc()==ON_1)
 505             ret=!ret;
 506         }
 507       else
 508         break;
 509     }
 510   return ret;
 511 }*/
 512
 513 bool ComposedEdge::getDirection() const
 514 {
 515   throw Exception("ComposedEdge::getDirection : no sense");
 516 }
 517
 518 bool ComposedEdge::intresincEqCoarse(const Edge *other) const
 519 {
 520   if(_sub_edges.size()!=1)
 521     return false;
 522   return _sub_edges.front()->intresincEqCoarse(other);
 523 }
 524
 525 void ComposedEdge::clearAll(std::list<ElementaryEdge *>::iterator startToDel)
 526 {
 527   for(std::list<ElementaryEdge *>::iterator iter=startToDel;iter!=_sub_edges.end();iter++)
 528     delete (*iter);
 529 }