From a143ec1de773b6f39c5aa0b7f556a21261b6f519 Mon Sep 17 00:00:00 2001 From: ageay Date: Fri, 23 Oct 2009 05:44:42 +0000 Subject: [PATCH] *** empty log message *** --- .../PointLocator3DIntersectorP0P0.hxx | 44 ++++++++ .../PointLocator3DIntersectorP0P0.txx | 77 +++++++++++++ .../PointLocator3DIntersectorP0P1.hxx | 44 ++++++++ .../PointLocator3DIntersectorP0P1.txx | 80 ++++++++++++++ .../PointLocator3DIntersectorP1P0.hxx | 44 ++++++++ .../PointLocator3DIntersectorP1P0.txx | 98 +++++++++++++++++ .../PointLocator3DIntersectorP1P1.hxx | 44 ++++++++ .../PointLocator3DIntersectorP1P1.txx | 101 ++++++++++++++++++ src/INTERP_KERNEL/PointLocatorAlgos.txx | 76 ++++++------- .../PolyhedronIntersectorP0P0.hxx | 63 +++++++++++ .../PolyhedronIntersectorP0P0.txx | 94 ++++++++++++++++ .../PolyhedronIntersectorP1P0Bary.txx | 2 +- 12 files changed, 730 insertions(+), 37 deletions(-) create mode 100644 src/INTERP_KERNEL/PointLocator3DIntersectorP0P0.hxx create mode 100644 src/INTERP_KERNEL/PointLocator3DIntersectorP0P0.txx create mode 100644 src/INTERP_KERNEL/PointLocator3DIntersectorP0P1.hxx create mode 100644 src/INTERP_KERNEL/PointLocator3DIntersectorP0P1.txx create mode 100644 src/INTERP_KERNEL/PointLocator3DIntersectorP1P0.hxx create mode 100644 src/INTERP_KERNEL/PointLocator3DIntersectorP1P0.txx create mode 100644 src/INTERP_KERNEL/PointLocator3DIntersectorP1P1.hxx create mode 100644 src/INTERP_KERNEL/PointLocator3DIntersectorP1P1.txx create mode 100644 src/INTERP_KERNEL/PolyhedronIntersectorP0P0.hxx create mode 100644 src/INTERP_KERNEL/PolyhedronIntersectorP0P0.txx diff --git a/src/INTERP_KERNEL/PointLocator3DIntersectorP0P0.hxx b/src/INTERP_KERNEL/PointLocator3DIntersectorP0P0.hxx new file mode 100644 index 000000000..f7fb73961 --- /dev/null +++ b/src/INTERP_KERNEL/PointLocator3DIntersectorP0P0.hxx @@ -0,0 +1,44 @@ +// 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 PointLocator3DIntersectorP0P0 : public Intersector3DP0P0 + { + 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& srcCells, MyMatrix& res); + protected: + double _precision; + }; +} + +#endif diff --git a/src/INTERP_KERNEL/PointLocator3DIntersectorP0P0.txx b/src/INTERP_KERNEL/PointLocator3DIntersectorP0P0.txx new file mode 100644 index 000000000..5de6bae9c --- /dev/null +++ b/src/INTERP_KERNEL/PointLocator3DIntersectorP0P0.txx @@ -0,0 +1,77 @@ +// 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 + PointLocator3DIntersectorP0P0::PointLocator3DIntersectorP0P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision): + Intersector3DP0P0(targetMesh,srcMesh),_precision(precision) + { + } + + template + PointLocator3DIntersectorP0P0::~PointLocator3DIntersectorP0P0() + { + } + + /** + * + * @param targetCell in C mode. + * @param srcCells in C mode. + * + */ + template + void PointLocator3DIntersectorP0P0::intersectCells(ConnType targetCell, const std::vector& srcCells, MyMatrix& res) + { + std::vector CoordsT; + Intersector3DP0P0::getRealTargetCoordinates(OTT::indFC(targetCell),CoordsT); + double bary[SPACEDIM]; + calculateBarycenterDyn2(&CoordsT[0],CoordsT.size()/SPACEDIM,bary); + typename MyMatrix::value_type& resRow=res[targetCell]; + const double *coordsS=Intersector3DP0P0::_src_mesh.getCoordinatesPtr(); + for(typename std::vector::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++) + { + NormalizedCellType tS=Intersector3DP0P0::_src_mesh.getTypeOfElement(OTT::indFC(*iterCellS)); + const CellModel& cmTypeS=CellModel::getCellModel(tS); + std::vector connOfCurCellS; + Intersector3DP0P0::getConnOfSourceCell(OTT::indFC(*iterCellS),connOfCurCellS); + if(PointLocatorAlgos::isElementContainsPointAlg3D(bary,&connOfCurCellS[0],coordsS,cmTypeS,_precision)) + { + resRow.insert(std::make_pair(OTT::indFC(*iterCellS),1)); + } + } + } +} + +#endif diff --git a/src/INTERP_KERNEL/PointLocator3DIntersectorP0P1.hxx b/src/INTERP_KERNEL/PointLocator3DIntersectorP0P1.hxx new file mode 100644 index 000000000..22d1abefd --- /dev/null +++ b/src/INTERP_KERNEL/PointLocator3DIntersectorP0P1.hxx @@ -0,0 +1,44 @@ +// 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 PointLocator3DIntersectorP0P1 : public Intersector3DP0P1 + { + 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& srcCells, MyMatrix& res); + protected: + double _precision; + }; +} + +#endif diff --git a/src/INTERP_KERNEL/PointLocator3DIntersectorP0P1.txx b/src/INTERP_KERNEL/PointLocator3DIntersectorP0P1.txx new file mode 100644 index 000000000..c59f8a962 --- /dev/null +++ b/src/INTERP_KERNEL/PointLocator3DIntersectorP0P1.txx @@ -0,0 +1,80 @@ +// 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 + PointLocator3DIntersectorP0P1::PointLocator3DIntersectorP0P1(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision): + Intersector3DP0P1(targetMesh,srcMesh),_precision(precision) + { + } + + template + PointLocator3DIntersectorP0P1::~PointLocator3DIntersectorP0P1() + { + } + + /** + * + * @param targetCell in C mode. + * @param srcCells in C mode. + * + */ + template + void PointLocator3DIntersectorP0P1::intersectCells(ConnType targetCell, const std::vector& srcCells, MyMatrix& res) + { + std::vector coordsTarget; + Intersector3DP0P1::getRealTargetCoordinates(OTT::indFC(targetCell),coordsTarget); + int nbNodesT=coordsTarget.size()/SPACEDIM; + const double *coordsS=Intersector3DP0P1::_src_mesh.getCoordinatesPtr(); + const ConnType *startOfCellNodeConnT=Intersector3DP0P1::getStartConnOfTargetCell(targetCell); + for(typename std::vector::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++) + { + NormalizedCellType tS=Intersector3DP0P1::_src_mesh.getTypeOfElement(OTT::indFC(*iterCellS)); + const CellModel& cmTypeS=CellModel::getCellModel(tS); + std::vector connOfCurCellS; + Intersector3DP0P1::getConnOfSourceCell(OTT::indFC(*iterCellS),connOfCurCellS); + for(int nodeIdT=0;nodeIdT::isElementContainsPointAlg3D(&coordsTarget[nodeIdT*SPACEDIM],&connOfCurCellS[0],coordsS,cmTypeS,_precision)) + { + ConnType curNodeTInCmode=OTT::coo2C(startOfCellNodeConnT[nodeIdT]); + typename MyMatrix::value_type& resRow=res[curNodeTInCmode]; + typename MyMatrix::value_type::const_iterator iterRes=resRow.find(OTT::indFC(*iterCellS)); + if(iterRes==resRow.end()) + resRow.insert(std::make_pair(OTT::indFC(*iterCellS),1.)); + } + } + } + } +} + +#endif diff --git a/src/INTERP_KERNEL/PointLocator3DIntersectorP1P0.hxx b/src/INTERP_KERNEL/PointLocator3DIntersectorP1P0.hxx new file mode 100644 index 000000000..921aaf8af --- /dev/null +++ b/src/INTERP_KERNEL/PointLocator3DIntersectorP1P0.hxx @@ -0,0 +1,44 @@ +// 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 PointLocator3DIntersectorP1P0 : public Intersector3DP1P0 + { + 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& srcCells, MyMatrix& res); + protected: + double _precision; + }; +} + +#endif diff --git a/src/INTERP_KERNEL/PointLocator3DIntersectorP1P0.txx b/src/INTERP_KERNEL/PointLocator3DIntersectorP1P0.txx new file mode 100644 index 000000000..e666f3a29 --- /dev/null +++ b/src/INTERP_KERNEL/PointLocator3DIntersectorP1P0.txx @@ -0,0 +1,98 @@ +// 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 + PointLocator3DIntersectorP1P0::PointLocator3DIntersectorP1P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision):Intersector3DP1P0(targetMesh,srcMesh),_precision(precision) + { + } + + template + PointLocator3DIntersectorP1P0::~PointLocator3DIntersectorP1P0() + { + } + + /** + * @param targetCell in C mode. + * @param srcCells in C mode. + * + * WARNING : for all methods on _split object source and target are switched ! + */ + template + void PointLocator3DIntersectorP1P0::intersectCells(ConnType targetCell, const std::vector& srcCells, MyMatrix& res) + { + std::vector CoordsT; + typename MyMatrix::value_type& resRow=res[targetCell]; + const double *coordsS=Intersector3DP1P0::_src_mesh.getCoordinatesPtr(); + Intersector3DP1P0::getRealTargetCoordinates(OTT::indFC(targetCell),CoordsT); + double baryT[SPACEDIM]; + calculateBarycenterDyn2(&CoordsT[0],CoordsT.size()/SPACEDIM,baryT); + for(typename std::vector::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++) + { + NormalizedCellType tS=Intersector3DP1P0::_src_mesh.getTypeOfElement(OTT::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 connOfCurCellS; + Intersector3DP1P0::getConnOfSourceCell(OTT::indFC(*iterCellS),connOfCurCellS); + if( PointLocatorAlgos::isElementContainsPointAlg3D(baryT,&connOfCurCellS[0],coordsS,cmTypeS,_precision) ) + { + double resLoc[4]; + std::vector srcCell; + Intersector3DP1P0::getRealSourceCoordinates(OTT::indFC(*iterCellS),srcCell); + std::vector 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::getStartConnOfSourceCell(*iterCellS); + for(int nodeIdS=0;nodeIdS<4;nodeIdS++) + { + if(fabs(resLoc[nodeIdS])>_precision) + { + ConnType curNodeSInCmode=OTT::coo2C(startOfCellNodeConn[nodeIdS]); + typename MyMatrix::value_type::const_iterator iterRes=resRow.find(OTT::indFC(curNodeSInCmode)); + if(iterRes==resRow.end()) + resRow.insert(std::make_pair(OTT::indFC(curNodeSInCmode),resLoc[nodeIdS])); + else + { + double val=(*iterRes).second+resLoc[nodeIdS]; + resRow.erase(OTT::indFC(curNodeSInCmode)); + resRow.insert(std::make_pair(OTT::indFC(curNodeSInCmode),val)); + } + } + } + } + } + } +} + +#endif diff --git a/src/INTERP_KERNEL/PointLocator3DIntersectorP1P1.hxx b/src/INTERP_KERNEL/PointLocator3DIntersectorP1P1.hxx new file mode 100644 index 000000000..b0444d137 --- /dev/null +++ b/src/INTERP_KERNEL/PointLocator3DIntersectorP1P1.hxx @@ -0,0 +1,44 @@ +// 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 PointLocator3DIntersectorP1P1 : public Intersector3DP1P1 + { + 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& srcCells, MyMatrix& res); + protected: + double _precision; + }; +} + +#endif diff --git a/src/INTERP_KERNEL/PointLocator3DIntersectorP1P1.txx b/src/INTERP_KERNEL/PointLocator3DIntersectorP1P1.txx new file mode 100644 index 000000000..8b29fc250 --- /dev/null +++ b/src/INTERP_KERNEL/PointLocator3DIntersectorP1P1.txx @@ -0,0 +1,101 @@ +// 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 + PointLocator3DIntersectorP1P1::PointLocator3DIntersectorP1P1(const MyMeshType& targetMesh, const MyMeshType& srcMesh, double precision): + Intersector3DP1P1(targetMesh,srcMesh),_precision(precision) + { + } + + template + PointLocator3DIntersectorP1P1::~PointLocator3DIntersectorP1P1() + { + } + + /** + * @param targetCell in C mode. + * @param srcCells in C mode. + */ + template + void PointLocator3DIntersectorP1P1::intersectCells(ConnType targetCell, const std::vector& srcCells, MyMatrix& res) + { + std::vector CoordsT; + Intersector3DP1P1::getRealTargetCoordinates(OTT::indFC(targetCell),CoordsT); + int nbOfNodesT=CoordsT.size()/SPACEDIM; + const double *coordsS=Intersector3DP1P1::_src_mesh.getCoordinatesPtr(); + for(typename std::vector::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++) + { + NormalizedCellType tS=Intersector3DP1P1::_src_mesh.getTypeOfElement(OTT::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::getStartConnOfTargetCell(targetCell); + for(int nodeIdT=0;nodeIdT::ind2C(startOfCellNodeConnT[nodeIdT])]; + std::vector connOfCurCellS; + Intersector3DP1P1::getConnOfSourceCell(OTT::indFC(*iterCellS),connOfCurCellS); + if( PointLocatorAlgos::isElementContainsPointAlg3D(&CoordsT[nodeIdT*SPACEDIM],&connOfCurCellS[0],coordsS,cmTypeS,_precision) ) + { + double resLoc[4]; + std::vector localCoordsS; + Intersector3DP1P1::getRealSourceCoordinates(OTT::indFC(*iterCellS),localCoordsS); + std::vector 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::getStartConnOfSourceCell(*iterCellS); + for(int nodeIdS=0;nodeIdS<4;nodeIdS++) + { + if(fabs(resLoc[nodeIdS])>_precision) + { + ConnType curNodeSInCmode=OTT::coo2C(startOfCellNodeConnS[nodeIdS]); + typename MyMatrix::value_type::const_iterator iterRes=resRow.find(OTT::indFC(curNodeSInCmode)); + if(iterRes==resRow.end()) + resRow.insert(std::make_pair(OTT::indFC(curNodeSInCmode),resLoc[nodeIdS])); + else + { + double val=(*iterRes).second+resLoc[nodeIdS]; + resRow.erase(OTT::indFC(curNodeSInCmode)); + resRow.insert(std::make_pair(OTT::indFC(curNodeSInCmode),val)); + } + } + } + } + } + } + } +} + +#endif diff --git a/src/INTERP_KERNEL/PointLocatorAlgos.txx b/src/INTERP_KERNEL/PointLocatorAlgos.txx index c4b4cdf35..1db1d3c7f 100644 --- a/src/INTERP_KERNEL/PointLocatorAlgos.txx +++ b/src/INTERP_KERNEL/PointLocatorAlgos.txx @@ -100,18 +100,18 @@ namespace INTERP_KERNEL 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; iedge1e-12) + else if(a>eps) sign[iedge]=1; else sign[iedge]=0; @@ -131,9 +131,32 @@ namespace INTERP_KERNEL 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::ind2C(conn_elem[connface[0]])); + const double* BB=coords+SPACEDIM*(OTT::ind2C(conn_elem[connface[1]])); + const double* CC=coords+SPACEDIM*(OTT::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) @@ -153,33 +176,14 @@ namespace INTERP_KERNEL const double* a=coords+SPACEDIM*(OTT::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::ind2C(conn_elem[connface[0]])); - const double* BB=coords+SPACEDIM*(OTT::ind2C(conn_elem[connface[1]])); - const double* CC=coords+SPACEDIM*(OTT::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); } } diff --git a/src/INTERP_KERNEL/PolyhedronIntersectorP0P0.hxx b/src/INTERP_KERNEL/PolyhedronIntersectorP0P0.hxx new file mode 100644 index 000000000..3d49f4df5 --- /dev/null +++ b/src/INTERP_KERNEL/PolyhedronIntersectorP0P0.hxx @@ -0,0 +1,63 @@ +// 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 PolyhedronIntersectorP0P0 : public Intersector3DP0P0 + { + 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& 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* > _tetra; + + SplitterTetra2 _split; + + }; +} + +#endif diff --git a/src/INTERP_KERNEL/PolyhedronIntersectorP0P0.txx b/src/INTERP_KERNEL/PolyhedronIntersectorP0P0.txx new file mode 100644 index 000000000..1e3b0b36e --- /dev/null +++ b/src/INTERP_KERNEL/PolyhedronIntersectorP0P0.txx @@ -0,0 +1,94 @@ +// 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 + PolyhedronIntersectorP0P0::PolyhedronIntersectorP0P0(const MyMeshType& targetMesh, const MyMeshType& srcMesh, SplittingPolicy policy):Intersector3DP0P0(targetMesh,srcMesh),_split(targetMesh,srcMesh,policy) + { + } + + /** + * Destructor. + * Liberates the SplitterTetra objects and potential sub-node points that have been allocated. + * + */ + template + PolyhedronIntersectorP0P0::~PolyhedronIntersectorP0P0() + { + releaseArrays(); + } + + template + void PolyhedronIntersectorP0P0::releaseArrays() + { + for(typename std::vector< SplitterTetra* >::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 + void PolyhedronIntersectorP0P0::intersectCells(ConnType targetCell, const std::vector& srcCells, MyMatrix& res) + { + int nbOfNodesT=Intersector3D::_target_mesh.getNumberOfNodesOfElement(OTT::indFC(targetCell)); + releaseArrays(); + _split.splitTargetCell(targetCell,nbOfNodesT,_tetra); + for(typename std::vector::const_iterator iterCellS=srcCells.begin();iterCellS!=srcCells.end();iterCellS++) + { + double volume = 0.; + for(typename std::vector*>::iterator iter = _tetra.begin(); iter != _tetra.end(); ++iter) + volume += (*iter)->intersectSourceCell(*iterCellS); + if(volume!=0.) + res[targetCell].insert(std::make_pair(OTT::indFC(*iterCellS), volume)); + } + _split.releaseArrays(); + } +} + +#endif diff --git a/src/INTERP_KERNEL/PolyhedronIntersectorP1P0Bary.txx b/src/INTERP_KERNEL/PolyhedronIntersectorP1P0Bary.txx index 0c8f976c0..15c8f97ce 100644 --- a/src/INTERP_KERNEL/PolyhedronIntersectorP1P0Bary.txx +++ b/src/INTERP_KERNEL/PolyhedronIntersectorP1P0Bary.txx @@ -56,7 +56,7 @@ namespace INTERP_KERNEL const unsigned long numSrcElems = srcMesh.getNumberOfElements(); for(unsigned long i = 0 ; i < numSrcElems ; ++i) - if ( srcMesh.getTypeOfElement( i ) != NORM_TETRA4 ) + if ( srcMesh.getTypeOfElement( OTT::indFC(i) ) != NORM_TETRA4 ) throw INTERP_KERNEL::Exception("P1P0 barycentric algorithm works only with tetrahedral source meshes"); } -- 2.39.2