UnitTestBulkData_Generate.cpp

00001 /*------------------------------------------------------------------------*/
00002 /*                 Copyright 2010 Sandia Corporation.                     */
00003 /*  Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive   */
00004 /*  license for use of this work by or on behalf of the U.S. Government.  */
00005 /*  Export of this program may require a license from the                 */
00006 /*  United States Government.                                             */
00007 /*------------------------------------------------------------------------*/
00008 
00009 
00010 #include <sstream>
00011 
00012 #include <unit_tests/stk_utest_macros.hpp>
00013 
00014 #include <stk_util/parallel/Parallel.hpp>
00015 #include <stk_mesh/base/BulkData.hpp>
00016 #include <stk_mesh/base/GetEntities.hpp>
00017 #include <stk_mesh/base/Comm.hpp>
00018 #include <stk_mesh/base/EntityComm.hpp>
00019 #include <stk_mesh/fem/EntityTypes.hpp>
00020 
00021 #include <unit_tests/UnitTestBulkData.hpp>
00022 
00023 //----------------------------------------------------------------------
00024 //----------------------------------------------------------------------
00025 
00026 namespace stk {
00027 namespace mesh {
00028 
00029 namespace {
00030 
00031 /* Recursively split a box into ( up - ip ) sub-boxes */
00032 
00033 typedef int BOX[3][2] ;
00034 
00035 void box_partition( int ip , int up , int axis , 
00036                     const BOX box ,
00037                     BOX p_box[] )
00038 { 
00039   const int np = up - ip ;
00040   if ( 1 == np ) {
00041     p_box[ip][0][0] = box[0][0] ; p_box[ip][0][1] = box[0][1] ;
00042     p_box[ip][1][0] = box[1][0] ; p_box[ip][1][1] = box[1][1] ;
00043     p_box[ip][2][0] = box[2][0] ; p_box[ip][2][1] = box[2][1] ;
00044   } 
00045   else {  
00046     const int n = box[ axis ][1] - box[ axis ][0] ;
00047     const int np_low = np / 2 ;  /* Rounded down */
00048     const int np_upp = np - np_low ;
00049  
00050     const int n_upp = (int) (((double) n) * ( ((double)np_upp) / ((double)np)));
00051     const int n_low = n - n_upp ;
00052     const int next_axis = ( axis + 2 ) % 3 ;
00053  
00054     if ( np_low ) { /* P = [ip,ip+np_low) */
00055       BOX dbox ;
00056       dbox[0][0] = box[0][0] ; dbox[0][1] = box[0][1] ;
00057       dbox[1][0] = box[1][0] ; dbox[1][1] = box[1][1] ;
00058       dbox[2][0] = box[2][0] ; dbox[2][1] = box[2][1] ;
00059  
00060       dbox[ axis ][1] = dbox[ axis ][0] + n_low ;
00061  
00062       box_partition( ip, ip + np_low, next_axis,
00063                      (const int (*)[2]) dbox, p_box );
00064     }
00065  
00066     if ( np_upp ) { /* P = [ip+np_low,ip+np_low+np_upp) */
00067       BOX dbox ;
00068       dbox[0][0] = box[0][0] ; dbox[0][1] = box[0][1] ;
00069       dbox[1][0] = box[1][0] ; dbox[1][1] = box[1][1] ;
00070       dbox[2][0] = box[2][0] ; dbox[2][1] = box[2][1] ;
00071  
00072       ip += np_low ;
00073       dbox[ axis ][0] += n_low ;
00074       dbox[ axis ][1]  = dbox[ axis ][0] + n_upp ;
00075  
00076       box_partition( ip, ip + np_upp, next_axis,
00077                      (const int (*)[2]) dbox, p_box );
00078     }
00079   }
00080 }
00081 
00082 }
00083 
00084 void UnitTestBulkData::generate_boxes(
00085   BulkData  & mesh ,
00086   const bool  generate_aura ,
00087   const int   root_box[][2] ,
00088         int   local_box[][2] )
00089 {
00090   const unsigned p_rank = mesh.parallel_rank();
00091   const unsigned p_size = mesh.parallel_size();
00092   const unsigned ngx = root_box[0][1] - root_box[0][0] ;
00093   const unsigned ngy = root_box[1][1] - root_box[1][0] ;
00094   const unsigned ngz = root_box[2][1] - root_box[2][0] ;
00095 /*
00096   const unsigned e_global = ngx * ngy * ngz ;
00097   const unsigned n_global = ( ngx + 1 ) * ( ngy + 1 ) * ( ngz + 1 );
00098 */
00099 
00100   if ( 0 == p_rank ) {
00101     std::cout << "Global box = " << ngx << " x " << ngy << " x " << ngz 
00102               << std::endl ;
00103   }
00104 
00105   BOX * const p_box = new BOX[ p_size ];
00106 
00107   box_partition( 0 , p_size , 2 , root_box , & p_box[0] );
00108 
00109   local_box[0][0] = p_box[ p_rank ][0][0] ;
00110   local_box[0][1] = p_box[ p_rank ][0][1] ;
00111   local_box[1][0] = p_box[ p_rank ][1][0] ;
00112   local_box[1][1] = p_box[ p_rank ][1][1] ;
00113   local_box[2][0] = p_box[ p_rank ][2][0] ;
00114   local_box[2][1] = p_box[ p_rank ][2][1] ;
00115 
00116   const unsigned nx = local_box[0][1] - local_box[0][0] ;
00117   const unsigned ny = local_box[1][1] - local_box[1][0] ;
00118   const unsigned nz = local_box[2][1] - local_box[2][0] ;
00119 
00120   const unsigned e_local = nx * ny * nz ;
00121   const unsigned n_local = ( nx + 1 ) * ( ny + 1 ) * ( nz + 1 );
00122 
00123   // Create elements:
00124 
00125   std::vector<unsigned> local_count ;
00126 
00127   const PartVector no_parts ;
00128 
00129   for ( int k = local_box[2][0] ; k < local_box[2][1] ; ++k ) {
00130   for ( int j = local_box[1][0] ; j < local_box[1][1] ; ++j ) {
00131   for ( int i = local_box[0][0] ; i < local_box[0][1] ; ++i ) {
00132     const EntityId n0 = 1 + (i+0) + (j+0) * (ngx+1) + (k+0) * (ngx+1) * (ngy+1);
00133     const EntityId n1 = 1 + (i+1) + (j+0) * (ngx+1) + (k+0) * (ngx+1) * (ngy+1);
00134     const EntityId n2 = 1 + (i+1) + (j+1) * (ngx+1) + (k+0) * (ngx+1) * (ngy+1);
00135     const EntityId n3 = 1 + (i+0) + (j+1) * (ngx+1) + (k+0) * (ngx+1) * (ngy+1);
00136     const EntityId n4 = 1 + (i+0) + (j+0) * (ngx+1) + (k+1) * (ngx+1) * (ngy+1);
00137     const EntityId n5 = 1 + (i+1) + (j+0) * (ngx+1) + (k+1) * (ngx+1) * (ngy+1);
00138     const EntityId n6 = 1 + (i+1) + (j+1) * (ngx+1) + (k+1) * (ngx+1) * (ngy+1);
00139     const EntityId n7 = 1 + (i+0) + (j+1) * (ngx+1) + (k+1) * (ngx+1) * (ngy+1);
00140 
00141     const EntityId elem_id =  1 + i + j * ngx + k * ngx * ngy;
00142 
00143     Entity & node0 = mesh.declare_entity( 0 , n0 , no_parts );
00144     Entity & node1 = mesh.declare_entity( 0 , n1 , no_parts );
00145     Entity & node2 = mesh.declare_entity( 0 , n2 , no_parts );
00146     Entity & node3 = mesh.declare_entity( 0 , n3 , no_parts );
00147     Entity & node4 = mesh.declare_entity( 0 , n4 , no_parts );
00148     Entity & node5 = mesh.declare_entity( 0 , n5 , no_parts );
00149     Entity & node6 = mesh.declare_entity( 0 , n6 , no_parts );
00150     Entity & node7 = mesh.declare_entity( 0 , n7 , no_parts );
00151     Entity & elem  = mesh.declare_entity( 3 , elem_id , no_parts );
00152 
00153     mesh.declare_relation( elem , node0 , 0 );
00154     mesh.declare_relation( elem , node1 , 1 );
00155     mesh.declare_relation( elem , node2 , 2 );
00156     mesh.declare_relation( elem , node3 , 3 );
00157     mesh.declare_relation( elem , node4 , 4 );
00158     mesh.declare_relation( elem , node5 , 5 );
00159     mesh.declare_relation( elem , node6 , 6 );
00160     mesh.declare_relation( elem , node7 , 7 );
00161   }
00162   }
00163   }
00164 
00165   Selector select_owned( mesh.mesh_meta_data().locally_owned_part() );
00166   Selector select_used( mesh.mesh_meta_data().locally_used_part() );
00167   Selector select_all(  mesh.mesh_meta_data().universal_part() );
00168 
00169   count_entities( select_used , mesh , local_count );
00170   STKUNIT_ASSERT_EQUAL( e_local , local_count[3] );
00171   STKUNIT_ASSERT_EQUAL( 0u , local_count[2] );
00172   STKUNIT_ASSERT_EQUAL( 0u , local_count[1] );
00173   STKUNIT_ASSERT_EQUAL( n_local , local_count[0] );
00174 
00175   // Set up sharing:
00176   STKUNIT_ASSERT( mesh.internal_modification_end( generate_aura ) );
00177 
00178   // Verify declarations and sharing
00179 
00180   count_entities( select_used , mesh , local_count );
00181   STKUNIT_ASSERT( local_count[3] == e_local );
00182   STKUNIT_ASSERT( local_count[2] == 0 );
00183   STKUNIT_ASSERT( local_count[1] == 0 );
00184   STKUNIT_ASSERT( local_count[0] == n_local );
00185 
00186   for ( int k = local_box[2][0] ; k <= local_box[2][1] ; ++k ) {
00187   for ( int j = local_box[1][0] ; j <= local_box[1][1] ; ++j ) {
00188   for ( int i = local_box[0][0] ; i <= local_box[0][1] ; ++i ) {
00189     EntityType node_type = 0;
00190     EntityId node_id = 1 + i + j * (ngx+1) + k * (ngx+1) * (ngy+1);
00191     Entity * const node = mesh.get_entity( node_type , node_id );
00192     STKUNIT_ASSERT( node != NULL );
00193     // Shared if on a processor boundary.
00194     const bool shared =
00195       ( k == local_box[2][0] && k != root_box[2][0] ) ||
00196       ( k == local_box[2][1] && k != root_box[2][1] ) ||
00197       ( j == local_box[1][0] && j != root_box[1][0] ) ||
00198       ( j == local_box[1][1] && j != root_box[1][1] ) ||
00199       ( i == local_box[0][0] && i != root_box[0][0] ) ||
00200       ( i == local_box[0][1] && i != root_box[0][1] );
00201     if (mesh.parallel_size() > 1) {
00202       STKUNIT_ASSERT_EQUAL( shared , ! node->sharing().empty() );
00203     }
00204   }
00205   }
00206   }
00207 
00208   size_t count_shared_node_pairs = 0 ;
00209   for ( unsigned p = 0 ; p < p_size ; ++p ) if ( p != p_rank ) {
00210     for ( int k = p_box[p][2][0] ; k <= p_box[p][2][1] ; ++k )
00211     if ( local_box[2][0] <= k && k <= local_box[2][1] ) {
00212 
00213       for ( int j = p_box[p][1][0] ; j <= p_box[p][1][1] ; ++j )
00214       if ( local_box[1][0] <= j && j <= local_box[1][1] ) {
00215 
00216         for ( int i = p_box[p][0][0] ; i <= p_box[p][0][1] ; ++i )
00217         if ( local_box[0][0] <= i && i <= local_box[0][1] ) {
00218 
00219           EntityType node_type = 0;
00220           EntityId node_id = 1 + i + j * (ngx+1) + k * (ngx+1) * (ngy+1);
00221           Entity * const node = mesh.get_entity( node_type , node_id );
00222           STKUNIT_ASSERT( node != NULL );
00223           // Must be shared with 'p'
00224           PairIterEntityComm iter = node->sharing();
00225           for ( ; ! iter.empty() && iter->proc != p ; ++iter );
00226           STKUNIT_ASSERT( ! iter.empty() );
00227 
00228           ++count_shared_node_pairs ;
00229         }
00230       }
00231     }
00232   }
00233 
00234   size_t count_shared_entities = 0 ;
00235   for ( std::vector<Entity*>::const_iterator
00236         i = mesh.entity_comm().begin() ; i != mesh.entity_comm().end() ; ++i ) {
00237     const PairIterEntityComm ec = (**i).sharing();
00238     count_shared_entities += ec.size();
00239   }
00240   STKUNIT_ASSERT_EQUAL( count_shared_entities , count_shared_node_pairs );
00241 
00242   delete[] p_box ;
00243 }
00244 
00245 } // namespace mesh
00246 } // namespace stk
00247 

Generated on Tue Jul 13 09:27:32 2010 for Sierra Toolkit by  doxygen 1.4.7