Back to index

salome-med  6.5.0
UnitTetraIntersectionBaryTest.cxx
Go to the documentation of this file.
00001 // Copyright (C) 2007-2012  CEA/DEN, EDF R&D
00002 //
00003 // This library is free software; you can redistribute it and/or
00004 // modify it under the terms of the GNU Lesser General Public
00005 // License as published by the Free Software Foundation; either
00006 // version 2.1 of the License.
00007 //
00008 // This library is distributed in the hope that it will be useful,
00009 // but WITHOUT ANY WARRANTY; without even the implied warranty of
00010 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00011 // Lesser General Public License for more details.
00012 //
00013 // You should have received a copy of the GNU Lesser General Public
00014 // License along with this library; if not, write to the Free Software
00015 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
00016 //
00017 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
00018 //
00019 
00020 // File      : UnitTetraIntersectionBaryTest.cxx
00021 // Created   : Thu Dec 11 15:54:41 2008
00022 // Author    : Edward AGAPOV (eap)
00023 //
00024 #include "UnitTetraIntersectionBaryTest.hxx"
00025 
00026 #include "UnitTetraIntersectionBary.hxx"
00027 #include "TetraAffineTransform.hxx"
00028 #include "InterpolationUtils.hxx"
00029 #include "SplitterTetra.txx"
00030 
00031 #include <iostream>
00032 
00033 using namespace INTERP_KERNEL;
00034 
00035 namespace INTERP_TEST
00036 {
00037   void fill_UnitTetraIntersectionBary(UnitTetraIntersectionBary& bary, double nodes[][3])
00038   {
00039     int faceConn[4][3] = { { 0, 1, 2 },// inverse order
00040                            { 0, 3, 1 },
00041                            { 1, 3, 2 },
00042                            { 3, 0, 2 } };
00043 //     int faceConn[4][3] = { { 0, 2, 1 },
00044 //                            { 0, 1, 3 },
00045 //                            { 1, 2, 3 },
00046 //                            { 3, 2, 0 } };
00047     bary.init(true);
00048     for ( int i = 0; i < 4; ++i ) {
00049       int* faceNodes = faceConn[ i ];
00050       TransformedTriangle tri(nodes[faceNodes[0]], nodes[faceNodes[1]], nodes[faceNodes[2]]);
00051       tri.calculateIntersectionVolume();
00052       bary.addSide( tri );
00053     }
00054   }
00055   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_1()
00056   {
00057     // cutting tetra coincides with the unit one
00058     double nodes[4][3] = { { 0.0, 0.0, 0.0 },
00059                            { 1.0, 0.0, 0.0 },
00060                            { 0.0, 1.0, 0.0 },
00061                            { 0.0, 0.0, 1.0 } };
00062     UnitTetraIntersectionBary bary;
00063     fill_UnitTetraIntersectionBary(bary,nodes);
00064     double baryCenter[3];
00065     bool ok    = bary.getBary( baryCenter );
00066     double vol = bary.getVolume();
00067     CPPUNIT_ASSERT( ok );
00068     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.166667, vol, 1e-5);
00069     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[0], 1e-5);
00070     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[1], 1e-5);
00071     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[2], 1e-5);
00072   }
00073   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_2()
00074   {
00075     // cutting tetra fully include the unit one
00076     double nodes[4][3] = { {-0.1,-0.1,-0.1 },
00077                            { 1.5,-0.1,-0.1 },
00078                            {-0.1, 1.5,-0.1 },
00079                            {-0.1,-0.1, 1.5 } };
00080     UnitTetraIntersectionBary bary;
00081     fill_UnitTetraIntersectionBary(bary,nodes);
00082     double baryCenter[3];
00083     bool ok    = bary.getBary( baryCenter );
00084     double vol = bary.getVolume();
00085     CPPUNIT_ASSERT( ok );
00086     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.166667, vol, 1e-5);
00087     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[0], 1e-5);
00088     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[1], 1e-5);
00089     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[2], 1e-5);
00090   }
00091   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_3()
00092   {
00093     // cutting tetra is same as the unit one but moved up by 0.5
00094     double nodes[4][3] = { { 0.0, 0.0, 0.5 },
00095                            { 1.0, 0.0, 0.5 },
00096                            { 0.0, 1.0, 0.5 },
00097                            { 0.0, 0.0, 1.5 } };
00098     UnitTetraIntersectionBary bary;
00099     fill_UnitTetraIntersectionBary(bary,nodes);
00100     double baryCenter[3];
00101     bool ok    = bary.getBary( baryCenter );
00102     double vol = bary.getVolume();
00103     CPPUNIT_ASSERT( ok );
00104     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.020833333333333332, vol, 1e-5);
00105     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.125, baryCenter[0], 1e-5);
00106     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.125, baryCenter[1], 1e-5);
00107     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.625, baryCenter[2], 1e-5);
00108   }
00109   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_4()
00110   {
00111     // same as previous but no cutting sides lay on the sides of unit tetra
00112     double nodes[4][3] = { {-0.2,-0.2, 0.5 },
00113                            { 1.0, 0.0, 0.5 },
00114                            { 0.0, 1.0, 0.5 },
00115                            { 0.0, 0.0, 2.0 } };
00116     UnitTetraIntersectionBary bary;
00117     fill_UnitTetraIntersectionBary(bary,nodes);
00118     double baryCenter[3];
00119     bool ok    = bary.getBary( baryCenter );
00120     double vol = bary.getVolume();
00121     CPPUNIT_ASSERT( ok );
00122     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.020833333333333332, vol, 1e-5);
00123     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.125, baryCenter[0], 1e-5);
00124     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.125, baryCenter[1], 1e-5);
00125     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.625, baryCenter[2], 1e-5);
00126   }
00127   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_5()
00128   {
00129     // cutting tetra is similar and parallel to the UT but moved (-0.1,-0.1,-0.1)
00130     double nodes[4][3] = { {-0.1,-0.1,-0.1 },
00131                            { 1.1,-0.1,-0.1 },
00132                            {-0.1, 1.1,-0.1 },
00133                            {-0.1,-0.1, 1.1 } };
00134     UnitTetraIntersectionBary bary;
00135     fill_UnitTetraIntersectionBary(bary,nodes);
00136     double baryCenter[3];
00137     bool ok    = bary.getBary( baryCenter );
00138     double vol = bary.getVolume();
00139     CPPUNIT_ASSERT( ok );
00140     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.1215, vol, 1e-5);
00141     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.225, baryCenter[0], 1e-5);
00142     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.225, baryCenter[1], 1e-5);
00143     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.225, baryCenter[2], 1e-5);
00144   }
00145   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_6()
00146   {
00147     // cutting tetra is deeped into the UT with one corner
00148     double nodes[4][3] = { { 0.2, 0.2, 0.2 },
00149                            { 1.0, 0.2, 0.2 },
00150                            { 0.9, 1.0, 0.2 },
00151                            { 0.9, 9.0, 1.0 } };
00152     UnitTetraIntersectionBary bary;
00153     fill_UnitTetraIntersectionBary(bary,nodes);
00154     double baryCenter[3];
00155     bool ok    = bary.getBary( baryCenter );
00156     double vol = bary.getVolume();
00157     CPPUNIT_ASSERT( ok );
00158     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.000441855, vol, 1e-5);
00159     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.353463 , baryCenter[0], 1e-5 );
00160     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.33877  , baryCenter[1], 1e-5 );
00161     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.207767 , baryCenter[2], 1e-5 );
00162   }
00163   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_7()
00164   {
00165     // cutting tetra passes through the UT with one corner
00166     double nodes[4][3] = { {-0.2, 0.2, 0.2 },
00167                            { 1.0, 0.2, 0.2 },
00168                            { 0.9, 1.0, 0.2 },
00169                            { 0.9, 0.9, 1.0 } };
00170     UnitTetraIntersectionBary bary;
00171     fill_UnitTetraIntersectionBary(bary,nodes);
00172     double baryCenter[3];
00173     bool ok    = bary.getBary( baryCenter );
00174     double vol = bary.getVolume();
00175     CPPUNIT_ASSERT( ok );
00176     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0103501, vol, 1e-5);
00177     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.215578 , baryCenter[0], 1e-5 );
00178     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.341363 , baryCenter[1], 1e-5 );
00179     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.263903 , baryCenter[2], 1e-5 );
00180   }
00181   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_8()
00182   {
00183     // cutting tetra passes through the UT with one edge
00184     double nodes[4][3] = { { 0.5, 0.2, -0.2 }, // O
00185                            {-0.5,-0.2, -0.2 }, // OX
00186                            { 1.0,-0.5, -0.2 }, // OY
00187                            { 0.5, 0.2,  1.5 } };//OZ
00188     UnitTetraIntersectionBary bary;
00189     fill_UnitTetraIntersectionBary(bary,nodes);
00190     double baryCenter[3];
00191     bool ok    = bary.getBary( baryCenter );
00192     double vol = bary.getVolume();
00193     CPPUNIT_ASSERT( ok );
00194     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0349217, vol, 1e-5);
00195     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.332275  , baryCenter[0], 1e-2 );
00196     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0565892 , baryCenter[1], 1e-3 );
00197     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.308713  , baryCenter[2], 1e-2 );
00198   }
00199   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_9()
00200   {
00201     // cutting tetra touches the UT at an edge, intersection volume == 0
00202     double nodes[4][3] = { { 1.0, 0.0, 0.0 }, // 0
00203                            {-1.0, 2.0, 2.0 }, // OX
00204                            {-1.0,-2.0, 2.0 }, // OY
00205                            { 1.0, 0.0, 2.0 } };//OZ
00206     UnitTetraIntersectionBary bary;
00207     fill_UnitTetraIntersectionBary(bary,nodes);
00208     double baryCenter[3];
00209     bool ok    = bary.getBary( baryCenter );
00210     double vol = bary.getVolume();
00211     CPPUNIT_ASSERT( !ok );
00212     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, vol, 1e-15);
00213     CPPUNIT_ASSERT_DOUBLES_EQUAL( -1. , baryCenter[0], 1e-5 );
00214     CPPUNIT_ASSERT_DOUBLES_EQUAL( -1. , baryCenter[1], 1e-5 );
00215     CPPUNIT_ASSERT_DOUBLES_EQUAL( -1. , baryCenter[2], 1e-5 );
00216   }
00217   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_10()
00218   {
00219     // cutting tetra fully includes theUT and touches it at an edge
00220     double nodes[4][3] = { { 1.0, 0.0, 0.0 }, // 0
00221                            {-1.0,-4.0, 2.0 }, // OX
00222                            {-1.0, 4.0, 2.0 }, // OY
00223                            { 1.0, 0.0,-2.0 } };//OZ
00224     UnitTetraIntersectionBary bary;
00225     fill_UnitTetraIntersectionBary(bary,nodes);
00226     double baryCenter[3];
00227     bool ok    = bary.getBary( baryCenter );
00228     double vol = bary.getVolume();
00229     CPPUNIT_ASSERT( ok );
00230     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.166667, vol, 1e-5);
00231     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[0], 1e-5);
00232     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[1], 1e-5);
00233     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.25, baryCenter[2], 1e-5);
00234   }
00235   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_11()
00236   {
00237     // cutting tetra intersects the UT and touches it at an edge
00238     double nodes[4][3] = { { 1.0, 0.0, 0.0 }, // 0
00239                            {-1.0,-4.0, 2.0 }, // OX
00240                            {-1.0, 4.0, 2.0 }, // OY
00241                            {-1.0, 0.0,-1.0 } };//OZ
00242     UnitTetraIntersectionBary bary;
00243     fill_UnitTetraIntersectionBary(bary,nodes);
00244     double baryCenter[3];
00245     bool ok    = bary.getBary( baryCenter );
00246     double vol = bary.getVolume();
00247     CPPUNIT_ASSERT( ok );
00248     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.15873 , vol, 1e-5);
00249     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.250000, baryCenter[0], 1e-5);
00250     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.230952, baryCenter[1], 1e-5);
00251     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.260714, baryCenter[2], 1e-5);
00252   }
00253   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_12()
00254   {
00255     // cutting tetra has one corner inside the UT and one its side passes through an UT edge
00256     double nodes[4][3] = { { 0.25, 0.25, 0.25 }, // 0
00257                            { 1.75,-0.25,-0.25 }, // OX
00258                            { 0.5 , 0.25, 0.25 }, // OY
00259                            { 0.5 , 0   , 0.5  } };//OZ
00260     UnitTetraIntersectionBary bary;
00261     fill_UnitTetraIntersectionBary(bary,nodes);
00262     double baryCenter[3];
00263     bool ok    = bary.getBary( baryCenter );
00264     double vol = bary.getVolume();
00265     CPPUNIT_ASSERT( ok );
00266     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.005208 , vol, 1e-5);
00267     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.562500, baryCenter[0], 1e-5);
00268     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.125000, baryCenter[1], 1e-5);
00269     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.250000, baryCenter[2], 1e-5);
00270   }
00271 
00272   struct __MESH_DUMMY
00273   {
00274     typedef int MyConnType;
00275   };
00276 
00277   void UnitTetraIntersectionBaryTest::test_UnitTetraIntersectionBary_13()
00278   {
00279     double T[] = {
00280       66.6666666666666714,133.333333333333343,66.6666666666666714,
00281       100,200,100,
00282       100,100,100,
00283       200,200,0 };
00284     
00285     double S[] = {
00286       100,166.666666666666657,66.6666666666666714,
00287       100,150,50,
00288       75,150,75,
00289       100,100,100};
00290 
00291     int conn[4] = { 0,1,2,3 };
00292     
00293     const double* tnodes[4]={ T, T+3, T+6, T+9 };
00294     const double* snodes[4]={ S, S+3, S+6, S+9 };
00295     
00296     __MESH_DUMMY dummyMesh;
00297     SplitterTetra<__MESH_DUMMY> src( dummyMesh, snodes, conn );
00298     double volume = src.intersectTetra( tnodes );
00299     CPPUNIT_ASSERT_DOUBLES_EQUAL(6944.4444444444443,volume,1e-9);
00300   }
00301 
00302   void UnitTetraIntersectionBaryTest::test_TetraAffineTransform_reverseApply()
00303   {
00304     double nodes[4][3] = { {-4.0, 9.0, 3.0 }, 
00305                            {11.0, 0.0, 2.0 }, 
00306                            { 0.0, 0.0, 0.0 }, 
00307                            { 2.0, 1.0,10.0 }};
00308     //    double pSrc[3] = { -4.0, 9.0, 3.0 };
00309     double pSrc[3] = { 40., -20., 100. };
00310     double pDest[] = {1,1,1};
00311     const double* n[4] = { &nodes[0][0], &nodes[1][0], &nodes[2][0], &nodes[3][0] };
00312     TetraAffineTransform a(&n[0]);
00313     a.apply( pDest, pSrc );
00314     a.reverseApply( pDest, pDest );
00315     CPPUNIT_ASSERT_DOUBLES_EQUAL( pSrc[0], pDest[0], 1e-12);
00316     CPPUNIT_ASSERT_DOUBLES_EQUAL( pSrc[1], pDest[1], 1e-12);
00317     CPPUNIT_ASSERT_DOUBLES_EQUAL( pSrc[2], pDest[2], 1e-12);
00318   }
00319 
00320   void UnitTetraIntersectionBaryTest::test_barycentric_coords()
00321   {
00322     // compute barycentric coordinates
00323     double nodes[4][3] = { {11.0, 0.0, 2.0 },
00324                            {-4.0, 9.0, 3.0 },
00325                            { 0.0, 0.0, 0.0 }, 
00326                            { 6.0, 1.0,10.0 }};
00327     std::vector<const double*> n (4);
00328     n[0] = &nodes[0][0];
00329     n[1] = &nodes[1][0];
00330     n[2] = &nodes[2][0];
00331     n[3] = &nodes[3][0];
00332     double p  [3] = { 2, 2, 5 }, bc[4];
00333     barycentric_coords(n, p, bc);
00334     double bcSum = 0;
00335     double p2 [3] = { 0,0,0 };
00336     for ( int i = 0; i < 4; ++i ) {
00337       bcSum += bc[i];
00338       p2[0] += bc[i]*n[i][0];
00339       p2[1] += bc[i]*n[i][1];
00340       p2[2] += bc[i]*n[i][2];
00341     }
00342     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1., bcSum, 1e-12);
00343     CPPUNIT_ASSERT_DOUBLES_EQUAL( p[0], p2[0], 1e-12);
00344     CPPUNIT_ASSERT_DOUBLES_EQUAL( p[1], p2[1], 1e-12);
00345     CPPUNIT_ASSERT_DOUBLES_EQUAL( p[2], p2[2], 1e-12);
00346   }  
00347 }