-// Copyright (C) 2007-2008 CEA/DEN, EDF R&D
+// Copyright (C) 2007-2013 CEA/DEN, EDF R&D
//
-// This library is free software; you can redistribute it and/or
-// modify it under the terms of the GNU Lesser General Public
-// License as published by the Free Software Foundation; either
-// version 2.1 of the License.
+// This library is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 2.1 of the License.
//
-// This library is distributed in the hope that it will be useful,
-// but WITHOUT ANY WARRANTY; without even the implied warranty of
-// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
-// Lesser General Public License for more details.
+// This library is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+// Lesser General Public License for more details.
//
-// You should have received a copy of the GNU Lesser General Public
-// License along with this library; if not, write to the Free Software
-// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+// You should have received a copy of the GNU Lesser General Public
+// License along with this library; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
//
-// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
//
+// Author : Anthony Geay (CEA/DEN)
#ifndef __GEOMETRIC2DINTERSECTOR_TXX__
#define __GEOMETRIC2DINTERSECTOR_TXX__
#include "Geometric2DIntersector.hxx"
#include "PlanarIntersectorP0P0.txx"
+#include "Planar2D1DIntersectorP0P0.txx"
#include "PlanarIntersectorP0P1.txx"
#include "PlanarIntersectorP1P0.txx"
+#include "PlanarIntersectorP1P1.txx"
+#include "PlanarIntersectorP1P0Bary.txx"
#include "CellModel.hxx"
-#include "QuadraticPolygon.hxx"
-#include "EdgeArcCircle.hxx"
-#include "EdgeLin.hxx"
-#include "Node.hxx"
+#include "InterpKernelGeo2DQuadraticPolygon.hxx"
+#include "InterpKernelGeo2DEdgeArcCircle.hxx"
+#include "InterpKernelGeo2DEdgeLin.hxx"
+#include "InterpKernelGeo2DNode.hxx"
+
+#define GEO2D_INTERSECTOR Geometric2DIntersector<MyMeshType,MyMatrix,InterpType>
+#define INTERSECTOR_TEMPLATE template<class MyMeshType, class MyMatrix, template <class MeshType, class TheMatrix, class ThisIntersector> class InterpType>
namespace INTERP_KERNEL
{
- template<class MyMeshType, class MyMatrix, template <class MeshType, class TheMatrix, class ThisIntersector> class InterpType>
- Geometric2DIntersector<MyMeshType,MyMatrix,InterpType>::Geometric2DIntersector(const MyMeshType& meshT, const MyMeshType& meshS,
- double dimCaracteristic, double medianPlane, double precision, int orientation):
- InterpType<MyMeshType,MyMatrix,Geometric2DIntersector<MyMeshType,MyMatrix,InterpType> >(meshT,meshS,dimCaracteristic, precision, medianPlane, true, orientation, 0)
+ INTERSECTOR_TEMPLATE
+ GEO2D_INTERSECTOR::Geometric2DIntersector(const MyMeshType& meshT, const MyMeshType& meshS,
+ double dimCaracteristic, double md3DSurf, double medianPlane,
+ double precision, int orientation):
+ InterpType<MyMeshType,MyMatrix,GEO2D_INTERSECTOR >(meshT,meshS,dimCaracteristic, precision, md3DSurf, medianPlane, true, orientation, 0)
{
- QUADRATIC_PLANAR::_precision=dimCaracteristic*precision;
+ QUADRATIC_PLANAR::_precision=precision;
}
- template<class MyMeshType, class MyMatrix, template <class MeshType, class TheMatrix, class ThisIntersector> class InterpType>
- double Geometric2DIntersector<MyMeshType,MyMatrix,InterpType>::intersectGeometry(ConnType icellT, ConnType icellS, ConnType nbNodesT, ConnType nbNodesS)
+ INTERSECTOR_TEMPLATE
+ double GEO2D_INTERSECTOR::intersectGeometry(ConnType icellT, ConnType icellS,
+ ConnType nbNodesT, ConnType nbNodesS)
{
int orientation = 1;
std::vector<double> CoordsT;
NormalizedCellType tS=PlanarIntersector<MyMeshType,MyMatrix>::_meshS.getTypeOfElement(icellS);
QuadraticPolygon *p1=buildPolygonFrom(CoordsT,tT);
QuadraticPolygon *p2=buildPolygonFrom(CoordsS,tS);
- double ret=p1->intersectWith(*p2);
+ double ret=p1->intersectWithAbs(*p2);
+ delete p1; delete p2;
+ return ret;
+ }
+
+ INTERSECTOR_TEMPLATE
+ double GEO2D_INTERSECTOR::intersectGeometry1D(ConnType icellT, ConnType icellS,
+ ConnType nbNodesT, ConnType nbNodesS,
+ bool& isColinear)
+ {
+ int orientation = 1;
+ std::vector<double> CoordsT;
+ std::vector<double> CoordsS;
+ PlanarIntersector<MyMeshType,MyMatrix>::getRealCoordinates(icellT,icellS,nbNodesT,nbNodesS,CoordsT,CoordsS,orientation);
+ NormalizedCellType tT=PlanarIntersector<MyMeshType,MyMatrix>::_meshT.getTypeOfElement(icellT);
+ NormalizedCellType tS=PlanarIntersector<MyMeshType,MyMatrix>::_meshS.getTypeOfElement(icellS);
+ QuadraticPolygon *p1=buildPolygonFrom(CoordsT,tT);
+ QuadraticPolygon *p2=buildPolygonOfOneEdgeFrom(CoordsS,tS);
+ double ret=p1->intersectWithAbs1D(*p2, isColinear);
delete p1; delete p2;
return ret;
}
- template<class MyMeshType, class MyMatrix, template <class MeshType, class TheMatrix, class ThisIntersector> class InterpType>
- double Geometric2DIntersector<MyMeshType,MyMatrix,InterpType>::intersectGeometryWithQuadrangle(const double *quadrangle, const std::vector<double>& sourceCoords, bool isSourceQuad)
+ INTERSECTOR_TEMPLATE
+ double GEO2D_INTERSECTOR::intersectGeometryWithQuadrangle(const double * quadrangle,
+ const std::vector<double>& sourceCoords,
+ bool isSourceQuad)
{
std::vector<Node *> nodes(4);
nodes[0]=new Node(quadrangle[0],quadrangle[1]);
std::vector<Node *> nodes2(nbOfSourceNodes);
for(int i=0;i<nbOfSourceNodes;i++)
nodes2[i]=new Node(sourceCoords[i*SPACEDIM],sourceCoords[i*SPACEDIM+1]);
- QuadraticPolygon *p1=QuadraticPolygon::buildLinearPolygon(nodes);
+ QuadraticPolygon *p1=QuadraticPolygon::BuildLinearPolygon(nodes);
QuadraticPolygon *p2;
if(!isSourceQuad)
- p2=QuadraticPolygon::buildLinearPolygon(nodes2);
+ p2=QuadraticPolygon::BuildLinearPolygon(nodes2);
else
- p2=QuadraticPolygon::buildArcCirclePolygon(nodes2);
- double ret=p1->intersectWith(*p2);
+ p2=QuadraticPolygon::BuildArcCirclePolygon(nodes2);
+ double ret=p1->intersectWithAbs(*p2);
delete p1; delete p2;
return ret;
}
- template<class MyMeshType, class MyMatrix, template <class MeshType, class TheMatrix, class ThisIntersector> class InterpType>
- QuadraticPolygon *Geometric2DIntersector<MyMeshType,MyMatrix,InterpType>::buildPolygonFrom(const std::vector<double>& coords, NormalizedCellType type)
+ INTERSECTOR_TEMPLATE
+ double GEO2D_INTERSECTOR::intersectGeometryGeneral(const std::vector<double>& targetCoords,
+ const std::vector<double>& sourceCoords)
+ {
+ int nbOfTargetNodes=targetCoords.size()/SPACEDIM;
+ std::vector<Node *> nodes(nbOfTargetNodes);
+ for(int i=0;i<nbOfTargetNodes;i++)
+ nodes[i]=new Node(targetCoords[i*SPACEDIM],targetCoords[i*SPACEDIM+1]);
+ int nbOfSourceNodes=sourceCoords.size()/SPACEDIM;
+ std::vector<Node *> nodes2(nbOfSourceNodes);
+ for(int i=0;i<nbOfSourceNodes;i++)
+ nodes2[i]=new Node(sourceCoords[i*SPACEDIM],sourceCoords[i*SPACEDIM+1]);
+ QuadraticPolygon *p1=QuadraticPolygon::BuildLinearPolygon(nodes);
+ QuadraticPolygon *p2=QuadraticPolygon::BuildLinearPolygon(nodes2);
+ double ret=p1->intersectWithAbs(*p2);
+ delete p1; delete p2;
+ return ret;
+ }
+
+ //================================================================================
+ /*!
+ * \brief Intersect a triangle and a polygon for P1P0 barycentric algorithm
+ * \param targetCell - list of coordinates of target polygon in full interlace
+ * \param targetCellQuadratic - specifies if target polygon is quadratic or not
+ * \param sourceTria - list of coordinates of source triangle
+ * \param res - coefficients a,b and c associated to nodes of sourceTria
+ */
+ //================================================================================
+
+ INTERSECTOR_TEMPLATE
+ double GEO2D_INTERSECTOR::intersectGeoBary(const std::vector<double>& targetCell,
+ bool targetCellQuadratic,
+ const double * sourceTria,
+ std::vector<double>& res)
+ {
+ std::vector<Node *> nodes(3);
+ nodes[0]=new Node(sourceTria[0*SPACEDIM],sourceTria[0*SPACEDIM+1]);
+ nodes[1]=new Node(sourceTria[1*SPACEDIM],sourceTria[1*SPACEDIM+1]);
+ nodes[2]=new Node(sourceTria[2*SPACEDIM],sourceTria[2*SPACEDIM+1]);
+ int nbOfTargetNodes=targetCell.size()/SPACEDIM;
+ std::vector<Node *> nodes2(nbOfTargetNodes);
+ for(int i=0;i<nbOfTargetNodes;i++)
+ nodes2[i]=new Node(targetCell[i*SPACEDIM],targetCell[i*SPACEDIM+1]);
+ QuadraticPolygon *p1=QuadraticPolygon::BuildLinearPolygon(nodes);
+ QuadraticPolygon *p2;
+ if(!targetCellQuadratic)
+ p2=QuadraticPolygon::BuildLinearPolygon(nodes2);
+ else
+ p2=QuadraticPolygon::BuildArcCirclePolygon(nodes2);
+ double barycenter[2];
+ double ret=p1->intersectWithAbs(*p2,barycenter);
+ delete p1; delete p2;
+ if ( ret > std::numeric_limits<double>::min() )
+ {
+ std::vector<const double* > sourceCell(3);
+ sourceCell[0] = &sourceTria[0];
+ sourceCell[1] = &sourceTria[SPACEDIM];
+ sourceCell[2] = &sourceTria[SPACEDIM*2];
+ res.resize(3);
+ barycentric_coords( sourceCell, barycenter, &res[0]);
+ res[0] *= ret;
+ res[1] *= ret;
+ res[2] *= ret;
+ }
+ else
+ {
+ ret = 0;
+ }
+ return ret;
+ }
+
+ INTERSECTOR_TEMPLATE
+ QuadraticPolygon *GEO2D_INTERSECTOR::buildPolygonFrom(const std::vector<double>& coords, NormalizedCellType type)
{
int nbNodes=coords.size()/SPACEDIM;
std::vector<Node *> nodes(nbNodes);
for(int i=0;i<nbNodes;i++)
nodes[i]=new Node(coords[i*SPACEDIM],coords[i*SPACEDIM+1]);
- if(!CellModel::getCellModel(type).isQuadratic())
- return QuadraticPolygon::buildLinearPolygon(nodes);
+ if(!CellModel::GetCellModel(type).isQuadratic())
+ return QuadraticPolygon::BuildLinearPolygon(nodes);
+ else
+ return QuadraticPolygon::BuildArcCirclePolygon(nodes);
+ }
+
+ INTERSECTOR_TEMPLATE
+ QuadraticPolygon *GEO2D_INTERSECTOR::buildPolygonOfOneEdgeFrom(const std::vector<double>& coords, NormalizedCellType type)
+ {
+ if(type==NORM_SEG2)
+ {
+ Node *node0=new Node(coords[0],coords[1]);
+ Node *node1=new Node(coords[SPACEDIM],coords[SPACEDIM+1]);
+ QuadraticPolygon *ret=new QuadraticPolygon;
+ ret->pushBack(new EdgeLin(node0,node1));
+ node0->decrRef(); node1->decrRef();
+ return ret;
+ }
+ else if(type==NORM_SEG3)
+ {
+ Node *nodeBg=new Node(coords[0],coords[1]);
+ Node *nodeEnd=new Node(coords[SPACEDIM],coords[SPACEDIM+1]);
+ Node *nodeMiddle=new Node(coords[2*SPACEDIM],coords[2*SPACEDIM+1]);
+ QuadraticPolygon *ret=new QuadraticPolygon;
+ ret->pushBack(new EdgeArcCircle(nodeBg,nodeMiddle,nodeEnd));
+ nodeBg->decrRef(); nodeEnd->decrRef(); nodeMiddle->decrRef();
+ return ret;
+ }
else
- return QuadraticPolygon::buildArcCirclePolygon(nodes);
+ throw INTERP_KERNEL::Exception("buildPolygonOfOneEdgeFrom : trying to build such non close QuadraticPolygon with 1D type !");
}
- template<class MyMeshType, class MyMatrix, template <class MeshType, class TheMatrix, class ThisIntersector> class InterpType>
- QuadraticPolygon *Geometric2DIntersector<MyMeshType,MyMatrix,InterpType>::buildPolygonAFrom(ConnType cell, int nbOfPoints, NormalizedCellType type)
+ INTERSECTOR_TEMPLATE
+ QuadraticPolygon *GEO2D_INTERSECTOR::buildPolygonAFrom(ConnType cell, int nbOfPoints, NormalizedCellType type)
{
const ConnType *startOfCellNodeConn=PlanarIntersector<MyMeshType,MyMatrix>::_connectT+OTT<ConnType,numPol>::conn2C(PlanarIntersector<MyMeshType,MyMatrix>::_connIndexT[OTT<ConnType,numPol>::ind2C(cell)]);
std::vector<Node *> nodes(nbOfPoints);
for(int i=0;i<nbOfPoints;i++)
nodes[i]=new Node(PlanarIntersector<MyMeshType,MyMatrix>::_coordsT+OTT<ConnType,numPol>::coo2C(startOfCellNodeConn[i])*SPACEDIM);
- if(CellModel::getCellModel(type).isQuadratic())
- return QuadraticPolygon::buildLinearPolygon(nodes);
+ if(CellModel::GetCellModel(type).isQuadratic())
+ return QuadraticPolygon::BuildLinearPolygon(nodes);
else
- return QuadraticPolygon::buildArcCirclePolygon(nodes);
+ return QuadraticPolygon::BuildArcCirclePolygon(nodes);
}
- template<class MyMeshType, class MyMatrix, template <class MeshType, class TheMatrix, class ThisIntersector> class InterpType>
- QuadraticPolygon *Geometric2DIntersector<MyMeshType,MyMatrix,InterpType>::buildPolygonBFrom(ConnType cell, int nbOfPoints, NormalizedCellType type)
+ INTERSECTOR_TEMPLATE
+ QuadraticPolygon *GEO2D_INTERSECTOR::buildPolygonBFrom(ConnType cell, int nbOfPoints, NormalizedCellType type)
{
const ConnType *startOfCellNodeConn=PlanarIntersector<MyMeshType,MyMatrix>::_connectS+OTT<ConnType,numPol>::conn2C(PlanarIntersector<MyMeshType,MyMatrix>::_connIndexS[OTT<ConnType,numPol>::ind2C(cell)]);
std::vector<Node *> nodes(nbOfPoints);
for(int i=0;i<nbOfPoints;i++)
nodes[i]=new Node(PlanarIntersector<MyMeshType,MyMatrix>::_coordsS+OTT<ConnType,numPol>::coo2C(startOfCellNodeConn[i])*SPACEDIM);
- if(type!=NORM_TRI6 && type!=NORM_QUAD8)
- return QuadraticPolygon::buildLinearPolygon(nodes);
+ const CellModel& cm=CellModel::GetCellModel(type);
+ if(!cm.isQuadratic())
+ return QuadraticPolygon::BuildLinearPolygon(nodes);
else
- return QuadraticPolygon::buildArcCirclePolygon(nodes);
+ return QuadraticPolygon::BuildArcCirclePolygon(nodes);
}
}