Sierra Toolkit Version of the Day
UseCase_Rebal_3.cpp
00001 /*------------------------------------------------------------------------*/
00002 /*                 Copyright 2010, 2011 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 #include <use_cases/UseCase_Rebal_3.hpp>
00010 
00011 #include <stk_util/parallel/Parallel.hpp>
00012 #include <stk_util/parallel/ParallelReduce.hpp>
00013 
00014 #include <stk_mesh/base/FieldData.hpp>
00015 #include <stk_mesh/base/GetEntities.hpp>
00016 #include <stk_mesh/base/MetaData.hpp>
00017 #include <stk_mesh/base/BulkData.hpp>
00018  
00019 #include <stk_mesh/fem/FEMMetaData.hpp>
00020 #include <stk_mesh/fem/FEMHelpers.hpp>
00021 #include <stk_mesh/fem/CoordinateSystems.hpp>
00022 #include <stk_mesh/fem/CreateAdjacentEntities.hpp>
00023 
00024 #include <stk_mesh/fixtures/HexFixture.hpp>
00025 
00026 #include <stk_rebalance/Rebalance.hpp>
00027 #include <stk_rebalance/Partition.hpp>
00028 #include <stk_rebalance/ZoltanPartition.hpp>
00029 
00030 #include <stk_rebalance_utils/RebalanceUtils.hpp>
00031 
00032 //----------------------------------------------------------------------
00033 
00034 using namespace stk::mesh::fixtures;
00035 
00036 typedef stk::mesh::Field<double> ScalarField ;
00037 
00038 namespace stk {
00039 namespace rebalance {
00040 namespace use_cases {
00041 
00042 class Test_Case_3_Partition : public stk::rebalance::Zoltan {
00043   public :
00044   explicit Test_Case_3_Partition(ParallelMachine pm,
00045                                  const unsigned ndim,
00046                                  Teuchos::ParameterList & rebal_region_parameters,
00047                                  std::string parameters_name=default_parameters_name());
00048   virtual ~Test_Case_3_Partition();
00049   virtual bool partition_dependents_needed() const; 
00050 };
00051 
00052 Test_Case_3_Partition::Test_Case_3_Partition(ParallelMachine pm,
00053                                              const unsigned ndim,
00054                                              Teuchos::ParameterList & rebal_region_parameters,
00055                                              std::string parameters_name) :
00056   stk::rebalance::Zoltan(pm, ndim, rebal_region_parameters, parameters_name) {}
00057 Test_Case_3_Partition::~Test_Case_3_Partition() {}
00058 bool Test_Case_3_Partition::partition_dependents_needed() const 
00059 {return false;} //Do NOT move dependent entities for this case
00060 
00061 enum { nx = 4, ny = 4 };
00062 
00063 bool test_contact_surfaces( stk::ParallelMachine comm )
00064 {
00065   unsigned spatial_dimension = 2;
00066   std::vector<std::string> rank_names = stk::mesh::fem::entity_rank_names(spatial_dimension);
00067   const stk::mesh::EntityRank constraint_rank = rank_names.size();
00068   rank_names.push_back("Constraint");
00069 
00070   stk::mesh::fem::FEMMetaData fem_meta;
00071   fem_meta.FEM_initialize(spatial_dimension, rank_names);
00072 
00073   stk::mesh::MetaData & meta_data = stk::mesh::fem::FEMMetaData::get_meta_data(fem_meta);
00074   stk::mesh::BulkData bulk_data( meta_data , comm , 100 );
00075   const stk::mesh::EntityRank element_rank    = fem_meta.element_rank();
00076   const stk::mesh::EntityRank side_rank       = fem_meta.side_rank();
00077   const stk::mesh::EntityRank node_rank       = fem_meta.node_rank();
00078 
00079   stk::mesh::fem::CellTopology quad_top(shards::getCellTopologyData<shards::Quadrilateral<4> >());
00080   stk::mesh::Part & quad_part( fem_meta.declare_part("quad", quad_top ) );
00081   stk::mesh::fem::CellTopology side_top(shards::getCellTopologyData<shards::Line<2> >());
00082   stk::mesh::Part & side_part( fem_meta.declare_part("line", side_top ) );
00083   VectorField & coord_field( fem_meta.declare_field< VectorField >( "coordinates" ) );
00084   ScalarField & weight_field( fem_meta.declare_field< ScalarField >( "element_weights" ) );
00085 
00086   stk::mesh::put_field( coord_field , node_rank , fem_meta.universal_part() );
00087   stk::mesh::put_field(weight_field , element_rank , fem_meta.universal_part() );
00088 
00089   fem_meta.commit();
00090 
00091   //const unsigned p_size = bulk_data.parallel_size();
00092   const unsigned p_rank = bulk_data.parallel_rank();
00093 
00094   bulk_data.modification_begin();
00095 
00096   if ( !p_rank ) { 
00097 
00098     std::vector<stk::mesh::EntityVector> quads(nx);
00099     for ( unsigned ix = 0 ; ix < nx ; ++ix ) quads[ix].resize(ny);
00100 
00101     const unsigned nnx = nx + 1 ; 
00102     const unsigned nny = ny + 1 ; 
00103     unsigned face_id   = 1;
00104     for ( unsigned iy = 0 ; iy < ny ; ++iy ) { 
00105       for ( unsigned ix = 0 ; ix < nx ; ++ix ) { 
00106         stk::mesh::EntityId elem = 1 + ix + iy * nx ;
00107         stk::mesh::EntityId nodes[4] ;
00108         nodes[0] = 1 + ix + iy * nnx ;
00109         nodes[1] = 2 + ix + iy * nnx ;
00110         nodes[2] = 2 + ix + ( iy + 1 ) * nnx ;
00111         nodes[3] = 1 + ix + ( iy + 1 ) * nnx ;
00112 
00113         stk::mesh::Entity &q = stk::mesh::fem::declare_element( bulk_data , quad_part , elem , nodes );
00114         if (0==ix) {
00115           stk::mesh::fem::declare_element_side( bulk_data, face_id, q, 0 /*local side id*/, &side_part);
00116           ++face_id;
00117         }
00118         quads[ix][iy] = &q; 
00119       }   
00120     }   
00121 
00122     for ( unsigned iy = 0 ; iy < ny ; ++iy ) { 
00123       for ( unsigned ix = 0 ; ix < nx ; ++ix ) { 
00124         stk::mesh::EntityId elem = 1 + ix + iy * nx ;
00125         stk::mesh::Entity * e = bulk_data.get_entity( element_rank, elem );
00126         double * const e_weight = stk::mesh::field_data( weight_field , *e );
00127         *e_weight = 1.0;
00128       }   
00129     }   
00130     for ( unsigned iy = 0 ; iy <= ny ; ++iy ) { 
00131       for ( unsigned ix = 0 ; ix <= nx ; ++ix ) { 
00132         stk::mesh::EntityId nid = 1 + ix + iy * nnx ;
00133         stk::mesh::Entity * n = bulk_data.get_entity( node_rank, nid );
00134         double * const coord = stk::mesh::field_data( coord_field , *n );
00135         coord[0] = .1*ix;
00136         coord[1] = .1*iy;
00137         coord[2] = 0;
00138       }
00139     }
00140 
00141     // Assign constraint relations between nodes at top and bottom of mesh
00142     {
00143       const unsigned iy_bottom  =  0;
00144       const unsigned iy_top = ny;
00145       stk::mesh::PartVector add(1, &fem_meta.locally_owned_part());
00146       for ( unsigned ix = 0 ; ix <= nx ; ++ix ) {
00147         stk::mesh::EntityId nid_bottom  = 1 + ix + iy_bottom  * nnx ;
00148         stk::mesh::EntityId nid_top = 1 + ix + iy_top * nnx ;
00149         stk::mesh::Entity * n_bottom  = bulk_data.get_entity( node_rank, nid_bottom  );
00150         stk::mesh::Entity * n_top = bulk_data.get_entity( node_rank, nid_top );
00151         const stk::mesh::EntityId constraint_entity_id =  1 + ix + nny * nnx;
00152         stk::mesh::Entity & c = bulk_data.declare_entity( constraint_rank, constraint_entity_id, add );
00153         bulk_data.declare_relation( c , *n_bottom  , 0 );
00154         bulk_data.declare_relation( c , *n_top , 1 );
00155       }
00156     } // end snippet
00157 
00158   }
00159 
00160   bulk_data.modification_end();
00161 
00162   stk::mesh::Selector selector(side_part);
00163   selector &=  fem_meta.locally_owned_part();
00164   // Coordinates are passed to support geometric-based load balancing algorithms
00165   // Zoltan partition is specialized form a virtual base class, stk::rebalance::Partition.
00166   // Other specializations are possible.
00167   Teuchos::ParameterList emptyList;
00168   stk::rebalance::use_cases::Test_Case_3_Partition zoltan_partition(comm, spatial_dimension, emptyList);
00169   stk::rebalance::rebalance(bulk_data, selector, &coord_field, NULL, zoltan_partition, side_rank);
00170 
00171   const double imbalance_threshold = 1.5;
00172   const bool do_rebal = imbalance_threshold < stk::rebalance::check_balance(bulk_data, NULL, side_rank, &selector);
00173 
00174   if( !p_rank )
00175     std::cerr << std::endl 
00176      << "imbalance_threshold after rebalance = " << imbalance_threshold <<", "<<do_rebal << std::endl;
00177 
00178   // Check that we satisfy our threshhold
00179   bool result = !do_rebal ;
00180 
00181   return result;
00182 }
00183 
00184 } //namespace use_cases
00185 } //namespace rebalance
00186 } //namespace stk
00187 
00188 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends