MOAB: Mesh Oriented datABase  (version 5.3.1)
SlaveBoundaryVerticesTest.cpp
Go to the documentation of this file.
00001 /* *****************************************************************
00002     MESQUITE -- The Mesh Quality Improvement Toolkit
00003 
00004     Copyright 2007 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 SlaveBoundaryVerticesTest.cpp
00028  *  \brief Test SlaveBoundaryVertices class
00029  *  \author Jason Kraftcheck
00030  */
00031 
00032 #include "Mesquite.hpp"
00033 #include "UnitUtil.hpp"
00034 #include "meshfiles.h"
00035 #include "Settings.hpp"
00036 #include "MsqError.hpp"
00037 #include "SlaveBoundaryVertices.hpp"
00038 
00039 #include "MsqVertex.hpp"
00040 #include "DomainClassifier.hpp"
00041 #include "PlanarDomain.hpp"
00042 #include "MeshDomain1D.hpp"
00043 #include "MeshImpl.hpp"
00044 
00045 #include <vector>
00046 #include <set>
00047 #include <algorithm>
00048 #include <fstream>
00049 #include <iostream>
00050 
00051 using namespace MBMesquite;
00052 
00053 class SlaveBoundaryVerticesTest : public CppUnit::TestFixture
00054 {
00055   private:
00056     CPPUNIT_TEST_SUITE( SlaveBoundaryVerticesTest );
00057     CPPUNIT_TEST( test_fail_if_slaves_not_calculated );
00058     CPPUNIT_TEST( test_one_depth_from_fixed );
00059     CPPUNIT_TEST( test_two_depth_from_fixed );
00060     CPPUNIT_TEST( test_zero_depth_from_surface );
00061     CPPUNIT_TEST( test_one_depth_from_surface );
00062     CPPUNIT_TEST( test_two_depth_from_surface );
00063     CPPUNIT_TEST( test_zero_depth_from_curve );
00064     CPPUNIT_TEST( test_one_depth_from_curve );
00065     CPPUNIT_TEST( test_two_depth_from_curve );
00066     CPPUNIT_TEST_SUITE_END();
00067 
00068     void test_slaved_common( unsigned depth, unsigned boundary );
00069 
00070     void make_mesh( MeshImpl& mesh, DomainClassifier& domain, const int intervals );
00071 
00072   public:
00073     void test_fail_if_slaves_not_calculated();
00074     void test_one_depth_from_fixed()
00075     {
00076         test_slaved_common( 1, 4 );
00077     }
00078     void test_two_depth_from_fixed()
00079     {
00080         test_slaved_common( 2, 4 );
00081     }
00082     void test_zero_depth_from_surface()
00083     {
00084         test_slaved_common( 0, 2 );
00085     }
00086     void test_one_depth_from_surface()
00087     {
00088         test_slaved_common( 1, 2 );
00089     }
00090     void test_two_depth_from_surface()
00091     {
00092         test_slaved_common( 2, 2 );
00093     }
00094     void test_zero_depth_from_curve()
00095     {
00096         test_slaved_common( 0, 1 );
00097     }
00098     void test_one_depth_from_curve()
00099     {
00100         test_slaved_common( 1, 1 );
00101     }
00102     void test_two_depth_from_curve()
00103     {
00104         test_slaved_common( 2, 1 );
00105     }
00106 };
00107 
00108 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION( SlaveBoundaryVerticesTest, "SlaveBoundaryVerticesTest" );
00109 CPPUNIT_TEST_SUITE_NAMED_REGISTRATION( SlaveBoundaryVerticesTest, "Unit" );
00110 
00111 void SlaveBoundaryVerticesTest::test_fail_if_slaves_not_calculated()
00112 {
00113     MsqError err;
00114 
00115     Settings settings;
00116     settings.set_slaved_ho_node_mode( Settings::SLAVE_ALL );
00117     SlaveBoundaryVertices tool( 1 );
00118 
00119     MeshImpl mesh;
00120     DomainClassifier domain;
00121     make_mesh( mesh, domain, 2 );
00122 
00123     MeshDomainAssoc mesh_and_domain = MeshDomainAssoc( &mesh, &domain );
00124     tool.loop_over_mesh( &mesh_and_domain, &settings, err );
00125     CPPUNIT_ASSERT( err );
00126     err.clear();
00127 }
00128 
00129 void SlaveBoundaryVerticesTest::make_mesh( MeshImpl& mesh, DomainClassifier& domain, const int intervals )
00130 {
00131     MsqPrintError err( std::cerr );
00132     const char input_file[] = MESH_FILES_DIR "3D/vtk/quadratic/6x6x6-hex20.vtk";
00133 
00134     const Vector3D min( -3, -3, -3 );
00135     const Vector3D max( 3, 3, 3 );
00136     const Vector3D space = ( max - min ) * 1.0 / ( intervals );
00137 
00138     mesh.clear();
00139     mesh.read_vtk( input_file, err );
00140     ASSERT_NO_ERROR( err );
00141 
00142     const Vector3D corners[8] = { Vector3D( min[0], min[1], min[2] ), Vector3D( max[0], min[1], min[2] ),
00143                                   Vector3D( max[0], max[1], min[2] ), Vector3D( min[0], max[1], min[2] ),
00144                                   Vector3D( min[0], min[1], max[2] ), Vector3D( max[0], min[1], max[2] ),
00145                                   Vector3D( max[0], max[1], max[2] ), Vector3D( min[0], max[1], max[2] ) };
00146 
00147     MeshDomain* subdomains[26] = { new PlanarDomain( PlanarDomain::XZ, min[1] ),
00148                                    new PlanarDomain( PlanarDomain::YZ, max[0] ),
00149                                    new PlanarDomain( PlanarDomain::XZ, max[1] ),
00150                                    new PlanarDomain( PlanarDomain::YZ, min[0] ),
00151                                    new PlanarDomain( PlanarDomain::XY, min[2] ),
00152                                    new PlanarDomain( PlanarDomain::XY, max[2] ),
00153                                    new LineDomain( corners[0], corners[1] - corners[0] ),
00154                                    new LineDomain( corners[1], corners[2] - corners[1] ),
00155                                    new LineDomain( corners[2], corners[3] - corners[2] ),
00156                                    new LineDomain( corners[3], corners[0] - corners[3] ),
00157                                    new LineDomain( corners[0], corners[4] - corners[0] ),
00158                                    new LineDomain( corners[1], corners[5] - corners[1] ),
00159                                    new LineDomain( corners[2], corners[6] - corners[0] ),
00160                                    new LineDomain( corners[3], corners[7] - corners[1] ),
00161                                    new LineDomain( corners[4], corners[5] - corners[4] ),
00162                                    new LineDomain( corners[5], corners[6] - corners[5] ),
00163                                    new LineDomain( corners[6], corners[7] - corners[6] ),
00164                                    new LineDomain( corners[7], corners[4] - corners[7] ),
00165                                    new PointDomain( corners[0] ),
00166                                    new PointDomain( corners[1] ),
00167                                    new PointDomain( corners[2] ),
00168                                    new PointDomain( corners[3] ),
00169                                    new PointDomain( corners[4] ),
00170                                    new PointDomain( corners[5] ),
00171                                    new PointDomain( corners[6] ),
00172                                    new PointDomain( corners[7] ) };
00173     const int subdims[26]      = { 2, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0 };
00174     DomainClassifier::classify_skin_geometrically( domain, &mesh, 1e-6, subdomains, subdims, 26, err );
00175     domain.delete_sub_domains( true );
00176     ASSERT_NO_ERROR( err );
00177 }
00178 
00179 void SlaveBoundaryVerticesTest::test_slaved_common( unsigned depth, unsigned boundary )
00180 {
00181     MeshImpl mesh;
00182     DomainClassifier domain;
00183     make_mesh( mesh, domain, 2 * depth + 2 );
00184 
00185     MsqPrintError err( std::cerr );
00186     std::vector< std::vector< Mesh::VertexHandle > > depths( depth + 1 );
00187     std::set< Mesh::VertexHandle > non_slave;
00188     std::set< Mesh::VertexHandle >::iterator p;
00189 
00190     // find boundary vertices
00191     std::vector< Mesh::VertexHandle > verts;
00192     mesh.get_all_vertices( verts, err );
00193     ASSERT_NO_ERROR( err );
00194     CPPUNIT_ASSERT( !verts.empty() );
00195     if( boundary >= 4 )
00196     {
00197         std::vector< bool > flags;
00198         mesh.vertices_get_fixed_flag( arrptr( verts ), flags, verts.size(), err );
00199         ASSERT_NO_ERROR( err );
00200         for( size_t i = 0; i < verts.size(); ++i )
00201             if( flags[i] )
00202             {
00203                 depths[0].push_back( verts[i] );
00204                 non_slave.insert( verts[i] );
00205             }
00206     }
00207     else
00208     {
00209         std::vector< unsigned short > dim( verts.size() );
00210         domain.domain_DoF( arrptr( verts ), arrptr( dim ), verts.size(), err );
00211         ASSERT_NO_ERROR( err );
00212         for( size_t i = 0; i < verts.size(); ++i )
00213             if( dim[i] <= boundary )
00214             {
00215                 depths[0].push_back( verts[i] );
00216                 non_slave.insert( verts[i] );
00217             }
00218     }
00219 
00220     // check that our input is usable for this test
00221     CPPUNIT_ASSERT( !verts.empty() );
00222 
00223     // find all vertices up to specified depth
00224     for( unsigned d = 0; d < depth; ++d )
00225     {
00226         for( size_t i = 0; i < depths[d].size(); ++i )
00227         {
00228             std::vector< Mesh::ElementHandle > adj;
00229             std::vector< size_t > junk;
00230             mesh.vertices_get_attached_elements( &depths[d][i], 1, adj, junk, err );
00231             ASSERT_NO_ERROR( err );
00232             for( size_t j = 0; j < adj.size(); ++j )
00233             {
00234                 junk.clear();
00235                 std::vector< Mesh::VertexHandle > conn;
00236                 mesh.elements_get_attached_vertices( &adj[j], 1, conn, junk, err );
00237                 ASSERT_NO_ERROR( err );
00238                 for( size_t k = 0; k < conn.size(); ++k )
00239                 {
00240                     p = non_slave.find( conn[k] );
00241                     if( p == non_slave.end() )
00242                     {
00243                         non_slave.insert( p, conn[k] );
00244                         depths[d + 1].push_back( conn[k] );
00245                     }
00246                 }
00247             }
00248         }
00249     }
00250 
00251     // Check that our input is usable for this test:
00252     // Should have some vertices that are not within the specified depth of
00253     // the boundary.
00254     CPPUNIT_ASSERT( non_slave.size() < verts.size() );
00255 
00256     // Now build a map of all higher-order nodes in the mesh
00257     std::set< Mesh::VertexHandle > higher_order;
00258     std::vector< Mesh::ElementHandle > elems;
00259     mesh.get_all_elements( elems, err );
00260     ASSERT_NO_ERROR( err );
00261     CPPUNIT_ASSERT( !elems.empty() );
00262     std::vector< EntityTopology > types( elems.size() );
00263     mesh.elements_get_topologies( arrptr( elems ), arrptr( types ), elems.size(), err );
00264     ASSERT_NO_ERROR( err );
00265     for( size_t i = 0; i < elems.size(); ++i )
00266     {
00267         std::vector< Mesh::VertexHandle > conn;
00268         std::vector< size_t > junk;
00269         mesh.elements_get_attached_vertices( &elems[i], 1, conn, junk, err );
00270         ASSERT_NO_ERROR( err );
00271         for( size_t j = TopologyInfo::corners( types[i] ); j < conn.size(); ++j )
00272             higher_order.insert( conn[j] );
00273     }
00274 
00275     // Check that our input is usable for this test:
00276     // Should have some higher-order vertices
00277     CPPUNIT_ASSERT( !higher_order.empty() );
00278 
00279     // Now build a map of all fixed vertices
00280     std::set< Mesh::VertexHandle > fixed_vertices;
00281     std::vector< bool > fixed;
00282     mesh.vertices_get_fixed_flag( arrptr( verts ), fixed, verts.size(), err );
00283     ASSERT_NO_ERROR( err );
00284     for( size_t i = 0; i < verts.size(); ++i )
00285         if( fixed[i] ) fixed_vertices.insert( verts[i] );
00286 
00287     // Now actually run the tool
00288     Settings settings;
00289     settings.set_slaved_ho_node_mode( Settings::SLAVE_CALCULATED );
00290     SlaveBoundaryVertices tool( depth, boundary );
00291     MeshDomainAssoc mesh_and_domain = MeshDomainAssoc( &mesh, &domain );
00292     tool.loop_over_mesh( &mesh_and_domain, &settings, err );
00293     ASSERT_NO_ERROR( err );
00294 
00295     // Now verify the results
00296     std::vector< unsigned char > bytes( verts.size() );
00297     mesh.vertices_get_byte( arrptr( verts ), arrptr( bytes ), verts.size(), err );
00298     ASSERT_NO_ERROR( err );
00299     for( size_t i = 0; i < verts.size(); ++i )
00300     {
00301         bool in_non_slave    = ( non_slave.find( verts[i] ) != non_slave.end() );
00302         bool in_fixed        = ( fixed_vertices.find( verts[i] ) != fixed_vertices.end() );
00303         bool in_higher_order = ( higher_order.find( verts[i] ) != higher_order.end() );
00304         if( bytes[i] & MsqVertex::MSQ_DEPENDENT )
00305         {  // if slave node
00306             // must not be within 'depth' of boundary
00307             CPPUNIT_ASSERT( !in_non_slave );
00308             // must be a higher-order vertex
00309             CPPUNIT_ASSERT( in_higher_order );
00310             // must not be fixed
00311             CPPUNIT_ASSERT( !in_fixed );
00312         }
00313         else
00314         {
00315             // there are three reasons that a vertex isn't slaved
00316             bool in_non_slave    = ( non_slave.find( verts[i] ) != non_slave.end() );
00317             bool in_fixed        = ( fixed_vertices.find( verts[i] ) != fixed_vertices.end() );
00318             bool in_higher_order = ( higher_order.find( verts[i] ) != higher_order.end() );
00319             CPPUNIT_ASSERT( in_fixed || !in_higher_order || in_non_slave );
00320         }
00321     }
00322 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines