--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POINTLOCATOR3DINTERSECTORP0P0_HXX__
+#define __POINTLOCATOR3DINTERSECTORP0P0_HXX__
+
+#include "Intersector3DP0P0.hxx"
+#include "NormalizedUnstructuredMesh.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class PointLocator3DIntersectorP0P0 : public Intersector3DP0P0<MyMeshType,MyMatrix>
+ {
+ public:
+ static const int SPACEDIM=MyMeshType::MY_SPACEDIM;
+ static const int MESHDIM=MyMeshType::MY_MESHDIM;
+ typedef typename MyMeshType::MyConnType ConnType;
+ static const NumberingPolicy numPol=MyMeshType::My_numPol;
+ public:
+ PointLocator3DIntersectorP0P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision);
+ ~PointLocator3DIntersectorP0P0();
+ void intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res);
+ protected:
+ double _precision;
+ };
+}
+
+#endif
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POINTLOCATOR3DINTERSECTORP0P0_TXX__
+#define __POINTLOCATOR3DINTERSECTORP0P0_TXX__
+
+#include "PointLocator3DIntersectorP0P0.hxx"
+#include "Intersector3DP0P0.txx"
+#include "MeshUtils.hxx"
+
+#include "SplitterTetra.txx"
+#include "PointLocatorAlgos.txx"
+
+namespace INTERP_KERNEL
+{
+
+ /**
+ * @param targetMesh mesh containing the target elements
+ * @param srcMesh mesh containing the source elements
+ * @param policy splitting policy to be used
+ */
+ template<class MyMeshType, class MyMatrix>
+ PointLocator3DIntersectorP0P0<MyMeshType,MyMatrix>::PointLocator3DIntersectorP0P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision):
+ Intersector3DP0P0<MyMeshType,MyMatrix>(targetMesh,srcMesh),_precision(precision)
+ {
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ PointLocator3DIntersectorP0P0<MyMeshType,MyMatrix>::~PointLocator3DIntersectorP0P0()
+ {
+ }
+
+ /**
+ *
+ * @param targetCell in C mode.
+ * @param srcCells in C mode.
+ *
+ */
+ template<class MyMeshType, class MyMatrix>
+ void PointLocator3DIntersectorP0P0<MyMeshType,MyMatrix>::intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res)
+ {
+ std::vector<double> CoordsT;
+ Intersector3DP0P0<MyMeshType,MyMatrix>::getRealTargetCoordinates(OTT<ConnType,numPol>::indFC(targetCell),CoordsT);
+ double bary[SPACEDIM];
+ calculateBarycenterDyn2<SPACEDIM>(&CoordsT[0],CoordsT.size()/SPACEDIM,bary);
+ typename MyMatrix::value_type& resRow=res[targetCell];
+ const double *coordsS=Intersector3DP0P0<MyMeshType,MyMatrix>::_src_mesh.getCoordinatesPtr();
+ for(typename std::vector<ConnType>::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++)
+ {
+ NormalizedCellType tS=Intersector3DP0P0<MyMeshType,MyMatrix>::_src_mesh.getTypeOfElement(OTT<ConnType,numPol>::indFC(*iterCellS));
+ const CellModel& cmTypeS=CellModel::getCellModel(tS);
+ std::vector<ConnType> connOfCurCellS;
+ Intersector3DP0P0<MyMeshType,MyMatrix>::getConnOfSourceCell(OTT<ConnType,numPol>::indFC(*iterCellS),connOfCurCellS);
+ if(PointLocatorAlgos<MyMeshType>::isElementContainsPointAlg3D(bary,&connOfCurCellS[0],coordsS,cmTypeS,_precision))
+ {
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(*iterCellS),1));
+ }
+ }
+ }
+}
+
+#endif
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POINTLOCATOR3DINTERSECTORP0P1_HXX__
+#define __POINTLOCATOR3DINTERSECTORP0P1_HXX__
+
+#include "Intersector3DP0P1.hxx"
+#include "NormalizedUnstructuredMesh.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class PointLocator3DIntersectorP0P1 : public Intersector3DP0P1<MyMeshType,MyMatrix>
+ {
+ public:
+ static const int SPACEDIM=MyMeshType::MY_SPACEDIM;
+ static const int MESHDIM=MyMeshType::MY_MESHDIM;
+ typedef typename MyMeshType::MyConnType ConnType;
+ static const NumberingPolicy numPol=MyMeshType::My_numPol;
+ public:
+ PointLocator3DIntersectorP0P1(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision);
+ ~PointLocator3DIntersectorP0P1();
+ void intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res);
+ protected:
+ double _precision;
+ };
+}
+
+#endif
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POINTLOCATOR3DINTERSECTORP0P1_TXX__
+#define __POINTLOCATOR3DINTERSECTORP0P1_TXX__
+
+#include "PointLocator3DIntersectorP0P1.hxx"
+#include "Intersector3DP0P1.txx"
+#include "MeshUtils.hxx"
+
+namespace INTERP_KERNEL
+{
+
+ /**
+ * @param targetMesh mesh containing the target elements
+ * @param srcMesh mesh containing the source elements
+ * @param policy splitting policy to be used
+ */
+ template<class MyMeshType, class MyMatrix>
+ PointLocator3DIntersectorP0P1<MyMeshType,MyMatrix>::PointLocator3DIntersectorP0P1(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision):
+ Intersector3DP0P1<MyMeshType,MyMatrix>(targetMesh,srcMesh),_precision(precision)
+ {
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ PointLocator3DIntersectorP0P1<MyMeshType,MyMatrix>::~PointLocator3DIntersectorP0P1()
+ {
+ }
+
+ /**
+ *
+ * @param targetCell in C mode.
+ * @param srcCells in C mode.
+ *
+ */
+ template<class MyMeshType, class MyMatrix>
+ void PointLocator3DIntersectorP0P1<MyMeshType,MyMatrix>::intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res)
+ {
+ std::vector<double> coordsTarget;
+ Intersector3DP0P1<MyMeshType,MyMatrix>::getRealTargetCoordinates(OTT<ConnType,numPol>::indFC(targetCell),coordsTarget);
+ int nbNodesT=coordsTarget.size()/SPACEDIM;
+ const double *coordsS=Intersector3DP0P1<MyMeshType,MyMatrix>::_src_mesh.getCoordinatesPtr();
+ const ConnType *startOfCellNodeConnT=Intersector3DP0P1<MyMeshType,MyMatrix>::getStartConnOfTargetCell(targetCell);
+ for(typename std::vector<ConnType>::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++)
+ {
+ NormalizedCellType tS=Intersector3DP0P1<MyMeshType,MyMatrix>::_src_mesh.getTypeOfElement(OTT<ConnType,numPol>::indFC(*iterCellS));
+ const CellModel& cmTypeS=CellModel::getCellModel(tS);
+ std::vector<ConnType> connOfCurCellS;
+ Intersector3DP0P1<MyMeshType,MyMatrix>::getConnOfSourceCell(OTT<ConnType,numPol>::indFC(*iterCellS),connOfCurCellS);
+ for(int nodeIdT=0;nodeIdT<nbNodesT;nodeIdT++)
+ {
+ if(PointLocatorAlgos<MyMeshType>::isElementContainsPointAlg3D(&coordsTarget[nodeIdT*SPACEDIM],&connOfCurCellS[0],coordsS,cmTypeS,_precision))
+ {
+ ConnType curNodeTInCmode=OTT<ConnType,numPol>::coo2C(startOfCellNodeConnT[nodeIdT]);
+ typename MyMatrix::value_type& resRow=res[curNodeTInCmode];
+ typename MyMatrix::value_type::const_iterator iterRes=resRow.find(OTT<ConnType,numPol>::indFC(*iterCellS));
+ if(iterRes==resRow.end())
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(*iterCellS),1.));
+ }
+ }
+ }
+ }
+}
+
+#endif
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POINTLOCATOR3DINTERSECTORP1P0_HXX__
+#define __POINTLOCATOR3DINTERSECTORP1P0_HXX__
+
+#include "Intersector3DP1P0.hxx"
+#include "NormalizedUnstructuredMesh.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class PointLocator3DIntersectorP1P0 : public Intersector3DP1P0<MyMeshType,MyMatrix>
+ {
+ public:
+ static const int SPACEDIM=MyMeshType::MY_SPACEDIM;
+ static const int MESHDIM=MyMeshType::MY_MESHDIM;
+ typedef typename MyMeshType::MyConnType ConnType;
+ static const NumberingPolicy numPol=MyMeshType::My_numPol;
+ public:
+ PointLocator3DIntersectorP1P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision);
+ ~PointLocator3DIntersectorP1P0();
+ void intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res);
+ protected:
+ double _precision;
+ };
+}
+
+#endif
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POINTLOCATOR3DINTERSECTORP1P0_TXX__
+#define __POINTLOCATOR3DINTERSECTORP1P0_TXX__
+
+#include "PointLocator3DIntersectorP1P0.hxx"
+#include "Intersector3DP1P0.txx"
+#include "MeshUtils.hxx"
+
+namespace INTERP_KERNEL
+{
+ /**
+ * @param targetMesh mesh containing the target elements
+ * @param srcMesh mesh containing the source elements
+ * @param policy splitting policy to be used
+ *
+ * WARNING : in _split attribute, sourceMesh and targetMesh are switched in order to fit intersectCells feature.
+ */
+ template<class MyMeshType, class MyMatrix>
+ PointLocator3DIntersectorP1P0<MyMeshType,MyMatrix>::PointLocator3DIntersectorP1P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision):Intersector3DP1P0<MyMeshType,MyMatrix>(targetMesh,srcMesh),_precision(precision)
+ {
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ PointLocator3DIntersectorP1P0<MyMeshType,MyMatrix>::~PointLocator3DIntersectorP1P0()
+ {
+ }
+
+ /**
+ * @param targetCell in C mode.
+ * @param srcCells in C mode.
+ *
+ * WARNING : for all methods on _split object source and target are switched !
+ */
+ template<class MyMeshType, class MyMatrix>
+ void PointLocator3DIntersectorP1P0<MyMeshType,MyMatrix>::intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res)
+ {
+ std::vector<double> CoordsT;
+ typename MyMatrix::value_type& resRow=res[targetCell];
+ const double *coordsS=Intersector3DP1P0<MyMeshType,MyMatrix>::_src_mesh.getCoordinatesPtr();
+ Intersector3DP1P0<MyMeshType,MyMatrix>::getRealTargetCoordinates(OTT<ConnType,numPol>::indFC(targetCell),CoordsT);
+ double baryT[SPACEDIM];
+ calculateBarycenterDyn2<SPACEDIM>(&CoordsT[0],CoordsT.size()/SPACEDIM,baryT);
+ for(typename std::vector<ConnType>::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++)
+ {
+ NormalizedCellType tS=Intersector3DP1P0<MyMeshType,MyMatrix>::_src_mesh.getTypeOfElement(OTT<ConnType,numPol>::indFC(*iterCellS));
+ if(tS!=NORM_TETRA4)
+ throw INTERP_KERNEL::Exception("Invalid source cell detected for meshdim==3. Only TETRA4 supported !");
+ const CellModel& cmTypeS=CellModel::getCellModel(tS);
+ std::vector<ConnType> connOfCurCellS;
+ Intersector3DP1P0<MyMeshType,MyMatrix>::getConnOfSourceCell(OTT<ConnType,numPol>::indFC(*iterCellS),connOfCurCellS);
+ if( PointLocatorAlgos<MyMeshType>::isElementContainsPointAlg3D(baryT,&connOfCurCellS[0],coordsS,cmTypeS,_precision) )
+ {
+ double resLoc[4];
+ std::vector<double> srcCell;
+ Intersector3DP1P0<MyMeshType,MyMatrix>::getRealSourceCoordinates(OTT<ConnType,numPol>::indFC(*iterCellS),srcCell);
+ std::vector<const double*> eap(4);
+ eap[0]=&srcCell[0]; eap[1]=&srcCell[3]; eap[2]=&srcCell[6]; eap[3]=&srcCell[9];
+ barycentric_coords(eap,baryT,resLoc);
+ const ConnType *startOfCellNodeConn=Intersector3DP1P0<MyMeshType,MyMatrix>::getStartConnOfSourceCell(*iterCellS);
+ for(int nodeIdS=0;nodeIdS<4;nodeIdS++)
+ {
+ if(fabs(resLoc[nodeIdS])>_precision)
+ {
+ ConnType curNodeSInCmode=OTT<ConnType,numPol>::coo2C(startOfCellNodeConn[nodeIdS]);
+ typename MyMatrix::value_type::const_iterator iterRes=resRow.find(OTT<ConnType,numPol>::indFC(curNodeSInCmode));
+ if(iterRes==resRow.end())
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(curNodeSInCmode),resLoc[nodeIdS]));
+ else
+ {
+ double val=(*iterRes).second+resLoc[nodeIdS];
+ resRow.erase(OTT<ConnType,numPol>::indFC(curNodeSInCmode));
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(curNodeSInCmode),val));
+ }
+ }
+ }
+ }
+ }
+ }
+}
+
+#endif
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POINTLOCATOR3DINTERSECTORP1P1_HXX__
+#define __POINTLOCATOR3DINTERSECTORP1P1_HXX__
+
+#include "Intersector3DP1P1.hxx"
+#include "NormalizedUnstructuredMesh.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class PointLocator3DIntersectorP1P1 : public Intersector3DP1P1<MyMeshType,MyMatrix>
+ {
+ public:
+ static const int SPACEDIM=MyMeshType::MY_SPACEDIM;
+ static const int MESHDIM=MyMeshType::MY_MESHDIM;
+ typedef typename MyMeshType::MyConnType ConnType;
+ static const NumberingPolicy numPol=MyMeshType::My_numPol;
+ public:
+ PointLocator3DIntersectorP1P1(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision);
+ ~PointLocator3DIntersectorP1P1();
+ void intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res);
+ protected:
+ double _precision;
+ };
+}
+
+#endif
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POINTLOCATOR3DINTERSECTORP1P1_TXX__
+#define __POINTLOCATOR3DINTERSECTORP1P1_TXX__
+
+#include "PointLocator3DIntersectorP1P1.hxx"
+#include "Intersector3DP1P1.txx"
+#include "MeshUtils.hxx"
+
+namespace INTERP_KERNEL
+{
+
+ /**
+ * Constructor creating object from target cell global number
+ *
+ * @param targetMesh mesh containing the target elements
+ * @param srcMesh mesh containing the source elements
+ * @param policy splitting policy to be used
+ */
+ template<class MyMeshType, class MyMatrix>
+ PointLocator3DIntersectorP1P1<MyMeshType,MyMatrix>::PointLocator3DIntersectorP1P1(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision):
+ Intersector3DP1P1<MyMeshType,MyMatrix>(targetMesh,srcMesh),_precision(precision)
+ {
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ PointLocator3DIntersectorP1P1<MyMeshType,MyMatrix>::~PointLocator3DIntersectorP1P1()
+ {
+ }
+
+ /**
+ * @param targetCell in C mode.
+ * @param srcCells in C mode.
+ */
+ template<class MyMeshType, class MyMatrix>
+ void PointLocator3DIntersectorP1P1<MyMeshType,MyMatrix>::intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res)
+ {
+ std::vector<double> CoordsT;
+ Intersector3DP1P1<MyMeshType,MyMatrix>::getRealTargetCoordinates(OTT<ConnType,numPol>::indFC(targetCell),CoordsT);
+ int nbOfNodesT=CoordsT.size()/SPACEDIM;
+ const double *coordsS=Intersector3DP1P1<MyMeshType,MyMatrix>::_src_mesh.getCoordinatesPtr();
+ for(typename std::vector<ConnType>::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++)
+ {
+ NormalizedCellType tS=Intersector3DP1P1<MyMeshType,MyMatrix>::_src_mesh.getTypeOfElement(OTT<ConnType,numPol>::indFC(*iterCellS));
+ if(tS!=NORM_TETRA4)
+ throw INTERP_KERNEL::Exception("Invalid source cell detected for meshdim==3. Only TETRA4 supported !");
+ const CellModel& cmTypeS=CellModel::getCellModel(tS);
+ const ConnType *startOfCellNodeConnT=Intersector3DP1P1<MyMeshType,MyMatrix>::getStartConnOfTargetCell(targetCell);
+ for(int nodeIdT=0;nodeIdT<nbOfNodesT;nodeIdT++)
+ {
+ typename MyMatrix::value_type& resRow=res[OTT<ConnType,numPol>::ind2C(startOfCellNodeConnT[nodeIdT])];
+ std::vector<ConnType> connOfCurCellS;
+ Intersector3DP1P1<MyMeshType,MyMatrix>::getConnOfSourceCell(OTT<ConnType,numPol>::indFC(*iterCellS),connOfCurCellS);
+ if( PointLocatorAlgos<MyMeshType>::isElementContainsPointAlg3D(&CoordsT[nodeIdT*SPACEDIM],&connOfCurCellS[0],coordsS,cmTypeS,_precision) )
+ {
+ double resLoc[4];
+ std::vector<double> localCoordsS;
+ Intersector3DP1P1<MyMeshType,MyMatrix>::getRealSourceCoordinates(OTT<ConnType,numPol>::indFC(*iterCellS),localCoordsS);
+ std::vector<const double*> eap(4);
+ eap[0]=&localCoordsS[0]; eap[1]=&localCoordsS[3]; eap[2]=&localCoordsS[6]; eap[3]=&localCoordsS[9];
+ barycentric_coords(eap,&CoordsT[nodeIdT*SPACEDIM],resLoc);
+ const ConnType *startOfCellNodeConnS=Intersector3DP1P1<MyMeshType,MyMatrix>::getStartConnOfSourceCell(*iterCellS);
+ for(int nodeIdS=0;nodeIdS<4;nodeIdS++)
+ {
+ if(fabs(resLoc[nodeIdS])>_precision)
+ {
+ ConnType curNodeSInCmode=OTT<ConnType,numPol>::coo2C(startOfCellNodeConnS[nodeIdS]);
+ typename MyMatrix::value_type::const_iterator iterRes=resRow.find(OTT<ConnType,numPol>::indFC(curNodeSInCmode));
+ if(iterRes==resRow.end())
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(curNodeSInCmode),resLoc[nodeIdS]));
+ else
+ {
+ double val=(*iterRes).second+resLoc[nodeIdS];
+ resRow.erase(OTT<ConnType,numPol>::indFC(curNodeSInCmode));
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(curNodeSInCmode),val));
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+}
+
+#endif
return retlist;
}
- static bool isElementContainsPointAlg2D(const double *ptToTest, const double *cellPts, int nbEdges)
+ static bool isElementContainsPointAlg2D(const double *ptToTest, const double *cellPts, int nbEdges, double eps)
{
- // with dimension 2, it suffices to check all the edges
- // and see if the sign of double products from the point
- //is always the same.
- // C
- // / \
- // / \
- // Xo / \
- // A-------B
- //
- //here XA^XC and XC^XB have different signs
+ /* with dimension 2, it suffices to check all the edges
+ and see if the sign of double products from the point
+ is always the same.
+ C
+ / \
+ / \
+ Xo / \
+ A-------B
+
+ here XA^XC and XC^XB have different signs*/
const int SPACEDIM=MyMeshType::MY_SPACEDIM;
int* sign = new int[nbEdges];
for (int iedge=0; iedge<nbEdges; iedge++)
const double* A=cellPts+SPACEDIM*iedge;
const double* B=cellPts+SPACEDIM*((iedge+1)%nbEdges);
double a=mon_determinant(ptToTest, A, B);
- if(a<-1e-12)
+ if(a<-eps)
sign[iedge]=-1;
- else if(a>1e-12)
+ else if(a>eps)
sign[iedge]=1;
else
sign[iedge]=0;
return ret;
}
- static bool isElementContainsPointAlg3D(const double *ptToTest, const double *cellPts, int nbEdges)
+ static bool isElementContainsPointAlg3D(const double *ptToTest, const int *conn_elem, const double *coords, const CellModel& cmType, double eps)
{
- return true;
+ const int SPACEDIM=MyMeshType::MY_SPACEDIM;
+ typedef typename MyMeshType::MyConnType ConnType;
+ const NumberingPolicy numPol=MyMeshType::My_numPol;
+
+ int nbfaces = cmType.getNumberOfSons();
+ int* sign = new int[nbfaces];
+ for (int iface=0; iface<nbfaces; iface++)
+ {
+ const unsigned* connface=cmType.getNodesConstituentTheSon(iface);
+ const double* AA=coords+SPACEDIM*(OTT<ConnType,numPol>::ind2C(conn_elem[connface[0]]));
+ const double* BB=coords+SPACEDIM*(OTT<ConnType,numPol>::ind2C(conn_elem[connface[1]]));
+ const double* CC=coords+SPACEDIM*(OTT<ConnType,numPol>::ind2C(conn_elem[connface[2]]));
+
+ double Vol=triple_product(AA,BB,CC,ptToTest);
+ if (Vol<-eps)
+ sign[iface]=-1;
+ else if (Vol>eps)
+ sign[iface]=1;
+ else
+ sign[iface]=0;
+ }
+ bool ret=decideFromSign(sign, nbfaces);
+ delete [] sign;
+ return ret;
}
static bool isElementContainsPoint(const double *ptToTest, NormalizedCellType type, const double *coords, const int *conn_elem)
const double* a=coords+SPACEDIM*(OTT<ConnType,numPol>::ind2C(conn_elem[iedge]));
std::copy(a,a+SPACEDIM,pts+iedge*SPACEDIM);
}
- bool ret=isElementContainsPointAlg2D(ptToTest,pts,nbEdges);
+ bool ret=isElementContainsPointAlg2D(ptToTest,pts,nbEdges,1e-12);
delete [] pts;
return ret;
}
if (SPACEDIM==3)
{
- int nbfaces = cmType.getNumberOfSons();
- int* sign = new int[nbfaces];
- for (int iface=0; iface<nbfaces; iface++)
- {
- const unsigned* connface=cmType.getNodesConstituentTheSon(iface);
- const double* AA=coords+SPACEDIM*(OTT<ConnType,numPol>::ind2C(conn_elem[connface[0]]));
- const double* BB=coords+SPACEDIM*(OTT<ConnType,numPol>::ind2C(conn_elem[connface[1]]));
- const double* CC=coords+SPACEDIM*(OTT<ConnType,numPol>::ind2C(conn_elem[connface[2]]));
-
- double Vol=triple_product(AA,BB,CC,ptToTest);
- if (Vol<-1e-12)
- sign[iface]=-1;
- else if (Vol>1e-12)
- sign[iface]=1;
- else
- sign[iface]=0;
- }
- bool ret=decideFromSign(sign, nbfaces);
- delete [] sign;
- return ret;
+ return isElementContainsPointAlg3D(ptToTest,conn_elem,coords,cmType,1e-12);
}
}
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POLYHEDRONINTERSECTORP0P0_HXX__
+#define __POLYHEDRONINTERSECTORP0P0_HXX__
+
+#include "Intersector3DP0P0.hxx"
+#include "SplitterTetra.hxx"
+#include "NormalizedUnstructuredMesh.hxx"
+
+namespace INTERP_KERNEL
+{
+
+
+ /**
+ * \brief Class responsible for calculating intersection between a hexahedron target element and
+ * the source elements.
+ *
+ */
+ template<class MyMeshType, class MyMatrix>
+ class PolyhedronIntersectorP0P0 : public Intersector3DP0P0<MyMeshType,MyMatrix>
+ {
+ public:
+ static const int SPACEDIM=MyMeshType::MY_SPACEDIM;
+ static const int MESHDIM=MyMeshType::MY_MESHDIM;
+ typedef typename MyMeshType::MyConnType ConnType;
+ static const NumberingPolicy numPol=MyMeshType::My_numPol;
+ public:
+
+ PolyhedronIntersectorP0P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, SplittingPolicy policy = GENERAL_24);
+
+ ~PolyhedronIntersectorP0P0();
+
+ void intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res);
+
+ private:
+ void releaseArrays();
+ private:
+ /// pointers to the SplitterTetra objects representing the tetrahedra
+ /// that result from the splitting of the hexahedron target cell
+ std::vector< SplitterTetra<MyMeshType>* > _tetra;
+
+ SplitterTetra2<MyMeshType> _split;
+
+ };
+}
+
+#endif
--- /dev/null
+// Copyright (C) 2007-2008 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 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
+//
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+//
+#ifndef __POLYHEDRONINTERSECTORP0P0_TXX__
+#define __POLYHEDRONINTERSECTORP0P0_TXX__
+
+#include "PolyhedronIntersectorP0P0.hxx"
+#include "Intersector3DP0P0.txx"
+#include "MeshUtils.hxx"
+
+#include "SplitterTetra.txx"
+
+namespace INTERP_KERNEL
+{
+
+ /**
+ * Constructor creating object from target cell global number
+ * The constructor first calculates the necessary nodes,
+ * (depending on the splitting policy) and then splits the hexahedron into
+ * tetrahedra, placing these in the internal vector _tetra.
+ *
+ * @param targetMesh mesh containing the target elements
+ * @param srcMesh mesh containing the source elements
+ * @param policy splitting policy to be used
+ */
+ template<class MyMeshType, class MyMatrix>
+ PolyhedronIntersectorP0P0<MyMeshType,MyMatrix>::PolyhedronIntersectorP0P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, SplittingPolicy policy):Intersector3DP0P0<MyMeshType,MyMatrix>(targetMesh,srcMesh),_split(targetMesh,srcMesh,policy)
+ {
+ }
+
+ /**
+ * Destructor.
+ * Liberates the SplitterTetra objects and potential sub-node points that have been allocated.
+ *
+ */
+ template<class MyMeshType, class MyMatrix>
+ PolyhedronIntersectorP0P0<MyMeshType,MyMatrix>::~PolyhedronIntersectorP0P0()
+ {
+ releaseArrays();
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ void PolyhedronIntersectorP0P0<MyMeshType,MyMatrix>::releaseArrays()
+ {
+ for(typename std::vector< SplitterTetra<MyMeshType>* >::iterator iter = _tetra.begin(); iter != _tetra.end(); ++iter)
+ delete *iter;
+ _split.releaseArrays();
+ _tetra.clear();
+ }
+
+ /**
+ * Calculates the volume of intersection of an element in the source mesh and the target element
+ * represented by the object.
+ * The calculation is performed by calling the corresponding method for
+ * each SplitterTetra object created by the splitting.
+ *
+ * @param targetCell in C mode.
+ * @param srcCells in C mode.
+ *
+ */
+ template<class MyMeshType, class MyMatrix>
+ void PolyhedronIntersectorP0P0<MyMeshType,MyMatrix>::intersectCells(ConnType targetCell, const std::vector<ConnType>& srcCells, MyMatrix& res)
+ {
+ int nbOfNodesT=Intersector3D<MyMeshType,MyMatrix>::_target_mesh.getNumberOfNodesOfElement(OTT<ConnType,numPol>::indFC(targetCell));
+ releaseArrays();
+ _split.splitTargetCell(targetCell,nbOfNodesT,_tetra);
+ for(typename std::vector<ConnType>::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++)
+ {
+ double volume = 0.;
+ for(typename std::vector<SplitterTetra<MyMeshType>*>::iterator iter = _tetra.begin(); iter != _tetra.end(); ++iter)
+ volume += (*iter)->intersectSourceCell(*iterCellS);
+ if(volume!=0.)
+ res[targetCell].insert(std::make_pair(OTT<ConnType,numPol>::indFC(*iterCellS), volume));
+ }
+ _split.releaseArrays();
+ }
+}
+
+#endif
const unsigned long numSrcElems = srcMesh.getNumberOfElements();
for(unsigned long i = 0 ; i < numSrcElems ; ++i)
- if ( srcMesh.getTypeOfElement( i ) != NORM_TETRA4 )
+ if ( srcMesh.getTypeOfElement( OTT<ConnType,numPol>::indFC(i) ) != NORM_TETRA4 )
throw INTERP_KERNEL::Exception("P1P0 barycentric algorithm works only with tetrahedral source meshes");
}
