--- /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 __CURVEINTERSECTOR_HXX__
+#define __CURVEINTERSECTOR_HXX__
+
+#include "TargetIntersector.hxx"
+#include "NormalizedUnstructuredMesh.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class CurveIntersector : public TargetIntersector<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:
+ CurveIntersector(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double adjustmentEpsAbs, double medianLine, int printLevel);
+ virtual ~CurveIntersector();
+ void createBoundingBoxes(const MyMeshType& mesh, std::vector<double>& bbox);
+ void adjustBoundingBoxes(std::vector<double>& bbox, double adjustmentEpsAbs);
+ static void getElemBB(double* bb, const MyMeshType& mesh, ConnType iP, ConnType nb_nodes);
+
+ protected :
+ bool getRealTargetCoordinates(ConnType icellT, std::vector<double>& coordsT);
+ bool getRealSourceCoordinates(ConnType icellS, std::vector<double>& coordsS);
+ double intersectSegments(double *Coords_T, double *Coords_S);
+
+ struct TDualSegment
+ {
+ std::vector<double> _coords;
+ int _nodeId; // in mesh mode
+ };
+ static void getDualSegments(ConnType icell,
+ const MyMeshType& mesh,
+ std::vector<TDualSegment>& segments);
+
+ protected:
+ const ConnType *_connectT;
+ const ConnType *_connectS;
+ const double *_coordsT;
+ const double *_coordsS;
+ const ConnType *_connIndexT;
+ const ConnType *_connIndexS;
+ const MyMeshType& _meshT;
+ const MyMeshType& _meshS;
+ double _tolerance;
+ double _precision;
+ double _median_line;
+ int _print_level;
+ };
+}
+
+#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 __CURVEINTERSECTOR_TXX__
+#define __CURVEINTERSECTOR_TXX__
+
+#include "CurveIntersector.hxx"
+
+#include <limits>
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ CurveIntersector<MyMeshType,MyMatrix>
+ ::CurveIntersector(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance, double medianLine, int printLevel):
+ _meshT(meshT),
+ _meshS(meshS),
+ _tolerance(tolerance),
+ _precision(precision),
+ _median_line(medianLine),
+ _print_level(printLevel)
+ {
+ if ( SPACEDIM != 1 && SPACEDIM != 2 )
+ throw Exception("CurveIntersector(): space dimension of mesh must be 1 or 2");
+ if ( MESHDIM != 1 )
+ throw Exception("CurveIntersector(): mesh dimension must be 1");
+
+ _connectT = meshT.getConnectivityPtr();
+ _connectS = meshS.getConnectivityPtr();
+ _connIndexT = meshT.getConnectivityIndexPtr();
+ _connIndexS = meshS.getConnectivityIndexPtr();
+ _coordsT = meshT.getCoordinatesPtr();
+ _coordsS = meshS.getCoordinatesPtr();
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ CurveIntersector<MyMeshType,MyMatrix>::~CurveIntersector()
+ {
+ }
+
+ //================================================================================
+ /*!
+ \brief creates the bounding boxes for all the cells of mesh \a mesh
+
+ \param mesh structure pointing to the mesh
+ \param bbox vector containing the bounding boxes
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ void CurveIntersector<MyMeshType,MyMatrix>::createBoundingBoxes (const MyMeshType& mesh,
+ std::vector<double>& bbox)
+ {
+ long nbelems = mesh.getNumberOfElements();
+ bbox.resize(2*SPACEDIM* nbelems);
+ const double* coords = mesh.getCoordinatesPtr();
+ const ConnType* conn = mesh.getConnectivityPtr();
+ const ConnType* conn_index = mesh.getConnectivityIndexPtr();
+ int ibox=0;
+ for(long icell=0; icell<nbelems; icell++)
+ {
+ int nb_nodes_per_elem = conn_index[icell+1]-conn_index[icell];
+ //initializing bounding box limits
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ bbox[2*SPACEDIM*ibox+2*idim] = std::numeric_limits<double>::max();
+ bbox[2*SPACEDIM*ibox+2*idim+1] = -std::numeric_limits<double>::max();
+ }
+ //updating the bounding box with each node of the element
+ for (int j=0; j<nb_nodes_per_elem; j++)
+ {
+ const double* coord_node = coords +
+ SPACEDIM*OTT<ConnType,numPol>
+ ::coo2C(conn[OTT<ConnType,numPol>::conn2C(conn_index[icell]+j)]);
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ double x = *(coord_node+idim);
+ bbox[ibox*2*SPACEDIM + 2*idim] =
+ ( bbox[ibox*2*SPACEDIM + 2*idim] < x ) ? bbox[ibox*2*SPACEDIM + 2*idim] : x;
+ bbox[ibox*2*SPACEDIM + 2*idim+1] =
+ ( bbox[ibox*2*SPACEDIM + 2*idim+1] > x ) ? bbox[ibox*2*SPACEDIM + 2*idim+1] : x;
+ }
+ }
+ ibox++;
+ }
+ }
+
+ //================================================================================
+ /*!
+ Computes the bouding box of a given element. iP in numPol mode.
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ void CurveIntersector<MyMeshType,MyMatrix>::getElemBB (double* bb,
+ const MyMeshType& mesh,
+ ConnType iP,
+ ConnType nb_nodes)
+ {
+ const double* coords = mesh.getCoordinatesPtr();
+ const ConnType* conn_index = mesh.getConnectivityIndexPtr();
+ const ConnType* conn = mesh.getConnectivityPtr();
+ //initializing bounding box limits
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ bb[2*idim ] = std::numeric_limits<double>::max();
+ bb[2*idim+1] = -std::numeric_limits<double>::max();
+ }
+
+ for (ConnType i=0; i<nb_nodes; i++)
+ {
+ //MN: iP= cell index, not node index, use of connectivity array ?
+ const double* coord_node = coords +
+ SPACEDIM*(OTT<ConnType,numPol>::coo2C(conn[OTT<ConnType,numPol>::conn2C(conn_index[OTT<ConnType,numPol>::ind2C(iP)]+i)]));
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ double x = *(coord_node+idim);
+ bb[2*idim ] = (x<bb[2*idim ]) ? x : bb[2*idim ];
+ bb[2*idim+1] = (x>bb[2*idim+1]) ? x : bb[2*idim+1];
+ }
+ }
+ }
+
+ //================================================================================
+ /*! Readjusts a set of bounding boxes so that they are extended
+ in all dimensions for avoiding missing interesting intersections
+
+ \param bbox vector containing the bounding boxes
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ void CurveIntersector<MyMeshType,MyMatrix>::adjustBoundingBoxes (std::vector<double>& bbox,
+ double adjustmentEpsAbs)
+ {
+ long size = bbox.size()/(2*SPACEDIM);
+ for (int i=0; i<size; i++)
+ {
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ bbox[i*2*SPACEDIM+2*idim ] -= adjustmentEpsAbs;
+ bbox[i*2*SPACEDIM+2*idim+1] += adjustmentEpsAbs;
+ }
+ }
+ }
+
+ //================================================================================
+ /*!
+ * @param icellT id in target mesh in format of MyMeshType.
+ * @param coordsT output val that stores coordinates of the target cell
+ * automatically resized to the right length.
+ * @return true if segment is quadratic and in this case coordinates of medium node
+ * are placed in the middle of coordsT
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ bool CurveIntersector<MyMeshType,MyMatrix>::getRealTargetCoordinates
+ (ConnType icellT, std::vector<double>& coordsT)
+ {
+ int nbNodesT = _connIndexT[OTT<ConnType,numPol>::ind2C(icellT)+1] -
+ _connIndexT[OTT<ConnType,numPol>::ind2C(icellT)];
+ coordsT.resize(SPACEDIM*nbNodesT);
+ for (ConnType iT=0; iT<nbNodesT; iT++)
+ {
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ coordsT[SPACEDIM*iT+idim] =
+ _coordsT[SPACEDIM*OTT<ConnType,numPol>::coo2C(_connectT[OTT<ConnType,numPol>::conn2C(_connIndexT[OTT<ConnType,numPol>::ind2C(icellT)]+iT)])+idim];
+ }
+ }
+ if ( nbNodesT > 2 )
+ {
+ for(int idim=0; idim<SPACEDIM; idim++)
+ std::swap( coordsT[SPACEDIM*1+idim], coordsT[SPACEDIM*2+idim]);
+ return true;
+ }
+ return false;
+ }
+
+ //================================================================================
+ /*!
+ * @param icellS id in source mesh in format of MyMeshType.
+ * @param coordsS output val that stores coordinates of the source cell automatically resized to the right length.
+ * @return true if segment is quadratic and in this case coordinates of medium node
+ * are placed in the middle of coordsS
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ bool CurveIntersector<MyMeshType,MyMatrix>::getRealSourceCoordinates
+ (ConnType icellS, std::vector<double>& coordsS)
+ {
+ int nbNodesS = _connIndexS[OTT<ConnType,numPol>::ind2C(icellS)+1] -
+ _connIndexS[OTT<ConnType,numPol>::ind2C(icellS)];
+ coordsS.resize(SPACEDIM*nbNodesS);
+ for(ConnType iS=0; iS<nbNodesS; iS++)
+ {
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ coordsS[SPACEDIM*iS+idim] =
+ _coordsS[SPACEDIM*OTT<ConnType,numPol>::coo2C(_connectS[OTT<ConnType,numPol>::conn2C(_connIndexS[OTT<ConnType,numPol>::ind2C(icellS)]+iS)])+idim];
+ }
+ }
+ if ( nbNodesS > 2 )
+ {
+ for(int idim=0; idim<SPACEDIM; idim++)
+ std::swap( coordsS[SPACEDIM*1+idim], coordsS[SPACEDIM*2+idim]);
+ return true;
+ }
+ return false;
+ }
+
+ //================================================================================
+ /*!
+ * \brief Return dual segments of given segment
+ * \param icell - given segment in C mode
+ * \param mesh - mesh
+ * \param segments - dual segments
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ void CurveIntersector<MyMeshType,MyMatrix>::getDualSegments(ConnType icell,
+ const MyMeshType& mesh,
+ std::vector<TDualSegment>& segments)
+ {
+ // get coordinates of cell nodes
+ int nbNodes;
+ std::vector<double> ncoords;
+ std::vector<int> nodeIds;
+ {
+ const ConnType *connect = mesh.getConnectivityPtr();
+ const ConnType *connIndex = mesh.getConnectivityIndexPtr();
+ const double *coords = mesh.getCoordinatesPtr();
+
+ nbNodes = connIndex[icell+1] - connIndex[icell];
+
+ ncoords.resize(SPACEDIM*nbNodes);
+ nodeIds.resize(nbNodes);
+
+ for(int i=0; i<nbNodes; i++)
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ nodeIds[i] = connect[OTT<ConnType,numPol>::conn2C(connIndex[OTT<ConnType,numPol>::ind2C(icell)]+i)];
+ ncoords[SPACEDIM*i+idim] = coords[SPACEDIM*OTT<ConnType,numPol>::coo2C(nodeIds[i])+idim];
+ }
+ if ( nbNodes > 2 ) // quadratic segment, put medium node in the middle
+ {
+ for(int idim=0; idim<SPACEDIM; idim++)
+ std::swap( ncoords[SPACEDIM*1+idim], ncoords[SPACEDIM*2+idim]);
+ std::swap( nodeIds[1], nodeIds[2] );
+ }
+ }
+
+ // fill segments
+ segments.clear();
+ segments.reserve( 2*nbNodes );
+ for(int i=0; i<nbNodes-1; i++)
+ {
+ segments.push_back(TDualSegment());
+ TDualSegment& seg1 = segments.back();
+ segments.push_back(TDualSegment());
+ TDualSegment& seg2 = segments.back();
+
+ seg1._nodeId = nodeIds[i];
+ seg2._nodeId = nodeIds[i+1];
+
+ seg1._coords.resize( SPACEDIM * 2 );
+ seg2._coords.resize( SPACEDIM * 2 );
+
+ for(int idim=0; idim<SPACEDIM; idim++)
+ {
+ double c1 = ncoords[SPACEDIM*i+idim];
+ double c2 = ncoords[SPACEDIM*(i+1)+idim];
+ double m = 0.5 * ( c1 + c2 );
+ seg1._coords[ idim ] = c1;
+ seg1._coords[ SPACEDIM + idim ] = m;
+ seg2._coords[ idim ] = m;
+ seg2._coords[ SPACEDIM + idim ] = c2;
+ }
+ }
+ }
+
+ //================================================================================
+ /*!
+ * \brief Return length of intersection of two segments
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ double CurveIntersector<MyMeshType,MyMatrix>::intersectSegments(double *Coords_T,
+ double *Coords_S)
+ {
+ double xt0 = Coords_T[0], xt1 = Coords_T[1];
+ double xs0 = Coords_S[0], xs1 = Coords_S[1];
+ if ( SPACEDIM == 2 )
+ {
+ // Pass 2D->1D
+
+ enum { X=0, Y };
+
+ // check if two segments overlap in 2D within tolerance
+
+ double* t0 = Coords_T;
+ double* t1 = Coords_T + 2;
+ double t01[2] = { t1[X]-t0[X], t1[Y]-t0[Y] }; // tgt segment direction
+ double tSize = sqrt( t01[X]*t01[X] + t01[Y]*t01[Y] ); // tgt segment size
+ if ( tSize < _precision ) return 0; // degenerated segment
+ t01[X] /= tSize, t01[Y] /= tSize; // normalize t01
+
+ double* s0 = Coords_S;
+ double* s1 = Coords_S + 2;
+ double t0s0[2] = { s0[X]-t0[X], s0[Y]-t0[Y] };
+ double t0s1[2] = { s1[X]-t0[X], s1[Y]-t0[Y] };
+ double nt01_x_t0s0 = t0s0[X] * t01[Y] - t0s0[Y] * t01[X]; // t0s0 dot norm of t01
+ double nt01_x_t0s1 = t0s1[X] * t01[Y] - t0s1[Y] * t01[X]; // t0s1 dot norm of t01
+ double dist_ts0 = fabs( nt01_x_t0s0 ); // dist from tgt seg to s0
+ double dist_ts1 = fabs( nt01_x_t0s1 ); // dist from tgt seg to s1
+ bool s0_out_of_tol = ( dist_ts0 > _tolerance );
+ bool s1_out_of_tol = ( dist_ts1 > _tolerance );
+ if ( nt01_x_t0s0 * nt01_x_t0s1 > 0 && ( s0_out_of_tol || s1_out_of_tol ))
+ return 0; // tgt segment is to far from src segment
+
+ double S0[2] = { s0[X], s0[Y] };
+ double S1[2] = { s1[X], s1[Y] };
+ if ( s0_out_of_tol ) // put s0 within tolerance
+ {
+ double t = _tolerance * nt01_x_t0s0 / dist_ts0; // signed tolerance
+ double r = ( nt01_x_t0s0 - t ) / ( nt01_x_t0s0 - nt01_x_t0s1 );
+ S0[X] = s0[X] * ( 1.-r ) + s1[X] * r;
+ S0[Y] = s0[Y] * ( 1.-r ) + s1[Y] * r;
+ }
+ if ( s1_out_of_tol ) // put s1 within tolerance
+ {
+ double t = _tolerance * nt01_x_t0s1 / dist_ts1; // signed tolerance
+ double r = ( nt01_x_t0s1 - t ) / ( nt01_x_t0s1 - nt01_x_t0s0 );
+ S1[X] = s1[X] * ( 1.-r ) + s0[X] * r;
+ S1[Y] = s1[Y] * ( 1.-r ) + s0[Y] * r;
+ }
+
+ // project tgt and src segments to median line
+
+ double s01[2] = { S1[X]-S0[X], S1[Y]-S0[Y] }; // src segment direction
+ double sSize = sqrt( s01[X]*s01[X] + s01[Y]*s01[Y] ); // src segment size
+ if ( sSize < _precision ) return 0; // degenerated segment
+ s01[X] /= sSize, s01[Y] /= sSize; // normalize s01
+
+ // make t01 and s01 codirected
+ double t01_x_s01 = t01[X] * s01[X] + t01[Y] * s01[Y]; // t01 dot s01
+ if ( t01_x_s01 < 0 )
+ s01[X] = -s01[X], s01[Y] = -s01[Y];
+
+ double medianDir[2] = {
+ t01[X] * ( 1.-_median_line) + s01[X] * _median_line,
+ t01[Y] * ( 1.-_median_line) + s01[Y] * _median_line
+ };
+ double medianSize = sqrt( medianDir[X]*medianDir[X] + medianDir[Y]*medianDir[Y] );
+ if ( medianSize < std::numeric_limits<double>::min() )
+ return 0; // strange...
+ medianDir[X] /= medianSize, medianDir[Y] /= medianSize;
+
+ xt0 = t0[X] * medianDir[X] + t0[Y] * medianDir[Y];
+ xt1 = t1[X] * medianDir[X] + t1[Y] * medianDir[Y];
+ xs0 = S0[X] * medianDir[X] + S0[Y] * medianDir[Y];
+ xs1 = S1[X] * medianDir[X] + S1[Y] * medianDir[Y];
+
+ } // if ( SPACEDIM == 2 )
+
+ if ( xt0 > xt1 ) std::swap( xt0, xt1 );
+ if ( xs0 > xs1 ) std::swap( xs0, xs1 );
+
+ double x0 = std::max( xt0, xs0 );
+ double x1 = std::min( xt1, xs1 );
+ return ( x0 < x1 ) ? ( x1 - x0 ) : 0.;
+ }
+
+}
+
+#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 __CURVEINTERSECTORP0P0_HXX__
+#define __CURVEINTERSECTORP0P0_HXX__
+
+#include "CurveIntersector.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class CurveIntersectorP0P0 : public CurveIntersector<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;
+
+ CurveIntersectorP0P0(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance,
+ double medianLine, int printLevel);
+ public:
+ int getNumberOfRowsOfResMatrix() const;
+ int getNumberOfColsOfResMatrix() const;
+ void intersectCells(ConnType icellT,
+ const std::vector<ConnType>& icellsS, MyMatrix& res);
+ };
+}
+
+#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 __CURVEINTERSECTORP0P0_TXX__
+#define __CURVEINTERSECTORP0P0_TXX__
+
+#include "CurveIntersectorP0P0.hxx"
+#include "CurveIntersector.txx"
+
+#define BASE_INTERSECTOR CurveIntersector<MyMeshType,MyMatrix>
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ CurveIntersectorP0P0<MyMeshType,MyMatrix>
+ ::CurveIntersectorP0P0(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance,
+ double medianLine, int printLevel):
+ BASE_INTERSECTOR(meshT, meshS, precision, tolerance, medianLine, printLevel)
+ {
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ int CurveIntersectorP0P0<MyMeshType,MyMatrix>
+ ::getNumberOfRowsOfResMatrix() const
+ {
+ return BASE_INTERSECTOR::_meshT.getNumberOfElements();
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ int CurveIntersectorP0P0<MyMeshType,MyMatrix>
+ ::getNumberOfColsOfResMatrix() const
+ {
+ return BASE_INTERSECTOR::_meshS.getNumberOfElements();
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ void CurveIntersectorP0P0<MyMeshType,MyMatrix>
+ ::intersectCells(ConnType icellT, const std::vector<ConnType>& icellsS, MyMatrix& res)
+ {
+ typename MyMatrix::value_type& resRow = res[icellT];
+ std::vector<double> coordsT;
+ int t, nbSegT = 1 + BASE_INTERSECTOR::getRealTargetCoordinates(icellT,coordsT);
+ for ( t = 0; t < nbSegT; ++t )
+ for(typename std::vector<ConnType>::const_iterator
+ iter=icellsS.begin(); iter!=icellsS.end(); iter++)
+ {
+ int iS = *iter;
+ std::vector<double> coordsS;
+ int s, nbSegS = 1 + BASE_INTERSECTOR::getRealSourceCoordinates(iS,coordsS);
+ for ( s = 0; s < nbSegS; ++s )
+ {
+ double surf = BASE_INTERSECTOR::intersectSegments(coordsT.data() + t*SPACEDIM,
+ coordsS.data() + s*SPACEDIM);
+ if(surf!=0.)
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(iS),surf));
+ }
+ }
+ }
+}
+#undef BASE_INTERSECTOR
+
+#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 __CURVEINTERSECTORP0P1_HXX__
+#define __CURVEINTERSECTORP0P1_HXX__
+
+#include "CurveIntersector.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class CurveIntersectorP0P1 : public CurveIntersector<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;
+
+ CurveIntersectorP0P1(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance,
+ double medianLine, int printLevel);
+ public:
+ void intersectCells(ConnType icellT,
+ const std::vector<ConnType>& icellsS, MyMatrix& res);
+ int getNumberOfRowsOfResMatrix() const;
+ int getNumberOfColsOfResMatrix() const;
+ };
+}
+
+#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 __CurveIntersectorP0P1_TXX__
+#define __CurveIntersectorP0P1_TXX__
+
+#include "CurveIntersectorP0P1.hxx"
+#include "CurveIntersector.txx"
+
+#define BASE_INTERSECTOR CurveIntersector<MyMeshType,MyMatrix>
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ CurveIntersectorP0P1<MyMeshType,MyMatrix>
+ ::CurveIntersectorP0P1(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance,
+ double medianLine, int printLevel):
+ BASE_INTERSECTOR(meshT, meshS, precision, tolerance, medianLine, printLevel)
+ {
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ int CurveIntersectorP0P1<MyMeshType,MyMatrix>
+ ::getNumberOfRowsOfResMatrix() const
+ {
+ return BASE_INTERSECTOR::_meshT.getNumberOfNodes();
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ int CurveIntersectorP0P1<MyMeshType,MyMatrix>
+ ::getNumberOfColsOfResMatrix() const
+ {
+ return BASE_INTERSECTOR::_meshS.getNumberOfElements();
+ }
+
+ //================================================================================
+ /*!
+ * \brief Project from segments to nodes
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ void CurveIntersectorP0P1<MyMeshType,MyMatrix>
+ ::intersectCells(ConnType icellT, const std::vector<ConnType>& icellsS, MyMatrix& res)
+ {
+ std::vector<typename BASE_INTERSECTOR::TDualSegment> segmentsT;
+ BASE_INTERSECTOR::getDualSegments( icellT, BASE_INTERSECTOR::_meshT, segmentsT);
+ for ( int t = 0; t < segmentsT.size(); ++t )
+ {
+ typename MyMatrix::value_type& resRow = res[ OTT<ConnType,numPol>::ind2C( segmentsT[t]._nodeId )];
+ for(typename std::vector<ConnType>::const_iterator
+ iter=icellsS.begin(); iter!=icellsS.end(); iter++)
+ {
+ int iS = *iter;
+ std::vector<double> coordsS;
+ int s, nbSegS = 1 + BASE_INTERSECTOR::getRealSourceCoordinates(iS,coordsS);
+ for ( s = 0; s < nbSegS; ++s )
+ {
+ double surf = BASE_INTERSECTOR::intersectSegments(segmentsT[t]._coords.data(),
+ coordsS.data() + s*SPACEDIM);
+ if(surf!=0.)
+ {
+ typename MyMatrix::value_type::const_iterator iterRes=resRow.find(OTT<ConnType,numPol>::indFC(iS));
+ if(iterRes==resRow.end())
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(iS),surf));
+ else
+ {
+ surf+=(*iterRes).second;
+ resRow.erase(OTT<ConnType,numPol>::indFC(iS));
+ resRow.insert(std::make_pair(OTT<ConnType,numPol>::indFC(iS),surf));
+ }
+ }
+ }
+ }
+ }
+ }
+}
+#undef BASE_INTERSECTOR
+
+#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 __CURVEINTERSECTORP1P0_HXX__
+#define __CURVEINTERSECTORP1P0_HXX__
+
+#include "CurveIntersector.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class CurveIntersectorP1P0 : public CurveIntersector<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;
+
+ CurveIntersectorP1P0(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance, double medianLine, int printLevel);
+ public:
+ void intersectCells(ConnType icellT, const std::vector<ConnType>& icellsS, MyMatrix& res);
+ int getNumberOfRowsOfResMatrix() const;
+ int getNumberOfColsOfResMatrix() const;
+ };
+}
+
+#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 __CurveIntersectorP1P0_TXX__
+#define __CurveIntersectorP1P0_TXX__
+
+#include "CurveIntersectorP1P0.hxx"
+#include "CurveIntersector.txx"
+
+#define BASE_INTERSECTOR CurveIntersector<MyMeshType,MyMatrix>
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ CurveIntersectorP1P0<MyMeshType,MyMatrix>
+ ::CurveIntersectorP1P0(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance,
+ double medianLine, int printLevel):
+ BASE_INTERSECTOR (meshT, meshS, precision, tolerance, medianLine, printLevel)
+ {
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ int CurveIntersectorP1P0<MyMeshType,MyMatrix>
+ ::getNumberOfRowsOfResMatrix() const
+ {
+ return BASE_INTERSECTOR::_meshT.getNumberOfElements();
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ int CurveIntersectorP1P0<MyMeshType,MyMatrix>
+ ::getNumberOfColsOfResMatrix() const
+ {
+ return BASE_INTERSECTOR::_meshS.getNumberOfNodes();
+ }
+
+ //================================================================================
+ /*!
+ * \brief Project from source nodes to target segments
+ */
+ //================================================================================
+
+ template<class MyMeshType, class MyMatrix>
+ void CurveIntersectorP1P0<MyMeshType,MyMatrix>
+ ::intersectCells(ConnType icellT, const std::vector<ConnType>& icellsS, MyMatrix& res)
+ {
+ typename MyMatrix::value_type& resRow = res[ icellT ];
+ std::vector<typename BASE_INTERSECTOR::TDualSegment> segmentsS;
+ std::vector<double> coordsT;
+ int t, nbSegT = 1 + BASE_INTERSECTOR::getRealTargetCoordinates(icellT,coordsT);
+ for ( t = 0; t < nbSegT; ++t )
+ for(typename std::vector<ConnType>::const_iterator
+ iter=icellsS.begin(); iter!=icellsS.end(); iter++)
+ {
+ int iS = *iter;
+ BASE_INTERSECTOR::getDualSegments( OTT<ConnType,numPol>::ind2C(iS),
+ BASE_INTERSECTOR::_meshS, segmentsS);
+ for ( int s = 0; s < segmentsS.size(); ++s )
+ {
+ double surf = BASE_INTERSECTOR::intersectSegments(segmentsS[s]._coords.data(),
+ coordsT.data() + t*SPACEDIM);
+ if(surf!=0.)
+ {
+ int nS = segmentsS[s]._nodeId;
+ typename MyMatrix::value_type::const_iterator iterRes=resRow.find(nS);
+ if(iterRes==resRow.end())
+ resRow.insert(std::make_pair(nS,surf));
+ else
+ {
+ surf+=(*iterRes).second;
+ resRow.erase(nS);
+ resRow.insert(std::make_pair(nS,surf));
+ }
+ }
+ }
+ }
+ }
+}
+#undef BASE_INTERSECTOR
+
+#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 __CURVEINTERSECTORP1P1_HXX__
+#define __CURVEINTERSECTORP1P1_HXX__
+
+#include "CurveIntersector.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ class CurveIntersectorP1P1 : public CurveIntersector<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;
+
+ CurveIntersectorP1P1(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance,
+ double medianLine, int printLevel);
+ public:
+ void intersectCells(ConnType icellT,
+ const std::vector<ConnType>& icellsS, MyMatrix& res);
+ int getNumberOfRowsOfResMatrix() const;
+ int getNumberOfColsOfResMatrix() const;
+ };
+}
+
+#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 __CurveIntersectorP1P1_TXX__
+#define __CurveIntersectorP1P1_TXX__
+
+#include "CurveIntersectorP1P1.hxx"
+#include "CurveIntersector.txx"
+
+#define BASE_INTERSECTOR CurveIntersector<MyMeshType,MyMatrix>
+
+namespace INTERP_KERNEL
+{
+ template<class MyMeshType, class MyMatrix>
+ CurveIntersectorP1P1<MyMeshType,MyMatrix>
+ ::CurveIntersectorP1P1(const MyMeshType& meshT, const MyMeshType& meshS,
+ double precision, double tolerance,
+ double medianLine, int printLevel):
+ BASE_INTERSECTOR (meshT, meshS, precision, tolerance, medianLine, printLevel)
+ {
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ int CurveIntersectorP1P1<MyMeshType,MyMatrix>
+ ::getNumberOfRowsOfResMatrix() const
+ {
+ return BASE_INTERSECTOR::_meshT.getNumberOfNodes();
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ int CurveIntersectorP1P1<MyMeshType,MyMatrix>
+ ::getNumberOfColsOfResMatrix() const
+ {
+ return BASE_INTERSECTOR::_meshS.getNumberOfNodes();
+ }
+
+ template<class MyMeshType, class MyMatrix>
+ void CurveIntersectorP1P1<MyMeshType,MyMatrix>
+ ::intersectCells(ConnType icellT, const std::vector<ConnType>& icellsS, MyMatrix& res)
+ {
+ std::vector<typename BASE_INTERSECTOR::TDualSegment> segmentsT, segmentsS;
+ BASE_INTERSECTOR::getDualSegments( icellT, BASE_INTERSECTOR::_meshT, segmentsT);
+ for ( int t = 0; t < segmentsT.size(); ++t )
+ {
+ typename MyMatrix::value_type& resRow = res[ OTT<ConnType,numPol>::ind2C( segmentsT[t]._nodeId )];
+ for(typename std::vector<ConnType>::const_iterator
+ iter=icellsS.begin(); iter!=icellsS.end(); iter++)
+ {
+ int iS = *iter;
+ BASE_INTERSECTOR::getDualSegments( OTT<ConnType,numPol>::ind2C(iS),
+ BASE_INTERSECTOR::_meshS, segmentsS);
+ for ( int s = 0; s < segmentsS.size(); ++s )
+ {
+ double surf = BASE_INTERSECTOR::intersectSegments(segmentsT[t]._coords.data(),
+ segmentsS[s]._coords.data());
+ if(surf!=0.)
+ {
+ int nS = segmentsS[s]._nodeId;
+ typename MyMatrix::value_type::const_iterator iterRes=resRow.find(nS);
+ if(iterRes==resRow.end())
+ resRow.insert(std::make_pair(nS,surf));
+ else
+ {
+ surf+=(*iterRes).second;
+ resRow.erase(nS);
+ resRow.insert(std::make_pair(nS,surf));
+ }
+ }
+ }
+ }
+ }
+ }
+}
+#undef BASE_INTERSECTOR
+
+#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 __INTERPOLATION1D_HXX__
+#define __INTERPOLATION1D_HXX__
+
+#include "InterpolationCurve.hxx"
+
+namespace INTERP_KERNEL
+{
+ class Interpolation1D : public InterpolationCurve<Interpolation1D>
+ {
+ public:
+ Interpolation1D() { }
+ Interpolation1D(const InterpolationOptions& io):InterpolationCurve<Interpolation1D>(io) {}
+ };
+}
+
+#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 __INTERPOLATION1D_TXX__
+#define __INTERPOLATION1D_TXX__
+
+#include "Interpolation1D.hxx"
+
+#include "InterpolationCurve.txx"
+
+#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 __INTERPOLATION2DCURVE_HXX__
+#define __INTERPOLATION2DCURVE_HXX__
+
+#include "InterpolationCurve.hxx"
+#include "InterpolationOptions.hxx"
+
+namespace INTERP_KERNEL
+{
+ class Interpolation2DCurve : public InterpolationCurve<Interpolation2DCurve>
+ {
+ public:
+ Interpolation2DCurve();
+ Interpolation2DCurve(const InterpolationOptions& io);
+ // geometric precision, intersection tolerance, coice of the median line,
+ void setOptions(double precision, double tolerance, double medianLine);
+ };
+}
+
+#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 __INTERPOLATION2DCURVE_TXX__
+#define __INTERPOLATION2DCURVE_TXX__
+
+#include "Interpolation2DCurve.hxx"
+#include "InterpolationCurve.txx"
+
+namespace INTERP_KERNEL
+{
+ Interpolation2DCurve::Interpolation2DCurve()
+ {
+ // to have non-zero default thickness of target element
+ InterpolationOptions::setBoundingBoxAdjustmentAbs( InterpolationOptions::getPrecision() );
+ }
+
+ Interpolation2DCurve::Interpolation2DCurve
+ (const InterpolationOptions& io):InterpolationCurve<Interpolation2DCurve>(io)
+ {
+ // to have non-zero default thickness of target element
+ InterpolationOptions::setBoundingBoxAdjustmentAbs( InterpolationOptions::getPrecision() );
+ }
+
+ /**
+ * \brief Function used to set the options for the intersection calculation
+ * \details The following options can be modified:
+ * -# Precision: Level of precision of the computations.
+ * - Values: positive real number.
+ * - Default: 1.0E-12.
+ * -# Tolerance: Thickness of target element.
+ * - Values: positive real number.
+ * - Default: 1.0E-12.
+ * -# Median line: Position of the median line where both segments will be projected.
+ * - Values: real number between 0.0 and 1.0.
+ * - Default: 0.5
+ */
+ void Interpolation2DCurve::setOptions (double precision,
+ double tolerance,
+ double medianLine)
+ {
+ InterpolationOptions::setPrecision(precision);
+ InterpolationOptions::setBoundingBoxAdjustmentAbs(tolerance);
+ InterpolationOptions::setMedianPlane(medianLine);
+ }
+}
+
+#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 __INTERPOLATIONCURVE_HXX__
+#define __INTERPOLATIONCURVE_HXX__
+
+#include "Interpolation.hxx"
+#include "InterpolationOptions.hxx"
+
+namespace INTERP_KERNEL
+{
+ template<class RealCurve>
+ class InterpolationCurve : public Interpolation< InterpolationCurve<RealCurve> >
+ {
+ public:
+ InterpolationCurve();
+ InterpolationCurve(const InterpolationOptions & io);
+
+ // Main function to interpolate
+ template<class MyMeshType, class MatrixType>
+ int interpolateMeshes(const MyMeshType& mesh1, const MyMeshType& mesh2,
+ MatrixType& result, const char *method);
+
+ };
+}
+
+#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 __INTERPOLATIONCURVE_TXX__
+#define __INTERPOLATIONCURVE_TXX__
+
+#include "InterpolationCurve.hxx"
+#include "InterpolationOptions.hxx"
+#include "CurveIntersectorP0P0.txx"
+#include "CurveIntersectorP1P0.txx"
+#include "CurveIntersectorP0P1.txx"
+#include "CurveIntersectorP1P1.txx"
+#include "BBTree.txx"
+
+#include <time.h>
+
+namespace INTERP_KERNEL
+{
+ /**
+ * \defgroup interpolationCurve InterpolationCurve
+ *
+ * \class InterpolationCurve
+ * \brief Class used to compute the coefficients of the interpolation matrix between
+ * two local meshes in two dimensions.
+ */
+ template<class RealCurve>
+ InterpolationCurve<RealCurve>::InterpolationCurve()
+ {
+ }
+
+ template<class RealCurve>
+ InterpolationCurve<RealCurve>::InterpolationCurve (const InterpolationOptions& io)
+ :Interpolation< InterpolationCurve<RealCurve> >(io)
+ {
+ }
+
+ /** \brief Main function to interpolate 1D meshes.
+ \details The algorithm proceeds in two steps: first a filtering process reduces the number of pairs of elements for which the
+ * calculation must be carried out by eliminating pairs that do not intersect based on their bounding boxes. Then, the
+ * volume of intersection is calculated by an object of type IntersectorPlanar for the remaining pairs, and entered into the
+ * intersection matrix.
+ *
+ * The matrix is partially sparse : it is a vector of maps of integer - double pairs.
+ * The length of the vector is equal to the number of target elements - for each target element there is a map, regardless
+ * of whether the element intersects any source elements or not. But in the maps there are only entries for those source elements
+ * which have a non-zero intersection volume with the target element. The vector has indices running from
+ * 0 to (#target elements - 1), meaning that the map for target element i is stored at index i - 1. In the maps, however,
+ * the indexing is more natural : the intersection volume of the target element i with source element j is found at matrix[i-1][j].
+ *
+
+ * @param myMeshS Planar source mesh
+ * @Param myMeshT Planar target mesh
+ * @return vector containing for each element i of the source mesh, a map giving for each element j
+ * of the target mesh which i intersects, the area of the intersection
+ *
+ */
+ template<class RealCurve>
+ template<class MyMeshType, class MatrixType>
+ int InterpolationCurve<RealCurve>::interpolateMeshes (const MyMeshType& myMeshS,
+ const MyMeshType& myMeshT,
+ MatrixType& result,
+ const char * method)
+ {
+ static const int SPACEDIM=MyMeshType::MY_SPACEDIM;
+ typedef typename MyMeshType::MyConnType ConnType;
+ static const NumberingPolicy numPol = MyMeshType::My_numPol;
+
+ long global_start = clock();
+ int counter=0;
+
+ long nbMailleS = myMeshS.getNumberOfElements();
+ long nbMailleT = myMeshT.getNumberOfElements();
+
+ CurveIntersector<MyMeshType,MatrixType>* intersector=0;
+ std::string meth(method);
+ if(meth=="P0P0") {
+ intersector = new CurveIntersectorP0P0<MyMeshType,MatrixType>
+ (myMeshT, myMeshS,
+ InterpolationOptions::getPrecision(),
+ InterpolationOptions::getBoundingBoxAdjustmentAbs(),
+ InterpolationOptions::getMedianPlane(),
+ InterpolationOptions::getPrintLevel());
+ }
+ else if(meth=="P0P1") {
+ intersector = new CurveIntersectorP0P1<MyMeshType,MatrixType>
+ (myMeshT, myMeshS,
+ InterpolationOptions::getPrecision(),
+ InterpolationOptions::getBoundingBoxAdjustmentAbs(),
+ InterpolationOptions::getMedianPlane(),
+ InterpolationOptions::getPrintLevel());
+ }
+ else if(meth=="P1P0") {
+ intersector = new CurveIntersectorP1P0<MyMeshType,MatrixType>
+ (myMeshT, myMeshS,
+ InterpolationOptions::getPrecision(),
+ InterpolationOptions::getBoundingBoxAdjustmentAbs(),
+ InterpolationOptions::getMedianPlane(),
+ InterpolationOptions::getPrintLevel());
+ }
+ else if(meth=="P1P1") {
+ intersector = new CurveIntersectorP1P1<MyMeshType,MatrixType>
+ (myMeshT, myMeshS,
+ InterpolationOptions::getPrecision(),
+ InterpolationOptions::getBoundingBoxAdjustmentAbs(),
+ InterpolationOptions::getMedianPlane(),
+ InterpolationOptions::getPrintLevel());
+ }
+ else
+ throw INTERP_KERNEL::Exception("Invalid method specified ! Must be in : \"P0P0\" \"P0P1\" \"P1P0\" or \"P1P1\"");
+ /****************************************************************/
+ /* Create a search tree based on the bounding boxes */
+ /* Instanciate the intersector and initialise the result vector */
+ /****************************************************************/
+
+ long start_filtering=clock();
+
+ std::vector<double> bbox;
+ intersector->createBoundingBoxes(myMeshS,bbox); // create the bounding boxes
+ intersector->adjustBoundingBoxes(bbox, InterpolationOptions::getBoundingBoxAdjustmentAbs());
+ BBTree<SPACEDIM,ConnType> my_tree(&bbox[0], 0, 0, nbMailleS);//creating the search structure
+
+ long end_filtering = clock();
+
+ result.resize(intersector->getNumberOfRowsOfResMatrix());//on initialise.
+
+ /****************************************************/
+ /* Loop on the target cells - core of the algorithm */
+ /****************************************************/
+ long start_intersection = clock();
+ const ConnType *connIndxT = myMeshT.getConnectivityIndexPtr();
+ for(int iT=0; iT<nbMailleT; iT++)
+ {
+ int nb_nodesT = connIndxT[iT+1] - connIndxT[iT];
+ std::vector<int> intersecting_elems;
+ double bb[2*SPACEDIM];
+ intersector->getElemBB(bb,myMeshT,OTT<ConnType,numPol>::indFC(iT),nb_nodesT);
+ my_tree.getIntersectingElems(bb, intersecting_elems);
+ intersector->intersectCells(iT,intersecting_elems,result);
+ counter += intersecting_elems.size();
+ }
+ int ret = intersector->getNumberOfColsOfResMatrix();
+ delete intersector;
+
+ if (InterpolationOptions::getPrintLevel() >= 1) {
+ long end_intersection=clock();
+ std::cout << "Filtering time= " << end_filtering-start_filtering << std::endl;
+ std::cout << "Intersection time= " << end_intersection-start_intersection << std::endl;
+ long global_end =clock();
+ std::cout << "Number of computed intersections = " << counter << std::endl;
+ std::cout << "Global time= " << global_end - global_start << std::endl;
+ }
+ return ret;
+ }
+}
+
+#endif
VectorUtils.hxx \
VolSurfFormulae.hxx \
VolSurfUser.hxx \
-VolSurfUser.txx
+VolSurfUser.txx \
+CurveIntersector.hxx \
+CurveIntersector.txx \
+CurveIntersectorP0P0.hxx \
+CurveIntersectorP0P0.txx \
+CurveIntersectorP0P1.hxx \
+CurveIntersectorP1P0.hxx \
+CurveIntersectorP1P1.hxx \
+Interpolation1D.hxx \
+Interpolation1D.txx \
+Interpolation2DCurve.hxx \
+Interpolation2DCurve.txx \
+InterpolationCurve.hxx \
+InterpolationCurve.txx
EXTRA_DIST += \
