src/INTERP_KERNELTest/ThreeDSurfProjectionTest.cxx

   1 // Copyright (C) 2007-2019  CEA/DEN, EDF R&D
   2 //
   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.
   7 //
   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.
  12 //
  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
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19
  20 #include "ThreeDSurfProjectionTest.hxx"
  21 #include "PlanarIntersector.txx"
  22
  23 class MyMeshType
  24 {
  25 public:
  26   static const int MY_SPACEDIM=3;
  27   static const int MY_MESHDIM=3;
  28   static const INTERP_KERNEL::NumberingPolicy My_numPol=INTERP_KERNEL::ALL_C_MODE;
  29   typedef int MyConnType;
  30 };
  31
  32 class MyMatrixType
  33 {
  34 };
  35
  36 void INTERP_TEST::ThreeDSurfProjectionTest::test1()
  37 {
  38   // Two triangles coo and coo2 are perfectly // each others with a distance equal to 1e-6.
  39   // A little rotation to make it more funny.
  40   //coo=DataArrayDouble([0.,0.,0.,1.,0.,0.,0.,1.,0.],3,3)
  41   //eps=1e-6
  42   //coo2=DataArrayDouble([0.,0.,eps,1.,0.,eps,0.,1.,eps],3,3)
  43   //MEDCouplingPointSet.Rotate3DAlg([0.,0.,0.],[2.,1.,3.],0.3,coo)
  44   //MEDCouplingPointSet.Rotate3DAlg([0.,0.,0.],[2.,1.,3.],0.3,coo2)
  45   const double coo[9]={0.,0.,0.,0.96809749223257568,0.24332379388106262,-0.059839592782071335,-0.23056279077409292,0.95852673990234838,0.16753294721527912};
  46   const double coo2[9]={9.8122602102980502e-08,-1.4839144255482456e-7,9.8404874611628791e-7,0.96809759035517784,0.24332364548962007,-0.059838608733325221,-0.23056269265149082,0.9585265915109058,0.16753393126402524};
  47   double *tmp0(new double[9]),*tmp1(new double[9]);
  48   int ret;
  49   //eps=1e-2. eps is a tolerance to detect that two points are the same or not in a same polygon.
  50   // here the max 3D distance is 1e-5 > 1e-6 so 1 is expected
  51   std::copy(coo,coo+9,tmp0);
  52   std::copy(coo2,coo2+9,tmp1);
  53   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-5/* <- */,-1.,0.5,true);
  54   CPPUNIT_ASSERT_EQUAL(1,ret);
  55   const double expected0[9]={0.,0.,0.,1.,0.,0.,0.,1.,0.};
  56   for(int i=0;i<9;i++)
  57     {
  58       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected0[i],tmp0[i],1e-15);
  59       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected0[i],tmp1[i],1e-15);
  60     }
  61   // here the max 3D distance is 1e-8 < 1e-6 so 0 is expected
  62   std::copy(coo,coo+9,tmp0);
  63   std::copy(coo2,coo2+9,tmp1);
  64   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-8/* <- */,-1.,0.5,true);
  65   CPPUNIT_ASSERT_EQUAL(0,ret);
  66   // here testing when max 3D distance is 1e-5 > 1e-6 with inverted cells
  67   std::copy(coo,coo+9,tmp0);
  68   std::copy(coo2,coo2+3,tmp1+6); std::copy(coo2+3,coo2+6,tmp1+3); std::copy(coo2+6,coo2+9,tmp1);
  69   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-5/* <- */,-1.,0.5,true);
  70   CPPUNIT_ASSERT_EQUAL(-1,ret);
  71   const double expected1[9]={-0.7071067811865476,-0.7071067811865476,0.,0.,-1.4142135623730951,0.,-1.4142135623730951,-1.4142135623730951,0.};
  72   const double expected2[9]={-1.4142135623730951,-1.4142135623730951,0.,0.,-1.4142135623730951,0.,-0.7071067811865476,-0.7071067811865476,0.};
  73   for(int i=0;i<9;i++)
  74     {
  75       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected1[i],tmp0[i],1e-14);
  76       CPPUNIT_ASSERT_DOUBLES_EQUAL(expected2[i],tmp1[i],1e-14);
  77     }
  78   //
  79   delete [] tmp0;
  80   delete [] tmp1;
  81 }
  82
  83 void INTERP_TEST::ThreeDSurfProjectionTest::test2()
  84 {// here the two triangles have their center of inertia very close (eps) but the angle between the two planes is "big"
  85   //coo=DataArrayDouble([0.,0.,0.,1.,0.,0.,0.,1.,0.],3,3)
  86   //coocpy=coo.deepCopy()
  87   //MEDCouplingPointSet.Rotate3DAlg([0.,0.,0.],[-1,-1.,0.],pi/3,coocpy)
  88   //coocpy+=[eps*sqrt(3)/2,eps/2,eps*0.]
  89   //
  90   const double coo[9]={0.,0.,0.,0.96809749223257568,0.24332379388106262,-0.059839592782071335,-0.23056279077409292,0.95852673990234838,0.16753294721527912};
  91   const double coo2[9]={7.2311562622637225e-07,6.8998795679738294e-07,3.1943866106249849e-08,0.72852072144314628,0.33125439126063028,0.5996079016637561,0.0090154262465889021,0.87059752249869415,-0.49191448334281612};
  92   double *tmp0(new double[9]),*tmp1(new double[9]);
  93   int ret;
  94   //eps=1e-2. eps is a tolerance to detect that two points are the same or not in a same polygon.
  95   // here the max 3D distance is 1e-5 > 1e-6 so 1 is expected
  96   std::copy(coo,coo+9,tmp0);
  97   std::copy(coo2,coo2+9,tmp1);
  98   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-5/* <- */,-1.,0.5,true);
  99   CPPUNIT_ASSERT_EQUAL(1,ret);
 100   // here the max 3D distance is 1e-8 < 1e-6 so 0 is expected
 101   std::copy(coo,coo+9,tmp0);
 102   std::copy(coo2,coo2+9,tmp1);
 103   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-8/* <- */,-1.,0.5,true);
 104   CPPUNIT_ASSERT_EQUAL(0,ret);
 105   // again max 3D distance is 1e-5 > 1e-6 so 1 is expected
 106   std::copy(coo,coo+9,tmp0);
 107   std::copy(coo2,coo2+9,tmp1);
 108   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-5/* <- */,-1.,0.5,true);
 109   CPPUNIT_ASSERT_EQUAL(1,ret);
 110   // again max 3D distance is 1e-5 > 1e-6 but minDot set to 0.8. 0 expected. because the angle is pi/4 so cos(pi/3) > 0.8
 111   std::copy(coo,coo+9,tmp0);
 112   std::copy(coo2,coo2+9,tmp1);
 113   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-5/* <- */,0.8/* <- */,0.5,true);
 114   CPPUNIT_ASSERT_EQUAL(0,ret);
 115   // again max 3D distance is 1e-5 > 1e-6 but minDot set to 0.7. 1 expected. because the angle is pi/4 so cos(pi/3) < 0.49
 116   std::copy(coo,coo+9,tmp0);
 117   std::copy(coo2,coo2+9,tmp1);
 118   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-5/* <- */,0.49/* <- */,0.5,true);
 119   CPPUNIT_ASSERT_EQUAL(1,ret);
 120   // again max 3D distance is 1e-5 > 1e-6 but minDot set to 0.7. 0 expected. because the angle is pi/4 so cos(pi/3) > 0.51
 121   std::copy(coo,coo+9,tmp0);
 122   std::copy(coo2,coo2+9,tmp1);
 123   ret=INTERP_KERNEL::PlanarIntersector<MyMeshType,MyMatrixType>::Projection(tmp0,tmp1,3,3,1e-2,1e-5/* <- */,0.51/* <- */,0.5,true);
 124   CPPUNIT_ASSERT_EQUAL(0,ret);
 125   //
 126   delete [] tmp0;
 127   delete [] tmp1;
 128 }