MOAB: Mesh Oriented datABase  (version 5.4.0)
IdealTargetTest.cpp
Go to the documentation of this file.
00001 /* *****************************************************************
00002     MESQUITE -- The Mesh Quality Improvement Toolkit
00003 
00004     Copyright 2006 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     retian 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     (2006) kraftche@cae.wisc.edu
00024 
00025   ***************************************************************** */
00026 
00027 /** \file IdealTargetTest.cpp
00028  *  \brief Test the IdealShapeTarget
00029  *  \author Jason Kraftcheck
00030  */
00031 
00032 #include "Mesquite.hpp"
00033 #include "IdealShapeTarget.hpp"
00034 #include "MsqError.hpp"
00035 #include "Vector3D.hpp"
00036 #include "PatchData.hpp"
00037 #include "PlanarDomain.hpp"
00038 #include "MappingFunction.hpp"
00039 #include "Settings.hpp"
00040 #include "IdealElements.hpp"
00041 #include "ElemSampleQM.hpp"
00042 #include "cppunit/extensions/HelperMacros.h"
00043 #include "UnitUtil.hpp"
00044 
00045 using namespace MBMesquite;
00046 
00047 class IdealTargetTest : public CppUnit::TestFixture
00048 {
00049   private:
00050     CPPUNIT_TEST_SUITE( IdealTargetTest );
00051     CPPUNIT_TEST( test_tri_corner );
00052     CPPUNIT_TEST( test_tri_edge );
00053     CPPUNIT_TEST( test_tri_center );
00054     CPPUNIT_TEST( test_hex_corner );
00055     CPPUNIT_TEST( test_hex_edge );
00056     CPPUNIT_TEST( test_hex_face );
00057     CPPUNIT_TEST( test_hex_center );
00058     CPPUNIT_TEST_SUITE_END();
00059 
00060   public:
00061     void test_tri_corner();
00062     void test_tri_edge();
00063     void test_tri_center();
00064     void test_hex_corner();
00065     void test_hex_edge();
00066     void test_hex_face();
00067     void test_hex_center();
00068 
00069   private:
00070     void get_calc_target( EntityTopology type, Sample sample, MsqMatrix< 3, 3 >&, MsqMatrix< 2, 2 >& );
00071 
00072     void get_ideal_target( EntityTopology type, Sample sample, MsqMatrix< 3, 3 >&, MsqMatrix< 2, 2 >& );
00073 
00074     void do_test( EntityTopology type, Sample location );
00075 
00076     Settings settings;
00077 };
00078 
00079 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION( IdealTargetTest, "IdealTargetTest" );
00080 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION( IdealTargetTest, "Unit" );
00081 
00082 void IdealTargetTest::test_tri_corner()
00083 {
00084     do_test( TRIANGLE, Sample( 0, 0 ) );
00085     do_test( TRIANGLE, Sample( 0, 1 ) );
00086     do_test( TRIANGLE, Sample( 0, 2 ) );
00087 }
00088 
00089 void IdealTargetTest::test_tri_edge()
00090 {
00091     do_test( TRIANGLE, Sample( 1, 0 ) );
00092     do_test( TRIANGLE, Sample( 1, 1 ) );
00093     do_test( TRIANGLE, Sample( 1, 2 ) );
00094 }
00095 
00096 void IdealTargetTest::test_tri_center()
00097 {
00098     do_test( TRIANGLE, Sample( 2, 0 ) );
00099 }
00100 
00101 void IdealTargetTest::test_hex_corner()
00102 {
00103     do_test( HEXAHEDRON, Sample( 0, 0 ) );
00104     do_test( HEXAHEDRON, Sample( 0, 1 ) );
00105     do_test( HEXAHEDRON, Sample( 0, 6 ) );
00106     do_test( HEXAHEDRON, Sample( 0, 7 ) );
00107 }
00108 
00109 void IdealTargetTest::test_hex_edge()
00110 {
00111     do_test( HEXAHEDRON, Sample( 1, 1 ) );
00112     do_test( HEXAHEDRON, Sample( 1, 4 ) );
00113     do_test( HEXAHEDRON, Sample( 1, 11 ) );
00114 }
00115 
00116 void IdealTargetTest::test_hex_face()
00117 {
00118     do_test( HEXAHEDRON, Sample( 2, 0 ) );
00119     do_test( HEXAHEDRON, Sample( 2, 1 ) );
00120     do_test( HEXAHEDRON, Sample( 2, 2 ) );
00121     do_test( HEXAHEDRON, Sample( 2, 3 ) );
00122     do_test( HEXAHEDRON, Sample( 2, 4 ) );
00123     do_test( HEXAHEDRON, Sample( 2, 5 ) );
00124 }
00125 
00126 void IdealTargetTest::test_hex_center()
00127 {
00128     do_test( HEXAHEDRON, Sample( 3, 0 ) );
00129 }
00130 
00131 void IdealTargetTest::get_calc_target( EntityTopology type,
00132                                        Sample location,
00133                                        MsqMatrix< 3, 3 >& w3,
00134                                        MsqMatrix< 2, 2 >& w2 )
00135 {
00136     MsqPrintError err( std::cout );
00137     const int elem_dim = TopologyInfo::dimension( type );
00138 
00139     // create a patch -- actual coords and such don't really matter
00140     std::vector< double > coords( 24, 0.0 );
00141     const size_t conn[] = { 0, 1, 2, 3, 4, 5, 6, 7 };
00142     PatchData pd;
00143     pd.fill( 8, arrptr( coords ), 1, type, conn, 0, err );
00144     CPPUNIT_ASSERT( !MSQ_CHKERR( err ) );
00145     pd.attach_settings( &settings );
00146 
00147     IdealShapeTarget tc;
00148     if( elem_dim == 2 )
00149         tc.get_2D_target( pd, 0, location, w2, err );
00150     else
00151         tc.get_3D_target( pd, 0, location, w3, err );
00152     CPPUNIT_ASSERT( !MSQ_CHKERR( err ) );
00153 }
00154 
00155 void IdealTargetTest::get_ideal_target( EntityTopology type,
00156                                         Sample location,
00157                                         MsqMatrix< 3, 3 >& w3,
00158                                         MsqMatrix< 2, 2 >& w2 )
00159 {
00160     MsqPrintError err( std::cout );
00161     const unsigned elem_dim = TopologyInfo::dimension( type );
00162 
00163     // get the target matrix for an ideal element
00164     size_t indices[100];
00165     size_t num_vtx;
00166     const Vector3D* coords = unit_element( type );
00167     Vector3D c[3];
00168     if( elem_dim == 2 )
00169     {
00170         MsqVector< 2 > derivs[100];
00171         const MappingFunction2D* func = settings.get_mapping_function_2D( type );
00172         func->derivatives( location, NodeSet(), indices, derivs, num_vtx, err );
00173         CPPUNIT_ASSERT( !MSQ_CHKERR( err ) );
00174 
00175         MsqMatrix< 3, 2 > J;
00176         for( size_t i = 0; i < num_vtx; ++i )
00177             for( unsigned j = 0; j < 2; ++j )
00178                 c[j] += derivs[i][j] * coords[indices[i]];
00179 
00180         for( unsigned i = 0; i < 3; ++i )
00181             for( unsigned j = 0; j < 2; ++j )
00182                 J( i, j ) = c[j][i];
00183 
00184         w2 = TargetCalculator::skew( J );
00185     }
00186     else
00187     {
00188         MsqVector< 3 > derivs[100];
00189         const MappingFunction3D* func = settings.get_mapping_function_3D( type );
00190         func->derivatives( location, NodeSet(), indices, derivs, num_vtx, err );
00191         CPPUNIT_ASSERT( !MSQ_CHKERR( err ) );
00192 
00193         for( size_t i = 0; i < num_vtx; ++i )
00194             for( unsigned j = 0; j < 3; ++j )
00195                 c[j] += derivs[i][j] * coords[indices[i]];
00196 
00197         for( unsigned i = 0; i < 3; ++i )
00198             for( unsigned j = 0; j < 3; ++j )
00199                 w3( i, j ) = c[j][i];
00200 
00201         w3 = TargetCalculator::skew( w3 );
00202     }
00203 }
00204 
00205 void IdealTargetTest::do_test( EntityTopology type, Sample location )
00206 {
00207     MsqMatrix< 3, 3 > w3_calc, w3_exp;
00208     MsqMatrix< 2, 2 > w2_calc, w2_exp;
00209     get_calc_target( type, location, w3_calc, w2_calc );
00210     get_ideal_target( type, location, w3_exp, w2_exp );
00211     if( TopologyInfo::dimension( type ) == 2 )
00212     {
00213         CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, det( w2_calc ), 1e-6 );
00214         CPPUNIT_ASSERT_DOUBLES_EQUAL( 1.0, det( w2_exp ), 1e-6 );
00215         // a rotation of the expected matrix is acceptable.
00216         MsqMatrix< 2, 2 > R = inverse( w2_calc ) * w2_exp;
00217         ASSERT_MATRICES_EQUAL( transpose( R ), inverse( R ), 1e-6 );
00218     }
00219     else
00220     {
00221         MsqMatrix< 3, 3 > R = inverse( w3_calc ) * w3_exp;
00222         ASSERT_MATRICES_EQUAL( transpose( R ), inverse( R ), 1e-6 );
00223     }
00224 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines