1 // Copyright (C) 2007-2024 CEA, EDF
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.
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.
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
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 // Author : Anthony Geay (EDF R&D)
23 #include "Interpolation3D1D.hxx"
24 #include "Interpolation.txx"
25 #include "MeshElement.txx"
26 #include "PointLocator3DIntersectorP0P0.txx"
27 #include "PointLocator3DIntersectorP0P1.txx"
28 #include "PointLocator3DIntersectorP1P0.txx"
29 #include "PointLocator3DIntersectorP1P1.txx"
32 #include "InterpolationHelper.txx"
34 namespace INTERP_KERNEL
37 * Very similar to Interpolation3D::interpolateMeshes, except for the bounding boxes that can be
38 * adjusted in a similar fashion as in InterpolationPlanar::performAdjustmentOfBB()
40 template<class MyMeshType, class MatrixType>
41 typename MyMeshType::MyConnType Interpolation3D1D::interpolateMeshes(const MyMeshType& srcMesh, const MyMeshType& targetMesh, MatrixType& result, const std::string& method)
43 if(InterpolationOptions::getIntersectionType() != PointLocator)
44 INTERP_KERNEL::Exception("Invalid 3D/1D-0D intersection type specified : must be PointLocator.");
46 typedef typename MyMeshType::MyConnType ConnType;
47 // create MeshElement objects corresponding to each element of the two meshes
48 const ConnType numTargetElems = targetMesh.getNumberOfElements();
50 LOG(2, "Target mesh has " << numTargetElems << " elements ");
52 std::unique_ptr< Intersector3D<MyMeshType,MatrixType> > intersector;
53 std::string methC = InterpolationOptions::filterInterpolationMethod(method);
55 { intersector.reset( new PointLocator3DIntersectorP0P0<MyMeshType,MatrixType>(targetMesh, srcMesh, getPrecision()) );
57 else if(methC=="P0P1")
58 { intersector.reset( new PointLocator3DIntersectorP0P1<MyMeshType,MatrixType>(targetMesh, srcMesh, getPrecision()) );
60 else if(methC=="P1P0")
61 { intersector.reset( new PointLocator3DIntersectorP1P0<MyMeshType,MatrixType>(targetMesh, srcMesh, getPrecision()) );
63 else if(methC=="P1P1")
64 { intersector.reset( new PointLocator3DIntersectorP1P1<MyMeshType,MatrixType>(targetMesh, srcMesh, getPrecision()) );
67 throw Exception("Invalid method chosen must be in \"P0P0\", \"P0P1\", \"P1P0\" or \"P1P1\".");
68 // create empty maps for all source elements
69 result.resize(intersector->getNumberOfRowsOfResMatrix());
71 BBTreeStandAlone<3,ConnType> tree( BuildBBTreeWithAdjustment(srcMesh,[this](double *bbox, typename MyMeshType::MyConnType sz){ this->adjustBoundingBoxes(bbox,sz); }) );
73 // for each target element, get source elements with which to calculate intersection
74 // - calculate intersection by calling intersectCells
75 for(ConnType i = 0; i < numTargetElems; ++i)
77 MeshElement<ConnType> trgMeshElem(i, targetMesh);
79 const BoundingBox* box = trgMeshElem.getBoundingBox();
80 // get target bbox in right order
82 box->fillInXMinXmaxYminYmaxZminZmaxFormat(targetBox);
84 std::vector<ConnType> intersectElems;
86 tree.getIntersectingElems(targetBox, intersectElems);
88 if ( !intersectElems.empty() )
89 intersector->intersectCells(i,intersectElems,result);
91 return intersector->getNumberOfColsOfResMatrix();