src/INTERP_KERNEL/PointLocator2DIntersector.txx

   1 // Copyright (C) 2007-2021  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, or (at your option) any later version.
   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 #ifndef __POINTLOCATORINTERSECTOR_TXX__
  21 #define __POINTLOCATORINTERSECTOR_TXX__
  22
  23 #include "PointLocator2DIntersector.hxx"
  24 #include "PlanarIntersectorP0P0.txx"
  25 #include "PlanarIntersectorP0P1.txx"
  26 #include "PlanarIntersectorP1P0.txx"
  27 #include "PlanarIntersectorP1P1.txx"
  28 #include "PlanarIntersectorP1P0Bary.txx"
  29 #include "CellModel.hxx"
  30
  31 #include "InterpKernelGeo2DQuadraticPolygon.hxx"
  32 #include "PointLocatorAlgos.txx"
  33
  34 #define PTLOC2D_INTERSECTOR PointLocator2DIntersector<MyMeshType,MyMatrix,InterpType>
  35 #define INTERSECTOR_TEMPLATE template<class MyMeshType, class MyMatrix, template <class MeshType, class TheMatrix, class ThisIntersector> class InterpType>
  36
  37 namespace INTERP_KERNEL
  38 {
  39   INTERSECTOR_TEMPLATE
  40   PTLOC2D_INTERSECTOR::PointLocator2DIntersector(const MyMeshType& meshT, const MyMeshType& meshS,
  41                                                double dimCaracteristic, double md3DSurf, double minDot3DSurf, double medianPlane,
  42                                                double precision, int orientation):
  43     InterpType<MyMeshType,MyMatrix,PTLOC2D_INTERSECTOR >(meshT,meshS,dimCaracteristic, precision, md3DSurf, minDot3DSurf, medianPlane, true, orientation, 0)
  44   {
  45   }
  46
  47   INTERSECTOR_TEMPLATE
  48   double PTLOC2D_INTERSECTOR::intersectGeometry(ConnType icellT,   ConnType icellS,
  49                                                 ConnType nbNodesT, ConnType nbNodesS)
  50   {
  51     int orientation = 1;
  52     std::vector<double> CoordsT;
  53     std::vector<double> CoordsS;
  54     PlanarIntersector<MyMeshType,MyMatrix>::getRealCoordinates(icellT,icellS,nbNodesT,nbNodesS,CoordsT,CoordsS,orientation);
  55     NormalizedCellType tT=PlanarIntersector<MyMeshType,MyMatrix>::_meshT.getTypeOfElement(icellT);
  56     NormalizedCellType tS=PlanarIntersector<MyMeshType,MyMatrix>::_meshS.getTypeOfElement(icellS);
  57     QuadraticPolygon *pT=buildPolygonFrom(CoordsT,tT);
  58     double baryT[SPACEDIM];
  59     pT->getBarycenterGeneral(baryT);
  60     delete pT;
  61     bool ret;
  62     double eps = InterpType<MyMeshType,MyMatrix,PTLOC2D_INTERSECTOR >::_precision;
  63     if (tS != INTERP_KERNEL::NORM_POLYGON && !CellModel::GetCellModel(tS).isQuadratic())
  64       ret = PointLocatorAlgos<MyMeshType>::isElementContainsPointAlg2DSimple(baryT,&CoordsS[0],nbNodesS,eps);
  65     else
  66       {
  67         std::vector<ConnType> conn_elem(nbNodesS);
  68         std::iota(conn_elem.begin(), conn_elem.end(), 0);
  69         ret = PointLocatorAlgos<MyMeshType>::isElementContainsPointAlgo2DPolygon(baryT, tS, &CoordsS[0], &conn_elem[0], nbNodesS, eps);
  70       }
  71     return (ret ? 1. : 0.);
  72   }
  73
  74   INTERSECTOR_TEMPLATE
  75   double PTLOC2D_INTERSECTOR::intersectGeometryWithQuadrangle(const double             * quadrangle,
  76                                                               const std::vector<double>& sourceCoords,
  77                                                               bool                       isSourceQuad)
  78   {
  79     int nbOfSourceNodes=sourceCoords.size()/SPACEDIM;
  80     std::vector<Node *> nodes2(nbOfSourceNodes);
  81     for(int i=0;i<nbOfSourceNodes;i++)
  82       nodes2[i]=new Node(sourceCoords[i*SPACEDIM],sourceCoords[i*SPACEDIM+1]);
  83     QuadraticPolygon *p2;
  84     if(!isSourceQuad)
  85       p2=QuadraticPolygon::BuildLinearPolygon(nodes2);
  86     else
  87       p2=QuadraticPolygon::BuildArcCirclePolygon(nodes2);
  88     double bary[SPACEDIM];
  89     p2->getBarycenter(bary);
  90     delete p2;
  91     if( PointLocatorAlgos<MyMeshType>::isElementContainsPointAlg2D(bary,quadrangle,4) )
  92       return 1.;
  93     return 0.;
  94   }
  95
  96   INTERSECTOR_TEMPLATE
  97   double PTLOC2D_INTERSECTOR::intersectGeometryGeneral(const std::vector<double>& targetCoords,
  98                                                        const std::vector<double>& sourceCoords)
  99   {
 100     int nbOfTargetNodes=targetCoords.size()/SPACEDIM;
 101     int nbOfSourceNodes=sourceCoords.size()/SPACEDIM;
 102     std::vector<Node *> nodes2(nbOfSourceNodes);
 103     for(int i=0;i<nbOfSourceNodes;i++)
 104       nodes2[i]=new Node(sourceCoords[i*SPACEDIM],sourceCoords[i*SPACEDIM+1]);
 105     QuadraticPolygon *p=QuadraticPolygon::BuildLinearPolygon(nodes2);
 106     double bary[SPACEDIM];
 107     p->getBarycenterGeneral(bary);
 108     delete p;
 109     if( PointLocatorAlgos<MyMeshType>::isElementContainsPointAlg2D(bary,&targetCoords[0],nbOfTargetNodes) )
 110       return 1.;
 111     return 0.;
 112   }
 113
 114   //================================================================================
 115   /*!
 116    * \brief Intersect a triangle and a polygon for P1P0 barycentric algorithm
 117    *  \param targetCell - list of coordinates of target polygon in full interlace
 118    *  \param targetCellQuadratic - specifies if target polygon is quadratic or not
 119    *  \param sourceTria - list of coordinates of source triangle
 120    *  \param res - coefficients a,b and c associated to nodes of sourceTria
 121    */
 122   //================================================================================
 123
 124   INTERSECTOR_TEMPLATE
 125   double PTLOC2D_INTERSECTOR::intersectGeoBary(const std::vector<double>& targetCell,
 126                                                bool                       targetCellQuadratic,
 127                                                const double *             sourceTria,
 128                                                std::vector<double>&       res)
 129   {
 130     throw INTERP_KERNEL::Exception("intersectGeoBary incompatible with PointLocator. Deactivate P1P0Bary to avoid the problem");
 131     return 0.;
 132   }
 133
 134   INTERSECTOR_TEMPLATE
 135   QuadraticPolygon *PTLOC2D_INTERSECTOR::buildPolygonFrom(const std::vector<double>& coords, NormalizedCellType type)
 136   {
 137     std::size_t nbNodes=coords.size()/SPACEDIM;
 138     std::vector<Node *> nodes(nbNodes);
 139     for(std::size_t i=0;i<nbNodes;i++)
 140       nodes[i]=new Node(coords[i*SPACEDIM],coords[i*SPACEDIM+1]);
 141     if(!CellModel::GetCellModel(type).isQuadratic())
 142       return QuadraticPolygon::BuildLinearPolygon(nodes);
 143     else
 144       return QuadraticPolygon::BuildArcCirclePolygon(nodes);
 145   }
 146
 147   INTERSECTOR_TEMPLATE
 148   QuadraticPolygon *PTLOC2D_INTERSECTOR::buildPolygonAFrom(ConnType cell, int nbOfPoints, NormalizedCellType type)
 149   {
 150     const ConnType *startOfCellNodeConn=PlanarIntersector<MyMeshType,MyMatrix>::_connectT+OTT<ConnType,numPol>::conn2C(PlanarIntersector<MyMeshType,MyMatrix>::_connIndexT[OTT<ConnType,numPol>::ind2C(cell)]);
 151     std::vector<Node *> nodes(nbOfPoints);
 152     for(int i=0;i<nbOfPoints;i++)
 153       nodes[i]=new Node(PlanarIntersector<MyMeshType,MyMatrix>::_coordsT+OTT<ConnType,numPol>::coo2C(startOfCellNodeConn[i])*SPACEDIM);
 154     if(CellModel::GetCellModel(type).isQuadratic())
 155       return QuadraticPolygon::BuildLinearPolygon(nodes);
 156     else
 157       return QuadraticPolygon::BuildArcCirclePolygon(nodes);
 158   }
 159
 160   INTERSECTOR_TEMPLATE
 161   QuadraticPolygon *PTLOC2D_INTERSECTOR::buildPolygonBFrom(ConnType cell, int nbOfPoints, NormalizedCellType type)
 162   {
 163     const ConnType *startOfCellNodeConn=PlanarIntersector<MyMeshType,MyMatrix>::_connectS+OTT<ConnType,numPol>::conn2C(PlanarIntersector<MyMeshType,MyMatrix>::_connIndexS[OTT<ConnType,numPol>::ind2C(cell)]);
 164     std::vector<Node *> nodes(nbOfPoints);
 165     for(int i=0;i<nbOfPoints;i++)
 166       nodes[i]=new Node(PlanarIntersector<MyMeshType,MyMatrix>::_coordsS+OTT<ConnType,numPol>::coo2C(startOfCellNodeConn[i])*SPACEDIM);
 167     if(type!=NORM_TRI6 && type!=NORM_QUAD8)
 168       return QuadraticPolygon::BuildLinearPolygon(nodes);
 169     else
 170       return QuadraticPolygon::BuildArcCirclePolygon(nodes);
 171   }
 172 }
 173
 174 #endif