MOAB: Mesh Oriented datABase  (version 5.4.0)
TargetCalculatorTest.cpp
Go to the documentation of this file.
00001 /* *****************************************************************
00002     MESQUITE -- The Mesh Quality Improvement Toolkit
00003 
00004     Copyright 2009 Sandia National Laboratories.  Developed at the
00005     University of Wisconsin--Madison under SNL contract number
00006     624796.  The U.S. Government and the University of Wisconsin
00007     retain certain rights to this software.
00008 
00009     This library is free software; you can redistribute it and/or
00010     modify it under the terms of the GNU Lesser General Public
00011     License as published by the Free Software Foundation; either
00012     version 2.1 of the License, or (at your option) any later version.
00013 
00014     This library is distributed in the hope that it will be useful,
00015     but WITHOUT ANY WARRANTY; without even the implied warranty of
00016     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00017     Lesser General Public License for more details.
00018 
00019     You should have received a copy of the GNU Lesser General Public License
00020     (lgpl.txt) along with this library; if not, write to the Free Software
00021     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00022 
00023     (2009) kraftche@cae.wisc.edu
00024 
00025   ***************************************************************** */
00026 
00027 /** \file TargetCalculatorTest.cpp
00028  *  \brief
00029  *  \author Jason Kraftcheck
00030  */
00031 
00032 #include "Mesquite.hpp"
00033 #include "UnitUtil.hpp"
00034 #include "TargetCalculator.hpp"
00035 #include "MsqError.hpp"
00036 #include "PatchData.hpp"
00037 #include "IdealElements.hpp"
00038 #include "LinearHexahedron.hpp"
00039 #include "LinearTriangle.hpp"
00040 #include "ReferenceMesh.hpp"
00041 #include <iostream>
00042 
00043 using namespace MBMesquite;
00044 
00045 double EPS    = 0.05;
00046 double EPSBIG = 0.5;
00047 
00048 class TargetCalculatorTest : public CppUnit::TestFixture
00049 {
00050   private:
00051     CPPUNIT_TEST_SUITE( TargetCalculatorTest );
00052 
00053     CPPUNIT_TEST( test_factor_2D );
00054     CPPUNIT_TEST( test_factor_surface );
00055     CPPUNIT_TEST( test_factor_3D );
00056     CPPUNIT_TEST( test_factor_2D_zero );
00057     CPPUNIT_TEST( test_factor_surface_zero );
00058     CPPUNIT_TEST( test_factor_3D_zero );
00059 
00060     CPPUNIT_TEST( test_size_2D );
00061     CPPUNIT_TEST( test_size_surface );
00062     CPPUNIT_TEST( test_size_3D );
00063 
00064     CPPUNIT_TEST( test_skew_2D );
00065     CPPUNIT_TEST( test_skew_surface );
00066     CPPUNIT_TEST( test_skew_3D );
00067 
00068     CPPUNIT_TEST( test_aspect_2D );
00069     CPPUNIT_TEST( test_aspect_surface );
00070     CPPUNIT_TEST( test_aspect_3D );
00071 
00072     CPPUNIT_TEST( test_shape_2D );
00073     CPPUNIT_TEST( test_shape_surface );
00074     CPPUNIT_TEST( test_shape_3D );
00075 
00076     CPPUNIT_TEST( test_ideal_skew_tri );
00077     CPPUNIT_TEST( test_ideal_skew_quad );
00078     CPPUNIT_TEST( test_ideal_skew_tet );
00079     CPPUNIT_TEST( test_ideal_skew_prism );
00080     CPPUNIT_TEST( test_ideal_skew_pyramid );
00081     CPPUNIT_TEST( test_ideal_skew_hex );
00082 
00083     CPPUNIT_TEST( test_ideal_shape_tri );
00084     CPPUNIT_TEST( test_ideal_shape_quad );
00085     CPPUNIT_TEST( test_ideal_shape_tet );
00086     CPPUNIT_TEST( test_ideal_shape_prism );
00087     CPPUNIT_TEST( test_ideal_shape_pyramid );
00088     CPPUNIT_TEST( test_ideal_shape_hex );
00089 
00090     CPPUNIT_TEST( test_new_orientatin_3D );
00091     CPPUNIT_TEST( test_new_orientatin_2D );
00092 
00093     CPPUNIT_TEST( test_new_aspect_3D );
00094     CPPUNIT_TEST( test_new_aspect_2D );
00095 
00096     CPPUNIT_TEST( test_jacobian_3D );
00097     CPPUNIT_TEST( test_jacobian_2D );
00098 
00099     CPPUNIT_TEST( test_get_refmesh_Jacobian_3D );
00100     CPPUNIT_TEST( test_get_refmesh_Jacobian_2D );
00101 
00102     CPPUNIT_TEST_SUITE_END();
00103 
00104     // define some matrices for use in testing.
00105     MsqMatrix< 3, 3 > V3D_Z45, V3D_X90, Q3D_45, D3D_123;
00106     MsqMatrix< 3, 2 > V2D_Z45, V2D_X90;
00107     MsqMatrix< 2, 2 > Q2D_45, D2D_21;
00108 
00109   public:
00110     TargetCalculatorTest();
00111     void setUp();
00112 
00113     void test_factor_2D();
00114     void test_factor_surface();
00115     void test_factor_3D();
00116     void test_factor_2D_zero();
00117     void test_factor_surface_zero();
00118     void test_factor_3D_zero();
00119 
00120     void test_size_2D();
00121     void test_size_surface();
00122     void test_size_3D();
00123 
00124     void test_skew_2D();
00125     void test_skew_surface();
00126     void test_skew_3D();
00127 
00128     void test_aspect_2D();
00129     void test_aspect_surface();
00130     void test_aspect_3D();
00131 
00132     void test_shape_2D();
00133     void test_shape_surface();
00134     void test_shape_3D();
00135 
00136     void test_ideal_skew_tri();
00137     void test_ideal_skew_quad();
00138     void test_ideal_skew_tet();
00139     void test_ideal_skew_prism();
00140     void test_ideal_skew_pyramid();
00141     void test_ideal_skew_hex();
00142 
00143     void test_ideal_shape_tri();
00144     void test_ideal_shape_quad();
00145     void test_ideal_shape_tet();
00146     void test_ideal_shape_prism();
00147     void test_ideal_shape_pyramid();
00148     void test_ideal_shape_hex();
00149 
00150     void test_new_orientatin_3D();
00151     void test_new_orientatin_2D();
00152 
00153     void test_new_aspect_3D();
00154     void test_new_aspect_2D();
00155 
00156     void test_jacobian_3D();
00157     void test_jacobian_2D();
00158 
00159     void test_get_refmesh_Jacobian_3D();
00160     void test_get_refmesh_Jacobian_2D();
00161 
00162     template < unsigned R, unsigned C >
00163     static inline void check_valid_V( MsqMatrix< R, C > V );
00164 
00165     template < unsigned D >
00166     static inline void check_valid_Q( MsqMatrix< D, D > Q );
00167 
00168     template < unsigned D >
00169     static inline void check_valid_delta( MsqMatrix< D, D > delta );
00170 };
00171 
00172 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION( TargetCalculatorTest, "TargetCalculatorTest" );
00173 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION( TargetCalculatorTest, "Unit" );
00174 
00175 TargetCalculatorTest::TargetCalculatorTest()
00176 {
00177     const double cos45              = MSQ_SQRT_TWO / 2.0;
00178     const double rotation_3D_Z45[9] = { cos45, -cos45, 0, cos45, cos45, 0, 0, 0, 1 };
00179     const double rotation_2D_Z45[6] = { cos45, -cos45, cos45, cos45, 0, 0 };
00180 
00181     const double rotation_3D_X90[9] = { 1, 0, 0, 0, 0, -1, 0, 1, 0 };
00182     const double rotation_2D_X90[6] = { 1, 0, 0, 0, 0, 1 };
00183 
00184     const double rc45          = sqrt( cos45 );
00185     const double skew_2D_45[4] = { 1 / rc45, rc45, 0, rc45 };
00186     const double skew_3D_45[9] = { 1, cos45, cos45, 0, cos45, 1 - cos45, 0, 0, sqrt( MSQ_SQRT_TWO - 1 ) };
00187 
00188     const double aspect_2D_2x[4]  = { MSQ_SQRT_TWO, 0, 0, MSQ_SQRT_TWO / 2 };
00189     const double r6               = MBMesquite::cbrt( 1.0 / 6.0 );
00190     const double aspect_3D_123[9] = { r6, 0, 0, 0, 2 * r6, 0, 0, 0, 3 * r6 };
00191 
00192     V3D_Z45 = MsqMatrix< 3, 3 >( rotation_3D_Z45 );
00193     V3D_X90 = MsqMatrix< 3, 3 >( rotation_3D_X90 );
00194     Q3D_45  = MsqMatrix< 3, 3 >( skew_3D_45 );
00195     Q3D_45 *= 1 / Mesquite::cbrt( det( Q3D_45 ) );
00196     D3D_123 = MsqMatrix< 3, 3 >( aspect_3D_123 );
00197 
00198     V2D_Z45 = MsqMatrix< 3, 2 >( rotation_2D_Z45 );
00199     V2D_X90 = MsqMatrix< 3, 2 >( rotation_2D_X90 );
00200     Q2D_45  = MsqMatrix< 2, 2 >( skew_2D_45 );
00201     D2D_21  = MsqMatrix< 2, 2 >( aspect_2D_2x );
00202 }
00203 
00204 template < unsigned R, unsigned C >
00205 inline void TargetCalculatorTest::check_valid_V( MsqMatrix< R, C > V )
00206 {
00207     // check that it is a rotation
00208     MsqMatrix< C, C > I( 1.0 );
00209     ASSERT_MATRICES_EQUAL( I, transpose( V ) * V, EPS );
00210 }
00211 
00212 template < unsigned D >
00213 inline void TargetCalculatorTest::check_valid_Q( MsqMatrix< D, D > Q )
00214 {
00215     // must have unit determinant
00216     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, det( Q ), EPS );
00217 
00218     // must be upper-triangular
00219     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, Q( 1, 0 ), EPS );
00220     if( D == 3 )
00221     {
00222         CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, Q( 2, 0 ), EPS );
00223         CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, Q( 2, 1 ), EPS );
00224     }
00225 
00226     // columns must be of equal length
00227     CPPUNIT_ASSERT_DOUBLES_EQUAL( length( Q.column( 0 ) ), length( Q.column( 1 ) ), EPSBIG );
00228     if( D == 3 )
00229     {
00230         CPPUNIT_ASSERT_DOUBLES_EQUAL( length( Q.column( 0 ) ), length( Q.column( 2 ) ), EPSBIG );
00231     }
00232 
00233     // diagonal elements must be greater than zero
00234     CPPUNIT_ASSERT( Q( 0, 0 ) - EPS >= 0.0 );
00235     CPPUNIT_ASSERT( Q( 1, 1 ) - EPS >= 0.0 );
00236     if( D == 3 )
00237     {
00238         CPPUNIT_ASSERT( Q( 2, 2 ) - EPS >= 0.0 );
00239     }
00240 }
00241 
00242 template < unsigned D >
00243 inline void TargetCalculatorTest::check_valid_delta( MsqMatrix< D, D > delta )
00244 {
00245     // must have unit determinant
00246     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, det( delta ), EPS );
00247 
00248     // must be diagonal
00249     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, delta( 1, 0 ), EPS );
00250     CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, delta( 0, 1 ), EPS );
00251     if( D == 3 )
00252     {
00253         CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, delta( 2, 0 ), EPS );
00254         CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, delta( 0, 2 ), EPS );
00255         CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, delta( 1, 2 ), EPS );
00256         CPPUNIT_ASSERT_DOUBLES_EQUAL( 0.0, delta( 2, 1 ), EPS );
00257     }
00258 }
00259 
00260 void TargetCalculatorTest::setUp()
00261 {
00262     // test the test: make sure we're testing with
00263     // valid matrix factors.
00264 
00265     check_valid_V( V3D_Z45 );
00266     check_valid_V( V3D_X90 );
00267     check_valid_V( V2D_Z45 );
00268     check_valid_V( V2D_X90 );
00269 
00270     check_valid_Q( Q3D_45 );
00271     check_valid_Q( Q2D_45 );
00272 
00273     check_valid_delta( D3D_123 );
00274     check_valid_delta( D2D_21 );
00275 }
00276 
00277 void TargetCalculatorTest::test_factor_2D()
00278 {
00279     MsqPrintError err( std::cout );
00280     MsqMatrix< 2, 2 > I( 1.0 ), W, V;
00281     double lambda;
00282     MsqMatrix< 2, 2 > Q, delta;
00283     bool valid;
00284 
00285     // first test with I
00286     valid = TargetCalculator::factor_2D( I, lambda, V, Q, delta, err );
00287     CPPUNIT_ASSERT( valid && !err );
00288     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, lambda, EPS );
00289     check_valid_V( V );
00290     check_valid_Q( Q );
00291     check_valid_delta( delta );
00292     W = lambda * V * Q * delta;
00293     ASSERT_MATRICES_EQUAL( I, W, EPS );
00294 
00295     // now test with 2*I
00296     valid = TargetCalculator::factor_2D( 2 * I, lambda, V, Q, delta, err );
00297     CPPUNIT_ASSERT( valid && !err );
00298     CPPUNIT_ASSERT_DOUBLES_EQUAL( 2.0, lambda, EPS );
00299     check_valid_V( V );
00300     check_valid_Q( Q );
00301     check_valid_delta( delta );
00302     W = lambda * V * Q * delta;
00303     ASSERT_MATRICES_EQUAL( 2 * I, W, EPS );
00304 
00305     // now rotate the matrix about the X axis by 90 degrees
00306     MsqMatrix< 2, 2 > I_rot( 0.0 );
00307     I_rot( 0, 1 ) = 1.0;
00308     I_rot( 1, 0 ) = -1.0;
00309     valid         = TargetCalculator::factor_2D( I_rot, lambda, V, Q, delta, err );
00310     CPPUNIT_ASSERT( valid && !err );
00311     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, lambda, EPS );
00312     check_valid_V( V );
00313     check_valid_Q( Q );
00314     check_valid_delta( delta );
00315     W = lambda * V * Q * delta;
00316     ASSERT_MATRICES_EQUAL( I_rot, W, EPS );
00317 
00318     // change the aspect ratio
00319     MsqMatrix< 2, 2 > A( 1.0 );
00320     A( 1, 1 ) = 1.5;
00321     valid     = TargetCalculator::factor_2D( A, lambda, V, Q, delta, err );
00322     CPPUNIT_ASSERT( valid && !err );
00323     CPPUNIT_ASSERT_DOUBLES_EQUAL( pow( det( transpose( A ) * A ), 0.25 ), lambda, EPS );
00324     check_valid_V( V );
00325     check_valid_Q( Q );
00326     check_valid_delta( delta );
00327     W = lambda * V * Q * delta;
00328     ASSERT_MATRICES_EQUAL( A, W, EPS );
00329 
00330     // try an arbitrary matrix
00331     double w[4] = { 3, 1, 4, 2 };
00332     MsqMatrix< 2, 2 > W2( w );
00333     valid = TargetCalculator::factor_2D( W2, lambda, V, Q, delta, err );
00334     CPPUNIT_ASSERT( valid && !err );
00335     CPPUNIT_ASSERT_DOUBLES_EQUAL( pow( det( transpose( W2 ) * W2 ), 0.25 ), lambda, EPS );
00336     check_valid_V( V );
00337     check_valid_Q( Q );
00338     check_valid_delta( delta );
00339     W = lambda * V * Q * delta;
00340     ASSERT_MATRICES_EQUAL( W2, W, EPS );
00341 
00342     // try a more elaborate test
00343     const double e = exp( (double)1 );
00344     MsqMatrix< 2, 2 > Z45( V2D_Z45.data() );  // copy first two rows
00345     W2    = e * Z45 * Q2D_45 * D2D_21;
00346     valid = TargetCalculator::factor_2D( W2, lambda, V, Q, delta, err );
00347     CPPUNIT_ASSERT( valid && !err );
00348     CPPUNIT_ASSERT_DOUBLES_EQUAL( e, lambda, EPS );
00349     ASSERT_MATRICES_EQUAL( Z45, V, EPS );
00350     ASSERT_MATRICES_EQUAL( Q2D_45, Q, EPS );
00351     ASSERT_MATRICES_EQUAL( D2D_21, delta, EPS );
00352 }
00353 
00354 void TargetCalculatorTest::test_factor_surface()
00355 {
00356     MsqPrintError err( std::cout );
00357     MsqMatrix< 3, 2 > I( 1.0 ), W, V;
00358     double lambda;
00359     MsqMatrix< 2, 2 > Q, delta;
00360     bool valid;
00361 
00362     // first test with I
00363     valid = TargetCalculator::factor_surface( I, lambda, V, Q, delta, err );
00364     CPPUNIT_ASSERT( valid && !err );
00365     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, lambda, EPS );
00366     check_valid_V( V );
00367     check_valid_Q( Q );
00368     check_valid_delta( delta );
00369     W = lambda * V * Q * delta;
00370     ASSERT_MATRICES_EQUAL( I, W, EPS );
00371 
00372     // now test with 2*I
00373     valid = TargetCalculator::factor_surface( 2 * I, lambda, V, Q, delta, err );
00374     CPPUNIT_ASSERT( valid && !err );
00375     CPPUNIT_ASSERT_DOUBLES_EQUAL( 2.0, lambda, EPS );
00376     check_valid_V( V );
00377     check_valid_Q( Q );
00378     check_valid_delta( delta );
00379     W = lambda * V * Q * delta;
00380     ASSERT_MATRICES_EQUAL( 2 * I, W, EPS );
00381 
00382     // now rotate the matrix about the X axis by 90 degrees
00383     MsqMatrix< 3, 2 > I_rot( 1.0 );
00384     I_rot( 1, 1 ) = 0.0;
00385     I_rot( 2, 1 ) = 1.0;
00386     valid         = TargetCalculator::factor_surface( I_rot, lambda, V, Q, delta, err );
00387     CPPUNIT_ASSERT( valid && !err );
00388     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, lambda, EPS );
00389     check_valid_V( V );
00390     check_valid_Q( Q );
00391     check_valid_delta( delta );
00392     W = lambda * V * Q * delta;
00393     ASSERT_MATRICES_EQUAL( I_rot, W, EPS );
00394 
00395     // change the aspect ratio
00396     MsqMatrix< 3, 2 > A( 1.0 );
00397     A( 1, 1 ) = 1.5;
00398     valid     = TargetCalculator::factor_surface( A, lambda, V, Q, delta, err );
00399     CPPUNIT_ASSERT( valid && !err );
00400     CPPUNIT_ASSERT_DOUBLES_EQUAL( pow( det( transpose( A ) * A ), 0.25 ), lambda, EPS );
00401     check_valid_V( V );
00402     check_valid_Q( Q );
00403     check_valid_delta( delta );
00404     W = lambda * V * Q * delta;
00405     ASSERT_MATRICES_EQUAL( A, W, EPS );
00406 
00407     // try an arbitrary matrix
00408     double w[6] = { 3, 1, 4, 2, -1, 5 };
00409     MsqMatrix< 3, 2 > W2( w );
00410     valid = TargetCalculator::factor_surface( W2, lambda, V, Q, delta, err );
00411     CPPUNIT_ASSERT( valid && !err );
00412     CPPUNIT_ASSERT_DOUBLES_EQUAL( pow( det( transpose( W2 ) * W2 ), 0.25 ), lambda, EPS );
00413     check_valid_V( V );
00414     check_valid_Q( Q );
00415     check_valid_delta( delta );
00416     W = lambda * V * Q * delta;
00417     ASSERT_MATRICES_EQUAL( W2, W, EPS );
00418 
00419     // try a more elaborate test
00420     const double e = exp( (double)1 );
00421     W2             = e * V2D_Z45 * Q2D_45 * D2D_21;
00422     valid          = TargetCalculator::factor_surface( W2, lambda, V, Q, delta, err );
00423     CPPUNIT_ASSERT( valid && !err );
00424     CPPUNIT_ASSERT_DOUBLES_EQUAL( e, lambda, EPS );
00425     ASSERT_MATRICES_EQUAL( V2D_Z45, V, EPS );
00426     ASSERT_MATRICES_EQUAL( Q2D_45, Q, EPS );
00427     ASSERT_MATRICES_EQUAL( D2D_21, delta, EPS );
00428 }
00429 
00430 void TargetCalculatorTest::test_factor_3D()
00431 {
00432     MsqPrintError err( std::cout );
00433     MsqMatrix< 3, 3 > I( 1.0 ), W, V, Q, delta;
00434     double lambda;
00435     bool valid;
00436 
00437     // first test with I
00438     valid = TargetCalculator::factor_3D( I, lambda, V, Q, delta, err );
00439     CPPUNIT_ASSERT( valid && !err );
00440     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, lambda, EPS );
00441     check_valid_V( V );
00442     check_valid_Q( Q );
00443     check_valid_delta( delta );
00444     W = lambda * V * Q * delta;
00445     ASSERT_MATRICES_EQUAL( I, W, EPS );
00446 
00447     // now test with 2*I
00448     valid = TargetCalculator::factor_3D( 2 * I, lambda, V, Q, delta, err );
00449     CPPUNIT_ASSERT( valid && !err );
00450     CPPUNIT_ASSERT_DOUBLES_EQUAL( 2.0, lambda, EPS );
00451     check_valid_V( V );
00452     check_valid_Q( Q );
00453     check_valid_delta( delta );
00454     W = lambda * V * Q * delta;
00455     ASSERT_MATRICES_EQUAL( 2 * I, W, EPS );
00456 
00457     // now rotate the matrix about the X axis by 90 degrees
00458     MsqMatrix< 3, 3 > I_rot( 1.0 );
00459     I_rot( 1, 1 ) = 0.0;
00460     I_rot( 2, 1 ) = 1.0;
00461     I_rot( 1, 2 ) = -1.0;
00462     I_rot( 2, 2 ) = 0.0;
00463     valid         = TargetCalculator::factor_3D( I_rot, lambda, V, Q, delta, err );
00464     CPPUNIT_ASSERT( valid && !err );
00465     CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, lambda, EPS );
00466     check_valid_V( V );
00467     check_valid_Q( Q );
00468     check_valid_delta( delta );
00469     W = lambda * V * Q * delta;
00470     ASSERT_MATRICES_EQUAL( I_rot, W, EPS );
00471 
00472     // change the aspect ratio
00473     MsqMatrix< 3, 3 > A( 1.0 );
00474     A( 1, 1 ) = 1.5;
00475     valid     = TargetCalculator::factor_3D( A, lambda, V, Q, delta, err );
00476     CPPUNIT_ASSERT( valid && !err );
00477     CPPUNIT_ASSERT_DOUBLES_EQUAL( MBMesquite::cbrt( det( A ) ), lambda, EPS );
00478     check_valid_V( V );
00479     check_valid_Q( Q );
00480     check_valid_delta( delta );
00481     W = lambda * V * Q * delta;
00482     ASSERT_MATRICES_EQUAL( A, W, EPS );
00483 
00484     // try an arbitrary matrix
00485     double vals[] = { 9, 8, 7, 1, 5, 4, 3, 2, 6 };
00486     MsqMatrix< 3, 3 > W2( vals );
00487     valid = TargetCalculator::factor_3D( W2, lambda, V, Q, delta, err );
00488     CPPUNIT_ASSERT( valid && !err );
00489     CPPUNIT_ASSERT_DOUBLES_EQUAL( MBMesquite::cbrt( det( W2 ) ), lambda, EPS );
00490     check_valid_V( V );
00491     check_valid_Q( Q );
00492     check_valid_delta( delta );
00493     W = lambda * V * Q * delta;
00494     ASSERT_MATRICES_EQUAL( W2, W, EPSBIG );
00495 
00496     // try a more elaborate test
00497     const double e = exp( (double)1 );
00498     W2             = e * V3D_Z45 * Q3D_45 * D3D_123;
00499     valid          = TargetCalculator::factor_3D( W2, lambda, V, Q, delta, err );
00500     CPPUNIT_ASSERT( valid && !err );
00501     CPPUNIT_ASSERT_DOUBLES_EQUAL( e, lambda, EPS );
00502     ASSERT_MATRICES_EQUAL( V3D_Z45, V, EPS );
00503     ASSERT_MATRICES_EQUAL( Q3D_45, Q, EPS );
00504     ASSERT_MATRICES_EQUAL( D3D_123, delta, EPS );
00505 }
00506 
00507 void TargetCalculatorTest::test_factor_2D_zero()
00508 {
00509     MsqPrintError err( std::cout );
00510     double lambda;
00511     MsqMatrix< 2, 2 > Q, delta;
00512     MsqMatrix< 2, 2 > V, Z( 0.0 );
00513     Z( 0, 0 ) = 1.0;
00514 
00515     bool valid = TargetCalculator::factor_2D( Z, lambda, V, Q, delta, err );
00516     ASSERT_NO_ERROR( err );
00517     CPPUNIT_ASSERT( !valid );
00518 }
00519 
00520 void TargetCalculatorTest::test_factor_surface_zero()
00521 {
00522     MsqPrintError err( std::cout );
00523     double lambda;
00524     MsqMatrix< 2, 2 > Q, delta;
00525     MsqMatrix< 3, 2 > V, Z( 0.0 );
00526     Z( 0, 0 ) = 1.0;
00527 
00528     bool valid = TargetCalculator::factor_surface( Z, lambda, V, Q, delta, err );
00529     ASSERT_NO_ERROR( err );
00530     CPPUNIT_ASSERT( !valid );
00531 }
00532 
00533 void TargetCalculatorTest::test_factor_3D_zero()
00534 {
00535     MsqPrintError err( std::cout );
00536     double lambda;
00537     MsqMatrix< 3, 3 > V, Q, delta, Z( 0.0 );
00538     Z( 0, 0 ) = 1.0;
00539 
00540     bool valid = TargetCalculator::factor_3D( Z, lambda, V, Q, delta, err );
00541     ASSERT_NO_ERROR( err );
00542     CPPUNIT_ASSERT( !valid );
00543 }
00544 
00545 void TargetCalculatorTest::test_size_2D()
00546 {
00547     MsqPrintError err( std::cout );
00548     double lambda;
00549     MsqMatrix< 2, 2 > V, Q, delta;
00550 
00551     MsqMatrix< 2, 2 > I( 1.0 );
00552     TargetCalculator::factor_2D( I, lambda, V, Q, delta, err );
00553     ASSERT_NO_ERROR( err );
00554     CPPUNIT_ASSERT_DOUBLES_EQUAL( lambda, TargetCalculator::size( I ), EPS );
00555 
00556     double w[] = { 3, 1, 4, 2 };
00557     MsqMatrix< 2, 2 > W( w );
00558     TargetCalculator::factor_2D( W, lambda, V, Q, delta, err );
00559     ASSERT_NO_ERROR( err );
00560     CPPUNIT_ASSERT_DOUBLES_EQUAL( lambda, TargetCalculator::size( W ), EPS );
00561 }
00562 
00563 void TargetCalculatorTest::test_size_surface()
00564 {
00565     MsqPrintError err( std::cout );
00566     double lambda;
00567     MsqMatrix< 3, 2 > V;
00568     MsqMatrix< 2, 2 > Q, delta;
00569 
00570     MsqMatrix< 3, 2 > I( 1.0 );
00571     TargetCalculator::factor_surface( I, lambda, V, Q, delta, err );
00572     ASSERT_NO_ERROR( err );
00573     CPPUNIT_ASSERT_DOUBLES_EQUAL( lambda, TargetCalculator::size( I ), EPS );
00574 
00575     double w[] = { 3, 1, 4, 2, -1, 5 };
00576     MsqMatrix< 3, 2 > W( w );
00577     TargetCalculator::factor_surface( W, lambda, V, Q, delta, err );
00578     ASSERT_NO_ERROR( err );
00579     CPPUNIT_ASSERT_DOUBLES_EQUAL( lambda, TargetCalculator::size( W ), EPS );
00580 }
00581 
00582 void TargetCalculatorTest::test_size_3D()
00583 {
00584     MsqPrintError err( std::cout );
00585     double lambda;
00586     MsqMatrix< 3, 3 > V, Q, delta;
00587 
00588     MsqMatrix< 3, 3 > I( 1.0 );
00589     TargetCalculator::factor_3D( I, lambda, V, Q, delta, err );
00590     ASSERT_NO_ERROR( err );
00591     CPPUNIT_ASSERT_DOUBLES_EQUAL( lambda, TargetCalculator::size( I ), EPS );
00592 
00593     double w[] = { 9, 8, 7, 1, 5, 4, 3, 2, 6 };
00594     MsqMatrix< 3, 3 > W( w );
00595     TargetCalculator::factor_3D( W, lambda, V, Q, delta, err );
00596     ASSERT_NO_ERROR( err );
00597     CPPUNIT_ASSERT_DOUBLES_EQUAL( lambda, TargetCalculator::size( W ), EPS );
00598 }
00599 
00600 void TargetCalculatorTest::test_skew_2D()
00601 {
00602     MsqPrintError err( std::cout );
00603     double lambda;
00604     MsqMatrix< 2, 2 > V, Q, delta;
00605 
00606     MsqMatrix< 2, 2 > I( 1.0 );
00607     TargetCalculator::factor_2D( I, lambda, V, Q, delta, err );
00608     ASSERT_NO_ERROR( err );
00609     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( I ), EPS );
00610 
00611     double r[] = { 0, -1, 2, 0 };
00612     MsqMatrix< 2, 2 > R( r );
00613     TargetCalculator::factor_2D( R, lambda, V, Q, delta, err );
00614     ASSERT_NO_ERROR( err );
00615     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( R ), EPS );
00616 
00617     double w[] = { 3, 1, 4, 2 };
00618     MsqMatrix< 2, 2 > W( w );
00619     TargetCalculator::factor_2D( W, lambda, V, Q, delta, err );
00620     ASSERT_NO_ERROR( err );
00621     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( W ), EPS );
00622 
00623     W = MsqMatrix< 2, 2 >( 6 ) * Q2D_45;
00624     ASSERT_MATRICES_EQUAL( Q2D_45, TargetCalculator::skew( W ), EPS );
00625 }
00626 
00627 void TargetCalculatorTest::test_skew_surface()
00628 {
00629     MsqPrintError err( std::cout );
00630     double lambda;
00631     MsqMatrix< 3, 2 > V;
00632     MsqMatrix< 2, 2 > Q, delta;
00633 
00634     MsqMatrix< 3, 2 > I( 1.0 );
00635     TargetCalculator::factor_surface( I, lambda, V, Q, delta, err );
00636     ASSERT_NO_ERROR( err );
00637     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( I ), EPS );
00638 
00639     double r[] = { 1, 0, 0, 0, 0, 2 };
00640     MsqMatrix< 3, 2 > R( r );
00641     TargetCalculator::factor_surface( R, lambda, V, Q, delta, err );
00642     ASSERT_NO_ERROR( err );
00643     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( R ), EPS );
00644 
00645     double w[] = { 3, 1, 4, 2, -1, 5 };
00646     MsqMatrix< 3, 2 > W( w );
00647     TargetCalculator::factor_surface( W, lambda, V, Q, delta, err );
00648     ASSERT_NO_ERROR( err );
00649     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( W ), EPS );
00650 
00651     W = MsqMatrix< 3, 2 >( 6 ) * Q2D_45;
00652     ASSERT_MATRICES_EQUAL( Q2D_45, TargetCalculator::skew( W ), EPS );
00653 }
00654 
00655 void TargetCalculatorTest::test_skew_3D()
00656 {
00657     MsqPrintError err( std::cout );
00658     double lambda;
00659     MsqMatrix< 3, 3 > V, Q, delta;
00660 
00661     MsqMatrix< 3, 3 > I( 1.0 );
00662     TargetCalculator::factor_3D( I, lambda, V, Q, delta, err );
00663     ASSERT_NO_ERROR( err );
00664     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( I ), EPS );
00665 
00666     double r[] = { 0.8, 0.0, 0.0, 0.0, 0.0, 1.13, 0.0, -0.5, 0.0 };
00667     MsqMatrix< 3, 3 > R( r );
00668     TargetCalculator::factor_3D( R, lambda, V, Q, delta, err );
00669     ASSERT_NO_ERROR( err );
00670     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( R ), EPS );
00671 
00672     double w[] = { 9, 8, 7, 1, 5, 4, 3, 2, 6 };
00673     MsqMatrix< 3, 3 > W( w );
00674     TargetCalculator::factor_3D( W, lambda, V, Q, delta, err );
00675     ASSERT_NO_ERROR( err );
00676     ASSERT_MATRICES_EQUAL( Q, TargetCalculator::skew( W ), 2 * EPS );
00677 
00678     W = MsqMatrix< 3, 3 >( 2.1 ) * Q3D_45;
00679     ASSERT_MATRICES_EQUAL( Q3D_45, TargetCalculator::skew( W ), EPS );
00680 }
00681 
00682 void TargetCalculatorTest::test_aspect_2D()
00683 {
00684     MsqPrintError err( std::cout );
00685     double lambda;
00686     MsqMatrix< 2, 2 > V, Q, delta;
00687 
00688     MsqMatrix< 2, 2 > I( 1.0 );
00689     TargetCalculator::factor_2D( I, lambda, V, Q, delta, err );
00690     ASSERT_NO_ERROR( err );
00691     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( I ), EPS );
00692 
00693     double r[] = { 0, -1, 2, 0 };
00694     MsqMatrix< 2, 2 > R( r );
00695     TargetCalculator::factor_2D( R, lambda, V, Q, delta, err );
00696     ASSERT_NO_ERROR( err );
00697     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( R ), EPS );
00698 
00699     double w[] = { 3, 1, 4, 2 };
00700     MsqMatrix< 2, 2 > W( w );
00701     TargetCalculator::factor_2D( W, lambda, V, Q, delta, err );
00702     ASSERT_NO_ERROR( err );
00703     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( W ), EPS );
00704 
00705     W = MsqMatrix< 2, 2 >( 6 ) * Q2D_45;
00706     ASSERT_MATRICES_EQUAL( I, TargetCalculator::aspect( W ), EPS );
00707 }
00708 
00709 void TargetCalculatorTest::test_aspect_surface()
00710 {
00711     MsqPrintError err( std::cout );
00712     double lambda;
00713     MsqMatrix< 3, 2 > V;
00714     MsqMatrix< 2, 2 > Q, delta;
00715 
00716     MsqMatrix< 3, 2 > I( 1.0 );
00717     TargetCalculator::factor_surface( I, lambda, V, Q, delta, err );
00718     ASSERT_NO_ERROR( err );
00719     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( I ), EPS );
00720 
00721     double r[] = { 1, 0, 0, 0, 0, 2 };
00722     MsqMatrix< 3, 2 > R( r );
00723     TargetCalculator::factor_surface( R, lambda, V, Q, delta, err );
00724     ASSERT_NO_ERROR( err );
00725     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( R ), EPS );
00726 
00727     double w[] = { 3, 1, 4, 2, -1, 5 };
00728     MsqMatrix< 3, 2 > W( w );
00729     TargetCalculator::factor_surface( W, lambda, V, Q, delta, err );
00730     ASSERT_NO_ERROR( err );
00731     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( W ), EPS );
00732 
00733     W = MsqMatrix< 3, 2 >( 6 ) * Q2D_45;
00734     ASSERT_MATRICES_EQUAL( ( MsqMatrix< 2, 2 >( 1.0 ) ), TargetCalculator::aspect( W ), EPS );
00735 }
00736 
00737 void TargetCalculatorTest::test_aspect_3D()
00738 {
00739     MsqPrintError err( std::cout );
00740     double lambda;
00741     MsqMatrix< 3, 3 > V, Q, delta;
00742 
00743     MsqMatrix< 3, 3 > I( 1.0 );
00744     TargetCalculator::factor_3D( I, lambda, V, Q, delta, err );
00745     ASSERT_NO_ERROR( err );
00746     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( I ), EPS );
00747 
00748     double r[] = { 0.8, 0.0, 0.0, 0.0, 0.0, 1.13, 0.0, -0.5, 0.0 };
00749     MsqMatrix< 3, 3 > R( r );
00750     TargetCalculator::factor_3D( R, lambda, V, Q, delta, err );
00751     ASSERT_NO_ERROR( err );
00752     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( R ), EPS );
00753 
00754     double w[] = { 9, 8, 7, 1, 5, 4, 3, 2, 6 };
00755     MsqMatrix< 3, 3 > W( w );
00756     TargetCalculator::factor_3D( W, lambda, V, Q, delta, err );
00757     ASSERT_NO_ERROR( err );
00758     ASSERT_MATRICES_EQUAL( delta, TargetCalculator::aspect( W ), 2 * EPS );
00759 
00760     W = MsqMatrix< 3, 3 >( 2.1 ) * Q3D_45;
00761     ASSERT_MATRICES_EQUAL( I, TargetCalculator::aspect( W ), EPS );
00762 }
00763 
00764 void TargetCalculatorTest::test_shape_2D()
00765 {
00766     MsqPrintError err( std::cout );
00767     double lambda;
00768     MsqMatrix< 2, 2 > V, Q, delta;
00769 
00770     MsqMatrix< 2, 2 > I( 1.0 );
00771     TargetCalculator::factor_2D( I, lambda, V, Q, delta, err );
00772     ASSERT_NO_ERROR( err );
00773     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( I ), EPS );
00774 
00775     double r[] = { 0, -1, 2, 0 };
00776     MsqMatrix< 2, 2 > R( r );
00777     TargetCalculator::factor_2D( R, lambda, V, Q, delta, err );
00778     ASSERT_NO_ERROR( err );
00779     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( R ), EPS );
00780 
00781     double w[] = { 3, 1, 4, 2 };
00782     MsqMatrix< 2, 2 > W( w );
00783     TargetCalculator::factor_2D( W, lambda, V, Q, delta, err );
00784     ASSERT_NO_ERROR( err );
00785     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( W ), EPS );
00786 
00787     W = MsqMatrix< 2, 2 >( 6 ) * Q2D_45;
00788     ASSERT_MATRICES_EQUAL( Q2D_45, TargetCalculator::shape( W ), EPS );
00789 }
00790 
00791 void TargetCalculatorTest::test_shape_surface()
00792 {
00793     MsqPrintError err( std::cout );
00794     double lambda;
00795     MsqMatrix< 3, 2 > V;
00796     MsqMatrix< 2, 2 > Q, delta;
00797 
00798     MsqMatrix< 3, 2 > I( 1.0 );
00799     TargetCalculator::factor_surface( I, lambda, V, Q, delta, err );
00800     ASSERT_NO_ERROR( err );
00801     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( I ), EPS );
00802 
00803     double r[] = { 1, 0, 0, 0, 0, 2 };
00804     MsqMatrix< 3, 2 > R( r );
00805     TargetCalculator::factor_surface( R, lambda, V, Q, delta, err );
00806     ASSERT_NO_ERROR( err );
00807     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( R ), EPS );
00808 
00809     double w[] = { 3, 1, 4, 2, -1, 5 };
00810     MsqMatrix< 3, 2 > W( w );
00811     TargetCalculator::factor_surface( W, lambda, V, Q, delta, err );
00812     ASSERT_NO_ERROR( err );
00813     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( W ), EPS );
00814 
00815     W = MsqMatrix< 3, 2 >( 6 ) * Q2D_45;
00816     ASSERT_MATRICES_EQUAL( Q2D_45, TargetCalculator::shape( W ), EPS );
00817 }
00818 
00819 void TargetCalculatorTest::test_shape_3D()
00820 {
00821     MsqPrintError err( std::cout );
00822     double lambda;
00823     MsqMatrix< 3, 3 > V, Q, delta;
00824 
00825     MsqMatrix< 3, 3 > I( 1.0 );
00826     TargetCalculator::factor_3D( I, lambda, V, Q, delta, err );
00827     ASSERT_NO_ERROR( err );
00828     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( I ), EPS );
00829 
00830     double r[] = { 0.8, 0.0, 0.0, 0.0, 0.0, 1.13, 0.0, -0.5, 0.0 };
00831     MsqMatrix< 3, 3 > R( r );
00832     TargetCalculator::factor_3D( R, lambda, V, Q, delta, err );
00833     ASSERT_NO_ERROR( err );
00834     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( R ), EPS );
00835 
00836     double w[] = { 9, 8, 7, 1, 5, 4, 3, 2, 6 };
00837     MsqMatrix< 3, 3 > W( w );
00838     TargetCalculator::factor_3D( W, lambda, V, Q, delta, err );
00839     ASSERT_NO_ERROR( err );
00840     ASSERT_MATRICES_EQUAL( Q * delta, TargetCalculator::shape( W ), 2 * EPS );
00841 
00842     W = MsqMatrix< 3, 3 >( 2.1 ) * Q3D_45;
00843     ASSERT_MATRICES_EQUAL( Q3D_45, TargetCalculator::shape( W ), EPS );
00844 }
00845 
00846 void TargetCalculatorTest::test_ideal_skew_tri()
00847 {
00848     MsqError err;
00849     PatchData pd;
00850     MsqMatrix< 2, 2 > Q, Q2;
00851 
00852     TargetCalculator::ideal_skew_2D( TRIANGLE, Sample( 0, 0 ), pd, Q, err );
00853     ASSERT_NO_ERROR( err );
00854     check_valid_Q( Q );
00855 
00856     TargetCalculator::ideal_skew_2D( TRIANGLE, Sample( 2, 0 ), pd, Q2, err );
00857     ASSERT_NO_ERROR( err );
00858     ASSERT_MATRICES_EQUAL( Q, Q2, EPS );
00859 
00860     const Vector3D* coords = unit_edge_element( TRIANGLE );
00861     MsqMatrix< 3, 2 > W;
00862     TargetCalculator::jacobian_2D( pd, TRIANGLE, 3, Sample( 0, 0 ), coords, W, err );
00863     ASSERT_NO_ERROR( err );
00864     Q = TargetCalculator::skew( W );
00865     ASSERT_MATRICES_EQUAL( Q, Q2, EPS );
00866 }
00867 
00868 void TargetCalculatorTest::test_ideal_skew_quad()
00869 {
00870     MsqError err;
00871     PatchData pd;
00872     MsqMatrix< 2, 2 > Q;
00873 
00874     TargetCalculator::ideal_skew_2D( QUADRILATERAL, Sample( 0, 0 ), pd, Q, err );
00875     ASSERT_NO_ERROR( err );
00876     ASSERT_IDENTITY_MATRIX( Q );
00877 
00878     TargetCalculator::ideal_skew_2D( QUADRILATERAL, Sample( 2, 0 ), pd, Q, err );
00879     ASSERT_NO_ERROR( err );
00880     ASSERT_IDENTITY_MATRIX( Q );
00881 }
00882 
00883 void TargetCalculatorTest::test_ideal_skew_tet()
00884 {
00885     MsqError err;
00886     PatchData pd;
00887     MsqMatrix< 3, 3 > Q, Q2;
00888 
00889     TargetCalculator::ideal_skew_3D( TETRAHEDRON, Sample( 0, 0 ), pd, Q, err );
00890     ASSERT_NO_ERROR( err );
00891     check_valid_Q( Q );
00892 
00893     TargetCalculator::ideal_skew_3D( TETRAHEDRON, Sample( 3, 0 ), pd, Q2, err );
00894     ASSERT_NO_ERROR( err );
00895     ASSERT_MATRICES_EQUAL( Q, Q2, EPS );
00896 
00897     const Vector3D* coords = unit_edge_element( TETRAHEDRON );
00898     MsqMatrix< 3, 3 > W;
00899     TargetCalculator::jacobian_3D( pd, TETRAHEDRON, 4, Sample( 0, 0 ), coords, W, err );
00900     ASSERT_NO_ERROR( err );
00901     Q = TargetCalculator::skew( W );
00902     ASSERT_MATRICES_EQUAL( Q, Q2, EPS );
00903 }
00904 
00905 void TargetCalculatorTest::test_ideal_skew_prism()
00906 {
00907     Sample points[]   = { Sample( 0, 0 ), Sample( 0, 1 ), Sample( 0, 2 ),
00908                         Sample( 0, 3 ), Sample( 2, 0 ), Sample( 3, 0 ) };
00909     const int num_pts = sizeof( points ) / sizeof( points[0] );
00910 
00911     MsqError err;
00912     PatchData pd;
00913     MsqMatrix< 3, 3 > Q, W;
00914     const Vector3D* coords = unit_edge_element( PRISM );
00915 
00916     for( int i = 0; i < num_pts; ++i )
00917     {
00918 
00919         TargetCalculator::ideal_skew_3D( PRISM, points[i], pd, Q, err );
00920         ASSERT_NO_ERROR( err );
00921         check_valid_Q( Q );
00922 
00923         TargetCalculator::jacobian_3D( pd, PRISM, 6, points[i], coords, W, err );
00924         ASSERT_NO_ERROR( err );
00925 
00926         ASSERT_MATRICES_EQUAL( TargetCalculator::skew( W ), Q, EPS );
00927     }
00928 }
00929 
00930 void TargetCalculatorTest::test_ideal_skew_pyramid()
00931 {
00932     Sample points[]   = { Sample( 0, 0 ), Sample( 0, 1 ), Sample( 0, 2 ),
00933                         Sample( 0, 3 ), Sample( 2, 0 ), Sample( 3, 0 ) };
00934     const int num_pts = sizeof( points ) / sizeof( points[0] );
00935 
00936     MsqError err;
00937     PatchData pd;
00938     MsqMatrix< 3, 3 > Q, W;
00939     const Vector3D* coords = unit_edge_element( PYRAMID, true );
00940 
00941     for( int i = 0; i < num_pts; ++i )
00942     {
00943 
00944         TargetCalculator::ideal_skew_3D( PYRAMID, points[i], pd, Q, err );
00945         ASSERT_NO_ERROR( err );
00946         check_valid_Q( Q );
00947 
00948         TargetCalculator::jacobian_3D( pd, PYRAMID, 5, points[i], coords, W, err );
00949         ASSERT_NO_ERROR( err );
00950 
00951         ASSERT_MATRICES_EQUAL( TargetCalculator::skew( W ), Q, EPS );
00952     }
00953 }
00954 
00955 void TargetCalculatorTest::test_ideal_skew_hex()
00956 {
00957     MsqError err;
00958     PatchData pd;
00959     MsqMatrix< 3, 3 > Q;
00960 
00961     TargetCalculator::ideal_skew_3D( HEXAHEDRON, Sample( 0, 0 ), pd, Q, err );
00962     ASSERT_NO_ERROR( err );
00963     ASSERT_IDENTITY_MATRIX( Q );
00964 
00965     TargetCalculator::ideal_skew_3D( HEXAHEDRON, Sample( 3, 0 ), pd, Q, err );
00966     ASSERT_NO_ERROR( err );
00967     ASSERT_IDENTITY_MATRIX( Q );
00968 }
00969 
00970 void TargetCalculatorTest::test_ideal_shape_tri()
00971 {
00972     // Ideal triangle should have Aspect (i.e. delta) == identity,
00973     // so shape should be equal to skew.
00974 
00975     MsqError err;
00976     PatchData pd;
00977     MsqMatrix< 2, 2 > Q, Q2;
00978 
00979     TargetCalculator::ideal_shape_2D( TRIANGLE, Sample( 0, 0 ), pd, Q, err );
00980     ASSERT_NO_ERROR( err );
00981     TargetCalculator::ideal_skew_2D( TRIANGLE, Sample( 0, 0 ), pd, Q2, err );
00982     ASSERT_NO_ERROR( err );
00983     ASSERT_MATRICES_EQUAL( Q2, Q, EPS );
00984 
00985     TargetCalculator::ideal_shape_2D( TRIANGLE, Sample( 2, 0 ), pd, Q, err );
00986     ASSERT_NO_ERROR( err );
00987     TargetCalculator::ideal_skew_2D( TRIANGLE, Sample( 2, 0 ), pd, Q2, err );
00988     ASSERT_NO_ERROR( err );
00989     ASSERT_MATRICES_EQUAL( Q2, Q, EPS );
00990 }
00991 
00992 void TargetCalculatorTest::test_ideal_shape_quad()
00993 {
00994     MsqError err;
00995     PatchData pd;
00996     MsqMatrix< 2, 2 > Q;
00997 
00998     TargetCalculator::ideal_shape_2D( QUADRILATERAL, Sample( 0, 0 ), pd, Q, err );
00999     ASSERT_NO_ERROR( err );
01000     ASSERT_IDENTITY_MATRIX( Q );
01001 
01002     TargetCalculator::ideal_shape_2D( QUADRILATERAL, Sample( 2, 0 ), pd, Q, err );
01003     ASSERT_NO_ERROR( err );
01004     ASSERT_IDENTITY_MATRIX( Q );
01005 }
01006 
01007 void TargetCalculatorTest::test_ideal_shape_tet()
01008 {
01009     // Ideal tetrahedron should have Aspect (i.e. delta) == identity,
01010     // so shape should be equal to skew.
01011 
01012     MsqError err;
01013     PatchData pd;
01014     MsqMatrix< 3, 3 > Q, Q2;
01015 
01016     TargetCalculator::ideal_shape_3D( TETRAHEDRON, Sample( 0, 0 ), pd, Q, err );
01017     ASSERT_NO_ERROR( err );
01018     TargetCalculator::ideal_skew_3D( TETRAHEDRON, Sample( 0, 0 ), pd, Q2, err );
01019     ASSERT_NO_ERROR( err );
01020     ASSERT_MATRICES_EQUAL( Q2, Q, EPS );
01021 
01022     TargetCalculator::ideal_shape_3D( TETRAHEDRON, Sample( 3, 0 ), pd, Q, err );
01023     ASSERT_NO_ERROR( err );
01024     TargetCalculator::ideal_skew_3D( TETRAHEDRON, Sample( 3, 0 ), pd, Q2, err );
01025     ASSERT_NO_ERROR( err );
01026     ASSERT_MATRICES_EQUAL( Q2, Q, EPS );
01027 }
01028 
01029 void TargetCalculatorTest::test_ideal_shape_prism()
01030 {
01031     // Ideal wedge should have Aspect (i.e. delta) == identity,
01032     // so shape should be equal to skew.
01033 
01034     Sample points[]   = { Sample( 0, 0 ), Sample( 0, 1 ), Sample( 0, 2 ),
01035                         Sample( 0, 3 ), Sample( 2, 0 ), Sample( 3, 0 ) };
01036     const int num_pts = sizeof( points ) / sizeof( points[0] );
01037 
01038     MsqError err;
01039     PatchData pd;
01040     MsqMatrix< 3, 3 > Q, W;
01041 
01042     for( int i = 0; i < num_pts; ++i )
01043     {
01044         TargetCalculator::ideal_shape_3D( PRISM, points[i], pd, Q, err );
01045         ASSERT_NO_ERROR( err );
01046         TargetCalculator::ideal_skew_3D( PRISM, points[i], pd, W, err );
01047         ASSERT_NO_ERROR( err );
01048         ASSERT_MATRICES_EQUAL( W, Q, EPS );
01049     }
01050 }
01051 
01052 void TargetCalculatorTest::test_ideal_shape_pyramid()
01053 {
01054     Sample points[]   = { Sample( 0, 0 ), Sample( 0, 1 ), Sample( 0, 2 ),
01055                         Sample( 0, 3 ), Sample( 2, 0 ), Sample( 3, 0 ) };
01056     const int num_pts = sizeof( points ) / sizeof( points[0] );
01057 
01058     MsqError err;
01059     PatchData pd;
01060     MsqMatrix< 3, 3 > Q, W;
01061     const Vector3D* coords = unit_edge_element( PYRAMID, true );
01062 
01063     for( int i = 0; i < num_pts; ++i )
01064     {
01065 
01066         TargetCalculator::ideal_shape_3D( PYRAMID, points[i], pd, Q, err );
01067         ASSERT_NO_ERROR( err );
01068 
01069         TargetCalculator::jacobian_3D( pd, PYRAMID, 5, points[i], coords, W, err );
01070         ASSERT_NO_ERROR( err );
01071 
01072         ASSERT_MATRICES_EQUAL( TargetCalculator::shape( W ), Q, EPS );
01073     }
01074 }
01075 
01076 void TargetCalculatorTest::test_ideal_shape_hex()
01077 {
01078     MsqError err;
01079     PatchData pd;
01080     MsqMatrix< 3, 3 > Q;
01081 
01082     TargetCalculator::ideal_shape_3D( HEXAHEDRON, Sample( 0, 0 ), pd, Q, err );
01083     ASSERT_NO_ERROR( err );
01084     ASSERT_IDENTITY_MATRIX( Q );
01085 
01086     TargetCalculator::ideal_shape_3D( HEXAHEDRON, Sample( 3, 0 ), pd, Q, err );
01087     ASSERT_NO_ERROR( err );
01088     ASSERT_IDENTITY_MATRIX( Q );
01089 }
01090 
01091 void TargetCalculatorTest::test_new_orientatin_3D()
01092 {
01093     MsqVector< 3 > x( 0.0 ), y( 0.0 ), z( 0.0 );
01094     x[0] = y[1] = z[2] = 1;
01095 
01096     MsqMatrix< 3, 3 > V = TargetCalculator::new_orientation_3D( x, y );
01097     check_valid_V( V );
01098     ASSERT_MATRICES_EQUAL( z, V.column( 0 ) * V.column( 1 ), EPS );
01099     ASSERT_MATRICES_EQUAL( z, V.column( 2 ), EPS );
01100 }
01101 
01102 void TargetCalculatorTest::test_new_orientatin_2D()
01103 {
01104     MsqVector< 3 > x( 0.0 ), y( 0.0 ), z( 0.0 );
01105     x[0] = y[1] = z[2] = 1;
01106 
01107     MsqMatrix< 3, 2 > V = TargetCalculator::new_orientation_2D( x, y );
01108     check_valid_V( V );
01109     ASSERT_MATRICES_EQUAL( z, V.column( 0 ) * V.column( 1 ), EPS );
01110 }
01111 
01112 void TargetCalculatorTest::test_new_aspect_3D()
01113 {
01114     MsqVector< 3 > r;
01115     r[0]                    = 2;
01116     r[1]                    = 4;
01117     r[2]                    = 6;
01118     MsqMatrix< 3, 3 > delta = TargetCalculator::new_aspect_3D( r );
01119     check_valid_delta( delta );
01120     CPPUNIT_ASSERT_DOUBLES_EQUAL( r[0] / r[1], delta( 0, 0 ) / delta( 1, 1 ), EPS );
01121     CPPUNIT_ASSERT_DOUBLES_EQUAL( r[0] / r[2], delta( 0, 0 ) / delta( 2, 2 ), EPS );
01122 }
01123 
01124 void TargetCalculatorTest::test_new_aspect_2D()
01125 {
01126     MsqVector< 2 > r;
01127     r[0]                    = 7;
01128     r[1]                    = 3;
01129     MsqMatrix< 2, 2 > delta = TargetCalculator::new_aspect_2D( r );
01130     check_valid_delta( delta );
01131     CPPUNIT_ASSERT_DOUBLES_EQUAL( r[0] / r[1], delta( 0, 0 ) / delta( 1, 1 ), EPS );
01132 }
01133 
01134 void TargetCalculatorTest::test_jacobian_3D()
01135 {
01136     MsqError err;
01137     PatchData pd;
01138     LinearHexahedron map;
01139     const Vector3D hex_coords[]  = { Vector3D( 0, 0, 0 ), Vector3D( 2, 0, 0 ), Vector3D( 3, 2, 0 ),
01140                                     Vector3D( 1, 2, 0 ), Vector3D( 1, 0, 1 ), Vector3D( 3, 0, 2 ),
01141                                     Vector3D( 3, 2, 2 ), Vector3D( 2, 2, 2 ) };
01142     const MsqVector< 3 >* coords = reinterpret_cast< const MsqVector< 3 >* >( hex_coords );
01143 
01144     Sample pts[]      = { Sample( 0, 0 ), Sample( 1, 6 ), Sample( 3, 0 ) };
01145     const int num_pts = sizeof( pts ) / sizeof( pts[0] );
01146     for( int i = 0; i < num_pts; ++i )
01147     {
01148         MsqMatrix< 3, 3 > J( 0.0 ), W;
01149         size_t indices[8], n;
01150         MsqVector< 3 > derivs[8];
01151         map.derivatives( pts[i], NodeSet(), indices, derivs, n, err );
01152         ASSERT_NO_ERROR( err );
01153         for( size_t j = 0; j < n; ++j )
01154             J += outer( coords[indices[j]], derivs[j] );
01155 
01156         TargetCalculator::jacobian_3D( pd, HEXAHEDRON, 8, pts[i], hex_coords, W, err );
01157         ASSERT_NO_ERROR( err );
01158 
01159         ASSERT_MATRICES_EQUAL( J, W, EPS );
01160     }
01161 }
01162 
01163 void TargetCalculatorTest::test_jacobian_2D()
01164 {
01165     MsqError err;
01166     PatchData pd;
01167     LinearTriangle map;
01168     const Vector3D tri_coords[]  = { Vector3D( 0, 0, 0 ), Vector3D( 2, 0, 0 ), Vector3D( 0, 0, 1 ) };
01169     const MsqVector< 3 >* coords = reinterpret_cast< const MsqVector< 3 >* >( tri_coords );
01170 
01171     Sample pts[]      = { Sample( 0, 0 ), Sample( 1, 1 ), Sample( 2, 0 ) };
01172     const int num_pts = sizeof( pts ) / sizeof( pts[0] );
01173     for( int i = 0; i < num_pts; ++i )
01174     {
01175         MsqMatrix< 3, 2 > J( 0.0 ), W;
01176         size_t indices[3], n;
01177         MsqVector< 2 > derivs[3];
01178         map.derivatives( pts[i], NodeSet(), indices, derivs, n, err );
01179         ASSERT_NO_ERROR( err );
01180         for( size_t j = 0; j < n; ++j )
01181             J += outer( coords[indices[j]], derivs[j] );
01182 
01183         TargetCalculator::jacobian_2D( pd, TRIANGLE, 3, pts[i], tri_coords, W, err );
01184         ASSERT_NO_ERROR( err );
01185 
01186         ASSERT_MATRICES_EQUAL( J, W, EPS );
01187     }
01188 }
01189 
01190 class DummyRefMesh : public ReferenceMeshInterface
01191 {
01192     std::vector< Vector3D > mCoords;
01193 
01194   public:
01195     DummyRefMesh( const Vector3D* coords, int num_verts ) : mCoords( coords, coords + num_verts ) {}
01196 
01197     void get_reference_vertex_coordinates( const Mesh::VertexHandle* handles,
01198                                            const size_t num_vertices,
01199                                            Vector3D* coords,
01200                                            MsqError& err )
01201     {
01202         const size_t* indices = reinterpret_cast< const size_t* >( handles );
01203         for( size_t i = 0; i < num_vertices; ++i )
01204         {
01205             if( i >= mCoords.size() )
01206             {
01207                 MSQ_SETERR( err )( MsqError::INVALID_ARG );
01208                 return;
01209             }
01210             coords[i] = mCoords[indices[i]];
01211         }
01212     }
01213 };
01214 
01215 void TargetCalculatorTest::test_get_refmesh_Jacobian_3D()
01216 {
01217     MsqError err;
01218 
01219     const Vector3D hex_coords[] = { Vector3D( 0, 0, 0 ), Vector3D( 2, 0, 0 ), Vector3D( 3, 2, 0 ),
01220                                     Vector3D( 1, 2, 0 ), Vector3D( 1, 0, 1 ), Vector3D( 3, 0, 2 ),
01221                                     Vector3D( 3, 2, 2 ), Vector3D( 2, 2, 2 ) };
01222     DummyRefMesh ref_mesh( hex_coords, 8 );
01223 
01224     const Vector3D rect_coords[] = { Vector3D( 0, 0, 0 ), Vector3D( 1, 0, 0 ), Vector3D( 1, 1, 0 ),
01225                                      Vector3D( 0, 1, 0 ), Vector3D( 0, 0, 5 ), Vector3D( 1, 0, 5 ),
01226                                      Vector3D( 1, 1, 5 ), Vector3D( 0, 1, 5 ) };
01227     size_t conn[]                = { 0, 1, 2, 3, 4, 5, 6, 7 };
01228     bool fixed[8];
01229     std::fill( fixed, fixed + sizeof( fixed ) / sizeof( fixed[0] ), false );
01230     PatchData pd;
01231     pd.fill( 8, rect_coords[0].to_array(), 1, HEXAHEDRON, conn, fixed, err );
01232     ASSERT_NO_ERROR( err );
01233 
01234     Sample pts[]      = { Sample( 0, 0 ), Sample( 1, 6 ), Sample( 3, 0 ) };
01235     const int num_pts = sizeof( pts ) / sizeof( pts[0] );
01236     for( int i = 0; i < num_pts; ++i )
01237     {
01238         MsqMatrix< 3, 3 > J, W;
01239         TargetCalculator::jacobian_3D( pd, HEXAHEDRON, 8, pts[i], hex_coords, J, err );
01240         ASSERT_NO_ERROR( err );
01241 
01242         TargetCalculator::get_refmesh_Jacobian_3D( &ref_mesh, pd, 0, pts[i], W, err );
01243         ASSERT_NO_ERROR( err );
01244 
01245         ASSERT_MATRICES_EQUAL( J, W, EPS );
01246     }
01247 }
01248 
01249 void TargetCalculatorTest::test_get_refmesh_Jacobian_2D()
01250 {
01251     MsqError err;
01252 
01253     const Vector3D tri_coords[] = { Vector3D( 0, 0, 0 ), Vector3D( 3, 0, 1 ), Vector3D( 1, 0, 3 ) };
01254     DummyRefMesh ref_mesh( tri_coords, 3 );
01255 
01256     const Vector3D right_coords[] = { Vector3D( 0, 0, 0 ), Vector3D( 5, 0, 0 ), Vector3D( 5, 1, 0 ) };
01257     size_t conn[]                 = { 0, 1, 2 };
01258     bool fixed[3];
01259     std::fill( fixed, fixed + sizeof( fixed ) / sizeof( fixed[0] ), false );
01260     PatchData pd;
01261     pd.fill( 3, right_coords[0].to_array(), 1, TRIANGLE, conn, fixed, err );
01262     ASSERT_NO_ERROR( err );
01263 
01264     Sample pts[]      = { Sample( 0, 0 ), Sample( 1, 1 ), Sample( 2, 0 ) };
01265     const int num_pts = sizeof( pts ) / sizeof( pts[0] );
01266     for( int i = 0; i < num_pts; ++i )
01267     {
01268         MsqMatrix< 3, 2 > J, W;
01269         TargetCalculator::jacobian_2D( pd, TRIANGLE, 3, pts[i], tri_coords, J, err );
01270         ASSERT_NO_ERROR( err );
01271 
01272         TargetCalculator::get_refmesh_Jacobian_2D( &ref_mesh, pd, 0, pts[i], W, err );
01273         ASSERT_NO_ERROR( err );
01274 
01275         ASSERT_MATRICES_EQUAL( J, W, EPS );
01276     }
01277 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines