-// Copyright (C) 2007-2014 CEA/DEN, EDF R&D
+// Copyright (C) 2007-2016 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
template<class MyMeshType, class MyMatrix>
PlanarIntersector<MyMeshType,MyMatrix>::PlanarIntersector(const MyMeshType& meshT, const MyMeshType& meshS, double dimCaracteristic, double precision, double md3DSurf, double minDot3DSurf, double medianPlane, bool doRotate, int orientation, int printLevel):
_meshT(meshT),_meshS(meshS),
- _dim_caracteristic(dimCaracteristic),_max_distance_3Dsurf_intersect(md3DSurf),_precision(precision),_median_plane(medianPlane),
+ _dim_caracteristic(dimCaracteristic),_max_distance_3Dsurf_intersect(md3DSurf),_min_dot_btw_3Dsurf_intersect(minDot3DSurf),_precision(precision),_median_plane(medianPlane),
_do_rotate(doRotate),_orientation(orientation),_print_level(printLevel)
{
_connectT=meshT.getConnectivityPtr();
}
/*!
- Computes the bouding box of a given element. iP in numPol mode.
+ Computes the bounding box of a given element. iP in numPol mode.
*/
template<class MyMeshType, class MyMatrix>
void PlanarIntersector<MyMeshType,MyMatrix>::getElemBB(double* bb, const MyMeshType& mesh, ConnType iP, ConnType nb_nodes)
if(i_A2<nb_NodesA && i_B2<nb_NodesB)
{
//Build the normal of the median plane
- same_orientation=dotprod<SPACEDIM>(normal_A,normal_B)>=0;
+
+ double dotProd(dotprod<SPACEDIM>(normal_A,normal_B)/(normA*normB));
+
+ if(fabs(dotProd)<minDot3DSurf)
+ return 0;
+
+ same_orientation=(dotProd>=0);
if(!same_orientation)
for(int idim =0; idim< SPACEDIM; idim++)
Coords_B[SPACEDIM*i_B+idim] -= proj*linear_comb[idim];
}
- //Buid the matrix sending A into the Oxy plane and apply it to A and B
+ //Build the matrix sending A into the Oxy plane and apply it to A and B
if(do_rotate)
{
TranslationRotationMatrix rotation;
//rotate3DTriangle(Coords_A, &Coords_A[SPACEDIM*i_A1], &Coords_A[SPACEDIM*i_A2], rotation);
- rotate3DTriangle(Coords_B, Coords_B+SPACEDIM*i_B1, Coords_B+SPACEDIM*i_B2, rotation);
+ Rotate3DTriangle(Coords_B, Coords_B+SPACEDIM*i_B1, Coords_B+SPACEDIM*i_B2, rotation);
for (int i=0; i<nb_NodesA; i++)
rotation.transform_vector(Coords_A+SPACEDIM*i);
for (int i=0; i<nb_NodesB; i++)
}
template<class MyMeshType, class MyMatrix>
- void PlanarIntersector<MyMeshType,MyMatrix>::rotate3DTriangle(double* PP1, double*PP2, double*PP3,
+ void PlanarIntersector<MyMeshType,MyMatrix>::Rotate3DTriangle(double* PP1, double*PP2, double*PP3,
TranslationRotationMatrix& rotation_matrix)
{
//initializes