Sierra Toolkit Version of the Day
ZoltanPartition.cpp
00001 
00002 #include <vector>
00003 #include <string>
00004 #include <sstream>
00005 #include <cstdlib>
00006 
00007 #include <stk_rebalance/ZoltanPartition.hpp>
00008 #include <stk_mesh/base/Field.hpp>
00009 #include <stk_mesh/base/FieldData.hpp>
00010 #include <stk_mesh/base/Entity.hpp>
00011 #include <stk_mesh/base/Bucket.hpp>
00012 #include <stk_mesh/base/BulkData.hpp>
00013 #include <stk_mesh/fem/FEMMetaData.hpp>
00014 
00015 #include <stk_util/parallel/ParallelReduce.hpp>
00016 
00017 #include <zoltan.h>
00018 
00019 #include <Teuchos_ParameterList.hpp>
00020 
00021 using namespace std;
00022 using namespace stk;
00023 using namespace stk::rebalance;
00024 
00025 #define STK_GEOMDECOMP_DEBUG 0
00026 
00027 namespace {
00028 
00029 
00030 double static_zoltan_version(const double v=0) {
00031   static double version=0;
00032   if (v) { version=v;}
00033   return version;
00034 }
00035 
00036 inline unsigned num_gid_entries() {
00037   static const unsigned n=2;
00038   return n;
00039 }
00040 
00041 inline unsigned num_lid_entries() {
00042   static const unsigned n=2;
00043   return n;
00044 }
00045 inline unsigned wdim() {
00046   static const unsigned n=1;
00047   return n;
00048 }
00049 
00050 inline void convert_param_to_string(const Parameters &from,
00051                                     vector < pair<std::string, std::string> > &to)
00052 {
00053   Parameters::ConstIterator
00054     from_iter  = from.begin(),
00055     from_end   = from.end();
00056 
00057   for (; from_iter != from_end; ++from_iter) {
00058     std::string s;
00059     const std::string name  = from.name(from_iter);
00060     const std::string entry = from.entry(from_iter).getValue(&s);
00061     std::pair<std::string,std::string> p(name,entry);
00062 
00063     to.push_back(p);
00064   }
00065 }
00066 
00067 inline void fill_parameters (const char *T[][2],
00068                              const int  i,
00069                              Parameters &Entry)
00070 {
00071   for (int j=0; j<i; ++j) Entry.set(T[j][0], T[j][1]);
00072 }
00073 
00074 void fill_name_conversion( Parameters & name_conversion )
00075 {
00076 
00077   Parameters & general = name_conversion.sublist("General");
00078 
00079   general.set("LOAD BALANCING METHOD"      , "LB_METHOD");
00080   general.set("ZOLTAN DEBUG LEVEL"         , "DEBUG_LEVEL");
00081   general.set("DEBUG PROCESSOR NUMBER"     , "DEBUG_PROCESSOR");
00082   general.set("TIMER"                      , "TIMER");
00083   general.set("DETERMINISTIC DECOMPOSITION", "DETERMINISTIC");
00084   general.set("DEBUG MEMORY"               , "DEBUG_MEMORY");
00085   general.set("IMBALANCE TOLERANCE"        , "IMBALANCE_TOL");
00086   general.set("RENUMBER PARTITIONS"        , "REMAP");
00087   general.set("KEEP CUTS"                  , "KEEP_CUTS");
00088   general.set("REUSE CUTS"                 , "RCB_REUSE");
00089   general.set("RCB RECOMPUTE BOX"          , "RCB_RECOMPUTE_BOX");
00090   general.set("CHECK GEOMETRY"             , "CHECK_GEOM");
00091   general.set("LOCK RCB DIRECTIONS"        , "RCB_LOCK_DIRECTIONS");
00092   general.set("SET RCB DIRECTIONS"         , "RCB_SET_DIRECTIONS");
00093   general.set("RCB MAX ASPECT RATIO"       , "RCB_MAX_ASPECT_RATIO");
00094   general.set("RECTILINEAR RCB BLOCKS"     , "RCB_RECTILINEAR_BLOCKS");
00095   general.set("OCTREE DIMENSION"           , "OCT_DIM");
00096   general.set("OCTREE METHOD"              , "OCT_METHOD");
00097   general.set("OCTREE MIN ENTITIES"         , "OCT_MINOBJECTS");
00098   general.set("OCTREE MAX ENTITIES"         , "OCT_MAXOBJECTS");
00099 //  // These values are never changed, but must
00100   // be set so default values work correctly.
00101   general.set("NUMBER GLOBAL ID ENTRIES"   , "NUM_GID_ENTRIES");
00102   general.set("NUMBER LOCAL ID ENTRIES"    , "NUM_LID_ENTRIES");
00103   general.set("ENTITY WEIGHT DIMENSION"    , "OBJ_WEIGHT_DIM");
00104   general.set("RETURN LISTS"               , "RETURN_LISTS");
00105   general.set("AUTOMATIC MIGRATION"        , "AUTO_MIGRATE");
00106   general.set("DISTANCE"                   , "DISTANCE");
00107 
00108   Parameters & rcb = name_conversion.sublist("0");
00109   rcb.set("OVER ALLOCATE MEMORY"       , "RCB_OVERALLOC");
00110   rcb.set("ALGORITHM DEBUG LEVEL"      , "RCB_OUTPUT_LEVEL");
00111 
00112   Parameters & rib = name_conversion.sublist("1");
00113   rib.set("OVER ALLOCATE MEMORY"       , "RIB_OVERALLOC");
00114   rib.set("ALGORITHM DEBUG LEVEL"      , "RIB_OUTPUT_LEVEL");
00115 
00116   Parameters & hsfc = name_conversion.sublist("2");
00117   hsfc.set("OVER ALLOCATE MEMORY"       , "");
00118   hsfc.set("ALGORITHM DEBUG LEVEL"      , "" );
00119 
00120   Parameters & oct = name_conversion.sublist("3");
00121   oct.set("OVER ALLOCATE MEMORY"       , "");
00122   oct.set("ALGORITHM DEBUG LEVEL"      , "OCT_OUTPUT_LEVEL");
00123 
00124   Parameters & graph = name_conversion.sublist("4");
00125   graph.set("OVER ALLOCATE MEMORY"       , "");
00126   graph.set("ALGORITHM DEBUG LEVEL"      , "" );
00127 }
00128 
00129 
00130 void fill_value_conversion( Parameters & value_conversion )
00131 {
00132   Parameters & lb_method = value_conversion.sublist("LOAD BALANCING METHOD");
00133   lb_method.set("0"   , "RCB");
00134   lb_method.set("1"   , "RIB");
00135   lb_method.set("2"   , "HSFC");
00136   lb_method.set("3"   , "OCTPART");
00137   lb_method.set("4"   , "GRAPH");
00138 
00139   Parameters & timer = value_conversion.sublist("TIMER");
00140   timer.set("0"   , "WALL");
00141   timer.set("1"   , "CPU");
00142 
00143 }
00144 
00145 void fill_default_values( Parameters & values )
00146 {
00147   Parameters & default_values = values; //values.sublist("General");
00148 
00149   default_values.set("LOAD BALANCING METHOD"      , "0");
00150   default_values.set("RENUMBER PARTITIONS"        , "1");
00151   default_values.set("ZOLTAN DEBUG LEVEL"         , "0");
00152   default_values.set("TIMER"                      , "0");
00153   default_values.set("DETERMINISTIC DECOMPOSITION", "1");
00154   default_values.set("DEBUG MEMORY"               , "1");
00155   default_values.set("IMBALANCE TOLERANCE"        , "1.1");
00156   default_values.set("KEEP CUTS"                  , "1");
00157   default_values.set("REUSE CUTS"                 , "1");
00158   default_values.set("OVER ALLOCATE MEMORY"       , "1.1");
00159   default_values.set("ALGORITHM DEBUG LEVEL"      , "0");
00160 //  default_values.set("OCTREE MIN ENTITIES"        , "1");
00161 //  default_values.set("OCTREE MAX ENTITIES"        , "1");
00162   default_values.set("NUMBER GLOBAL ID ENTRIES"   , "2");
00163   default_values.set("NUMBER LOCAL ID ENTRIES"    , "2");
00164   default_values.set("ENTITY WEIGHT DIMENSION"    , "1");
00165   default_values.set("RETURN LISTS"               , "EXPORT");
00166 }
00167 
00168 
00169 
00170 #if STK_GEOMDECOMP_DEBUG>=2
00171 void debug_print_decomp_export(Zoltan   *zoltan,
00172                                Zoltan *zoltan_id )
00173 {
00174   int i;
00175   int           num_export   = zoltan->Num_Exported();
00176   ZOLTAN_ID_PTR export_lids  = zoltan->Export_Local_IDs();
00177   ZOLTAN_ID_PTR export_gids  = zoltan->Export_Global_IDs();
00178   int*          export_procs = zoltan->Export_Proc_IDs();
00179   int           Z_LID_SIZE   = zoltan->Num_Lid_Entries();
00180   int           Z_GID_SIZE   = zoltan->Num_Gid_Entries();
00181 
00182   if ( export_gids!=NULL && export_lids!=NULL && export_procs!=NULL ) {
00183     Env::output() << ": Zoltan RCB EXPORTS" << std::endl;
00184     for ( i = 0; i < num_export; i++ ) {
00185       Env::output()
00186         << "  " << i
00187         << ":  GID = "
00188         << "T" << zoltan_id->Type( &export_gids[i*Z_GID_SIZE]) << "  "
00189         << "I" << zoltan_id->Index(&export_gids[i*Z_GID_SIZE]) << "  "
00190         << "P" << zoltan_id->Proc( &export_gids[i*Z_GID_SIZE]) << "  "
00191         << "    LID = "
00192         << "T" << zoltan_id->Type( &export_lids[i*Z_LID_SIZE]) << "  "
00193         << "I" << zoltan_id->Index(&export_lids[i*Z_LID_SIZE]) << "  "
00194         << "  EXPORT_PROC_ID = "
00195         << export_procs[i]
00196         << std::endl;
00197     }
00198     Env::output_flush();
00199   }
00200 }
00201 #endif
00202 
00203 
00204 
00205 extern "C" {
00206   int Callback_Num_Elements( void *data, int *ierr );
00207   void Callback_Element_List( void *data,
00208                              int Num_gid_entries,
00209                              int Num_lid_entries,
00210                              ZOLTAN_ID_PTR global_ids,
00211                              ZOLTAN_ID_PTR local_ids,     //{Iterator or index}
00212                              int wdim,
00213                              float *weights,
00214                              int *ierr );
00215   int Callback_Num_Dimensions( void *data, int *ierr );
00216   void Callback_Centroid_Coord( void *data,
00217                                 int Num_gid_entries,
00218                                 int Num_lid_entries,
00219                                 ZOLTAN_ID_PTR global_id,
00220                                 ZOLTAN_ID_PTR local_id,  //{Iterator or index}
00221                                 double *geom,
00222                                 int *ierr );
00223   int Callback_Num_Edges( void *data,
00224                           int Num_gid_entries,
00225                           int Num_lid_entries,
00226                           ZOLTAN_ID_PTR global_id,
00227                           ZOLTAN_ID_PTR local_id,  //{Iterator or index}
00228                           int *ierr );
00229   void Callback_Edge_List( void *data,
00230                            int Num_gid_entries,
00231                            int Num_lid_entries,
00232                            ZOLTAN_ID_PTR global_id,
00233                            ZOLTAN_ID_PTR local_id,  //{Iterator or index}
00234                            ZOLTAN_ID_PTR nbor_global_id,
00235                            int *nbor_procs,
00236                            int wgt_dim,
00237                            float *ewgts,
00238                            int *ierr );
00239 }
00240 
00241 
00242 int Callback_Num_Elements( void *data, int *ierr )
00243 {
00244   if ( data == NULL ) {
00245     *ierr = ZOLTAN_FATAL;  // Set FATAL Zoltan error flag
00246     return 0;
00247   }
00248   stk::rebalance::GeomDecomp   *gdata =  static_cast<stk::rebalance::GeomDecomp*>  (data);
00249   stk::rebalance::Zoltan       *zdata = dynamic_cast<stk::rebalance::Zoltan*>(gdata);
00250 
00251   if (zdata == NULL ) {
00252     *ierr = ZOLTAN_FATAL;  // Set FATAL Zoltan error flag
00253     return 0;
00254   }
00255 
00256   int ne = zdata->num_elems();
00257 
00258   *ierr = ZOLTAN_OK;
00259   return ne;
00260 }
00261 
00262 void Callback_Element_List( void *data,
00263                             int Num_gid_entries,
00264                             int Num_lid_entries,
00265                             ZOLTAN_ID_PTR global_ids,
00266                             ZOLTAN_ID_PTR local_ids,     //{Iterator or index}
00267                             int weightdim,
00268                             float *weights,
00269                             int *ierr )
00270 {
00271   if (!data) {
00272     *ierr = ZOLTAN_FATAL;           // Set FATAL Zoltan error flag
00273     return;
00274   }
00275 
00276   stk::rebalance::GeomDecomp   *gdata =  static_cast<stk::rebalance::GeomDecomp*>  (data);
00277   stk::rebalance::Zoltan       *zdata = dynamic_cast<stk::rebalance::Zoltan*>     (gdata);
00278 
00279   if (!zdata) {
00280     *ierr = ZOLTAN_FATAL;           // Set FATAL Zoltan error flag
00281     return;
00282   }
00283 
00284   unsigned k=0, l=0;
00285   const unsigned num_local_ids = zdata->num_moid();
00286   for (unsigned j=0; j<num_local_ids; ++j) {
00287     local_ids [k] = j;
00288     global_ids[k] = 0; // region_id
00289     ++k;
00290     local_ids [k] = 0;
00291     global_ids[k] = static_cast<ZOLTAN_ID_TYPE>(zdata->globalID(j));
00292     ++k;
00293     if (weightdim) {
00294       weights[l++] = zdata->entity_weight(j);
00295     } else {
00296       ++l;
00297     }
00298   }
00299   *ierr = ZOLTAN_OK;
00300   return;
00301 
00302 }
00303 
00304 int Callback_Num_Dimensions( void *data, int *ierr )
00305 {
00306   if ( !data ) {
00307     *ierr = ZOLTAN_FATAL;
00308     return 0;
00309   }
00310   stk::rebalance::GeomDecomp *gdata = static_cast<stk::rebalance::GeomDecomp*>  (data);
00311   stk::rebalance::Zoltan     *zdata = dynamic_cast<stk::rebalance::Zoltan*>    (gdata);
00312 
00313   if ( !zdata ) {
00314     *ierr = ZOLTAN_FATAL;
00315     return 0;
00316   }
00317 
00318   const int  nd = zdata->spatial_dimension();
00319 
00320   *ierr = ZOLTAN_OK;
00321   return nd;
00322 
00323 }
00324 
00325 void Callback_Centroid_Coord( void *data,
00326                                      int Num_gid_entries,
00327                                      int Num_lid_entries,
00328                                      ZOLTAN_ID_PTR global_id,
00329                                      ZOLTAN_ID_PTR local_id,  //{Iterator or index}
00330                                      double *geom,
00331                                      int *ierr )
00332 {
00333   std::vector<double> temp(3,0.0);
00334 
00335   // (from include Fmwk_Sierra_Zoltan_Defines.h:)
00336   // ( Num_gid_entries = ZOLTAN_GID_SIZE 2 )
00337   // ( Num_lid_entries = ZOLTAN_LID_SIZE 2 )
00338 
00339   if ( !data ) {
00340     *ierr = ZOLTAN_FATAL;
00341     return ;
00342   }
00343   stk::rebalance::GeomDecomp *gdata = static_cast<stk::rebalance::GeomDecomp*>  (data);
00344   stk::rebalance::Zoltan     *zdata = dynamic_cast<stk::rebalance::Zoltan*>    (gdata);
00345 
00346   if ( !zdata ) {
00347     *ierr = ZOLTAN_FATAL;
00348     return ;
00349   }
00350 
00351   int lid = local_id[  0 ]; // Local Element ID
00352 
00353   const mesh::Entity & target_entity = * zdata->mesh_entity( lid );
00354   const GeomDecomp::VectorField & coor = * zdata->entity_coord_ref();
00355   const unsigned                        nd   =  zdata->spatial_dimension();
00356 
00357   /*
00358    * Obtain the centroid coordinates of the element by averaging all
00359    * the nodal coordinates of the nodes associated with the element.
00360    * Use GeomDecomp friend function, ent( , , )
00361    *
00362    */
00363 
00364   rebalance::GeomDecomp::entity_to_point( target_entity, coor, temp );
00365 
00366   for (size_t i=0 ; i < nd ; i++ ) geom[ i ] = (double) temp[ i ];
00367 
00368   *ierr = ZOLTAN_OK;
00369 }
00370 
00371 
00372 
00373 void getNeighbors( const mesh::Entity & entity,
00374                    std::set<const mesh::Entity*> & nodes )
00375 {
00376   stk::mesh::fem::FEMMetaData &fem_meta = stk::mesh::fem::FEMMetaData::get(entity);
00377   const stk::mesh::EntityRank element_rank = fem_meta.element_rank();
00378 
00379   nodes.clear();
00380 
00381   mesh::PairIterRelation iElem = entity.relations(element_rank);
00382 
00383   for ( ; iElem.first != iElem.second; ++iElem.first ) {
00384     mesh::Entity * elem = iElem.first->entity();
00385     mesh::PairIterRelation iNode = elem->relations(fem_meta.node_rank());
00386     for ( ; iNode.first != iNode.second; ++iNode.first ) {
00387       mesh::Entity * node = iNode.first->entity();
00388       if (&entity != node) nodes.insert( node );
00389     }
00390   }
00391 }
00392 
00393 int numEdges( const mesh::Entity & entity ) {
00394 
00395   std::set<const mesh::Entity*> nodes;
00396 
00397   getNeighbors( entity, nodes );
00398   return nodes.size();
00399 }
00400 
00401 int Callback_Num_Edges( void *data,
00402                         int Num_gid_entries,
00403                         int Num_lid_entries,
00404                         ZOLTAN_ID_PTR global_id,
00405                         ZOLTAN_ID_PTR local_id,  //{Iterator or index}
00406                         int *ierr )
00407 {
00408   // (from include Fmwk_Sierra_Zoltan_Defines.h:)
00409   // ( Num_gid_entries = ZOLTAN_GID_SIZE 2 )
00410   // ( Num_lid_entries = ZOLTAN_LID_SIZE 2 )
00411 
00412   *ierr = ZOLTAN_OK;
00413 
00414   if ( !data ) {
00415     *ierr = ZOLTAN_FATAL;
00416     return 0;
00417   }
00418   stk::rebalance::GeomDecomp *gdata = static_cast<stk::rebalance::GeomDecomp*>  (data);
00419   stk::rebalance::Zoltan     *zdata = dynamic_cast<stk::rebalance::Zoltan*>    (gdata);
00420 
00421   if ( !zdata ) {
00422     *ierr = ZOLTAN_FATAL;
00423     return 0;
00424   }
00425 
00426   int lid = local_id[  0 ]; // Local Element ID
00427 
00428   const mesh::Entity & target_entity = * zdata->mesh_entity( lid );
00429 
00430   return  numEdges( target_entity );
00431 
00432 }
00433 
00434 void Callback_Edge_List( void *data,
00435                          int Num_gid_entries,
00436                          int Num_lid_entries,
00437                          ZOLTAN_ID_PTR global_id,
00438                          ZOLTAN_ID_PTR local_id,  //{Iterator or index}
00439 
00440                          // Returned
00441                          ZOLTAN_ID_PTR nbor_global_id, //owning processor
00442                          int *nbor_procs,
00443                          int wgt_dim,
00444                          float *ewgts,
00445 
00446                          int *ierr )
00447 {
00448   // (from include Fmwk_Sierra_Zoltan_Defines.h:)
00449   // ( Num_gid_entries = ZOLTAN_GID_SIZE 2 )
00450   // ( Num_lid_entries = ZOLTAN_LID_SIZE 2 )
00451 
00452   *ierr = ZOLTAN_OK;
00453 
00454   if ( !data ) {
00455     *ierr = ZOLTAN_FATAL;
00456     return ;
00457   }
00458   stk::rebalance::GeomDecomp *gdata = static_cast<stk::rebalance::GeomDecomp*>  (data);
00459   stk::rebalance::Zoltan     *zdata = dynamic_cast<stk::rebalance::Zoltan*>    (gdata);
00460 
00461   if ( !zdata ) {
00462     *ierr = ZOLTAN_FATAL;
00463     return ;
00464   }
00465 
00466   int lid = local_id[  0 ]; // Local Node ID
00467 
00468   const mesh::Entity & target_entity = * zdata->mesh_entity( lid );
00469 
00470   std::set<const mesh::Entity*> nodes;
00471   getNeighbors( target_entity, nodes );
00472 
00473   int counter(0);
00474   for ( std::set<const mesh::Entity*>::iterator i = nodes.begin();
00475         i != nodes.end(); ++i ) {
00476     nbor_global_id[counter*2+0] = 0; // region_id
00477 
00478     const mesh::Entity & n = **i;
00479     nbor_global_id[counter*2+1] = n.key().id();
00480 
00481     nbor_procs[counter] = n.owner_rank();
00482 
00483     if ( wgt_dim ) {
00484       ewgts[counter] = 1;
00485     }
00486     ++counter;
00487   }
00488 }
00489 
00490 }
00491 
00492 
00493 
00494 const std::string Zoltan::m_zoltanparametersname_="Zoltan_Parameters";
00495 const std::string Zoltan::m_defaultparametersname_="DEFAULT";
00496 
00497 
00498 const std::string Zoltan::zoltan_parameters_name()
00499 {
00500   return m_zoltanparametersname_;
00501 }
00502 
00503 const std::string Zoltan::default_parameters_name()
00504 {
00505   return m_defaultparametersname_;
00506 }
00507 
00508 void Zoltan::init_default_parameters()
00509 {
00510   fill_default_values(m_default_parameters_);
00511 }
00512 
00513 
00514 static Parameters *Name_Conversion =NULL;
00515 static Parameters *Value_Conversion=NULL;
00516 
00517 double Zoltan::zoltan_version()   const { return static_zoltan_version();  }
00518 
00519 //: ===========
00520 //: Constructor
00521 //: ===========
00522 
00523 namespace {
00524 void merge_parameters(std::vector <std::pair<std::string, std::string> > &str_zoltan_params,
00525                       const Parameters &Zoltan_Params) {
00526   Parameters Merged_Zoltan_Params   ;
00527   Parameters Converted_Zoltan_Params;
00528 
00529   rebalance::Zoltan::merge_default_values (Zoltan_Params,
00530                                       Merged_Zoltan_Params);
00531 
00532   rebalance::Zoltan::convert_names_and_values(Merged_Zoltan_Params,
00533                                          Converted_Zoltan_Params);
00534 
00535   convert_param_to_string (Converted_Zoltan_Params,
00536                            str_zoltan_params);
00537   return;
00538 }
00539 }
00540 
00541 Zoltan::Zoltan(ParallelMachine pm, const unsigned ndim, Parameters & rebal_region_parameters, const std::string parameters_name) :
00542   GeomDecomp(pm),
00543   m_zoltan_id_(NULL),
00544   m_spatial_dimension_(ndim),
00545   m_total_number_entities_(0)
00546 {
00547   /* Determine if the default set of parameters already exists. */
00548   if( !rebal_region_parameters.isSublist(default_parameters_name()) )
00549   {
00550     init_default_parameters();
00551     rebal_region_parameters.sublist(default_parameters_name()) = m_default_parameters_;
00552   }
00553 
00554   /* If name is empty, use default values */
00555   std::string default_name =
00556     (parameters_name.empty()) ? default_parameters_name() : parameters_name ;
00557 
00558   if( !rebal_region_parameters.isSublist(default_name) ) {
00559     throw std::runtime_error("The Zoltan parameter set '" + default_name + "' does not exist.");
00560   }
00561   const Parameters & zoltan_params = rebal_region_parameters.sublist(default_name);
00562 
00563   std::vector <std::pair<std::string, std::string> > str_zoltan_params;
00564   merge_parameters(str_zoltan_params, zoltan_params);
00565 
00566   init(str_zoltan_params);
00567 }
00568 
00569 void
00570 Zoltan::set_mesh_info( const std::vector<mesh::Entity *> &mesh_entities,
00571                           const VectorField * nodal_coord_ref,
00572                           const ScalarField * elem_weight_ref)
00573 {
00574   MeshInfo mesh_info;
00575 
00576   /* Keep track of the total number of elements. */
00577   m_total_number_entities_ = mesh_entities.size();
00578 
00579   mesh_info.mesh_entities = mesh_entities;
00580   mesh_info.nodal_coord_ref = nodal_coord_ref;
00581   mesh_info.elem_weight_ref = elem_weight_ref;
00582 
00588   mesh_info.dest_proc_ids.assign(mesh_entities.size(), stk::parallel_machine_rank(comm_));
00589 
00590   m_mesh_information_ = mesh_info;
00591 }
00592 
00593 void Zoltan::init( const vector< pair<std::string,std::string> >
00594                    &dynamicLoadRebalancingParameters ) {
00595   if (0==static_zoltan_version()) {
00596     const double v = init_zoltan_library();
00597     static_zoltan_version(v);
00598   }
00599 
00604   m_zoltan_id_ = Zoltan_Create( comm_ );
00605   if ( m_zoltan_id_ == NULL ) {
00606     throw runtime_error ("(FATAL ERROR) Zoltan_Create() returned NULL");
00607   }
00608 
00613   vector<pair<std::string,std::string> >::const_iterator
00614     P  =  dynamicLoadRebalancingParameters.begin(),
00615     PE =  dynamicLoadRebalancingParameters.end();
00616 
00617   for ( ; PE != P ; P++ ) {
00618 
00619     char * label = const_cast<char*>( P->first.c_str() ) ;
00620     char * value = const_cast<char*>( P->second.c_str() ) ;
00621 
00622     if (ZOLTAN_OK != (Zoltan_Set_Param(m_zoltan_id_,label,value)))
00623     {
00624       throw runtime_error(": FATAL ERROR returned from Zoltan_Set_Param ");
00625     }
00626   }
00627 
00631   if ( ZOLTAN_OK != register_callbacks() )
00632     throw runtime_error ("zoltan->Register_Callbacks error. ");
00633 
00634 #if STK_GEOMDECOMP_DEBUG>=2
00635   {
00636     debug_print_decomp_export( zoltan, m_zoltan_id_ );
00637   }
00638 #endif
00639   return;
00640 
00641 }
00642 
00643 double Zoltan::init_zoltan_library () {
00644   float version = 0.0;
00645   if ( Zoltan_Initialize( 0, NULL, &version) != ZOLTAN_OK )
00646     throw std::runtime_error("Return code from Zoltan_Initialize() != ZOLTAN_OK ");
00647 
00648   static_zoltan_version(version);
00649   std::ostringstream s;
00650   s << version;
00651 
00652   //sierra::ProductRegistry::instance().addTPL("Zoltan", s.str());
00653 
00654   return version;
00655 }
00656 
00657 
00658 
00659 //: ==========
00660 //: Destructor
00661 //: ==========
00662 
00663 Zoltan::~Zoltan()
00664 {
00665   if ( m_zoltan_id_ != NULL ) {
00666     Zoltan_Destroy( &m_zoltan_id_ );
00667   }
00668   m_zoltan_id_ = NULL ;
00669   if(Name_Conversion) {
00670     delete Name_Conversion;
00671     Name_Conversion = NULL;
00672   }
00673   if(Value_Conversion) {
00674     delete Value_Conversion;
00675     Value_Conversion = NULL;
00676   }
00677 }
00678 
00679 void
00680 Zoltan::reset_dest_proc_data()
00681 {
00682   const int  proc = 0; //Env::parallel_rank();
00683   const unsigned size = m_mesh_information_.mesh_entities.size();
00684   m_mesh_information_.dest_proc_ids.assign(size, proc);
00685 }
00686 
00687 void
00688 Zoltan::set_destination_proc(const unsigned moid,
00689                              const unsigned proc )
00690 {
00691   m_mesh_information_.dest_proc_ids[ moid ] = proc;
00692 }
00693 
00694 unsigned
00695 Zoltan::destination_proc(const unsigned moid) const
00696 {
00697   return m_mesh_information_.dest_proc_ids[ moid ];
00698 }
00699 
00700 bool
00701 Zoltan::find_mesh_entity(const mesh::Entity * entity, unsigned & moid) const
00702 {
00703   unsigned len = m_mesh_information_.mesh_entities.size();
00704   for(moid = 0; moid < len; ++moid)
00705   {
00706     if(m_mesh_information_.mesh_entities[moid] == entity) return true;
00707   }
00708   return false;
00709 }
00710 
00711 int
00712 Zoltan::get_new_partition(stk::mesh::EntityProcVec &rebal_spec)
00713 {
00714   const unsigned entity_iter_len = m_mesh_information_.mesh_entities.size();
00715   for (unsigned entity_iter =0; entity_iter != entity_iter_len; ++entity_iter) {
00716     mesh::Entity * mesh_ent = mesh_entity(entity_iter);
00717     int proc = destination_proc(entity_iter);
00718     mesh::EntityProc et(mesh_ent, proc);
00719     rebal_spec.push_back(et);
00720   }
00721   return 0;
00722 }
00723 
00724 const VectorField *
00725 Zoltan::entity_coord_ref() const
00726 {
00727    return m_mesh_information_.nodal_coord_ref;
00728 }
00729 
00730 mesh::Entity *
00731 Zoltan::mesh_entity(const unsigned moid ) const
00732 {
00733   return m_mesh_information_.mesh_entities[ moid ];
00734 }
00735 
00736 double
00737 Zoltan::entity_weight(const unsigned moid ) const
00738 {
00739   double mo_weight = 1.0;
00740   if (entity_weight_ref()) {
00741     mo_weight = * static_cast<double *>
00742       ( mesh::field_data (*entity_weight_ref(), *m_mesh_information_.mesh_entities[ moid ]));
00743   }
00744   return mo_weight;
00745 }
00746 
00747 const ScalarField *
00748 Zoltan::entity_weight_ref() const
00749 {
00750   return m_mesh_information_.elem_weight_ref;
00751 }
00752 
00753 unsigned
00754 Zoltan::num_moid() const
00755 {
00756   return m_mesh_information_.mesh_entities.size() ;
00757 }
00758 
00759 const std::string &Zoltan::parameter_entry_name() const
00760 {
00761   return m_parameter_entry_name_;
00762 }
00763 
00764 
00765 //: Load Balance calls or Load "Re-partitioning" calls
00766 
00767 int Zoltan::register_callbacks()
00768 {
00769   /*
00770    * Register the Zoltan/SIERRA "call-back" (querry) functions.
00771    * Use ONLY THE BARE ESSENTIALS for decompositions:
00772    *
00773    *    Zoltan_Set_Num_Obj_Fn      
00774    *    Zoltan_Set_Obj_List_Fn     
00775    *    Zoltan_Set_Num_Geom_Fn
00776    *    Zoltan_Set_Geom_Fn
00777    */
00778 
00779   /*
00780    * flag for what data is to be registered with the zoltan callback
00781    * functions in combination with the "static_cast<CLASS>(data)"
00782    * statement used in the zoltan interface routine,
00783    * Fmwk_Zoltaninterface.C
00784    *
00785    */
00786 
00787   if ( Zoltan_Set_Num_Obj_Fn( m_zoltan_id_,
00788                               Callback_Num_Elements,
00789                               this )
00790        != ZOLTAN_OK ) {
00791     throw std::runtime_error("Zoltan_Set_Num_Obj_Fn using Callback_Num_Elements failed to register");
00792   }
00793   if ( Zoltan_Set_Obj_List_Fn( m_zoltan_id_, Callback_Element_List,this )
00794        != ZOLTAN_OK ) {
00795     throw std::runtime_error("Zoltan_Set_Obj_List_Fn using Callback_Element_List");
00796   }
00797   if ( Zoltan_Set_Num_Geom_Fn( m_zoltan_id_, Callback_Num_Dimensions,this )
00798        != ZOLTAN_OK ) {
00799     throw std::runtime_error("Zoltan_Set_Num_Geom_Fn using Callback_Num_Dimensions");
00800   }
00801   if ( Zoltan_Set_Geom_Fn( m_zoltan_id_, Callback_Centroid_Coord,this )
00802        != ZOLTAN_OK ) {
00803     throw std::runtime_error("Zoltan_Set_Geom_Fn using Callback_Centroid_Coord");
00804   }
00805   if ( Zoltan_Set_Num_Edges_Fn( m_zoltan_id_, Callback_Num_Edges,this )
00806        != ZOLTAN_OK ) {
00807     throw std::runtime_error("Zoltan_Set_Num_Edges_Fn using Callback_Num_Edges");
00808   }
00809   if ( Zoltan_Set_Edge_List_Fn( m_zoltan_id_, Callback_Edge_List,this )
00810        != ZOLTAN_OK ) {
00811     throw std::runtime_error("Zoltan_Set_Edge_List_Fn using Callback_Edge_List");
00812   }
00813 
00814   return 0;
00815 
00816 }
00817 
00818 
00819 int  Zoltan::evaluate( int    print_stats,
00820                        int*   nentity,
00821                        double*  entity_wgt,
00822                        int*   ncuts,
00823                        double*  cut_wgt,
00824                        int*   nboundary,
00825                        int*   nadj      )
00826 {
00827   int ierr        = 0;
00828 
00829   ZOLTAN_BALANCE_EVAL eval  = {0};
00830   ZOLTAN_GRAPH_EVAL   graph = {{0}};
00831   if (Zoltan_LB_Eval_Balance( m_zoltan_id_, print_stats, &eval)) ierr = 1;
00832   if (Zoltan_LB_Eval_Graph( m_zoltan_id_, print_stats, &graph) ) ierr = 1;
00833   *nentity      = (int)eval.nobj[0];
00834   *entity_wgt   =      eval.obj_wgt[0];
00835   *ncuts        = (int)graph.cuts[0];
00836   *cut_wgt      =      graph.cut_wgt[0];
00837   *nboundary    = (int)graph.num_boundary[0];
00838   *nadj         = (int)graph.nnborparts[0];
00839 
00840   return ierr;
00841 
00842 }
00843 
00844 void Zoltan::determine_new_partition (bool &RebalancingNeeded)
00845 {
00846   //: Transfer export global ID lists and owning processors
00847   //: to SIERRA Framework's data structures
00848 
00855   int length_gid, length_lid;
00856 
00857   int           new_decomp;
00858   int           num_imported;
00859   ZOLTAN_ID_PTR import_gids;
00860   ZOLTAN_ID_PTR import_lids;
00861   int          *import_procs=NULL;
00862   int           num_exported;
00863   ZOLTAN_ID_PTR export_gids;
00864   ZOLTAN_ID_PTR export_lids;
00865   int          *export_procs=NULL;
00866 
00867 
00868 
00924   int status = Zoltan_LB_Balance( m_zoltan_id_,        &new_decomp,
00925                                   &length_gid     , &length_lid,
00926                                   &num_imported,    &import_gids,
00927                                   &import_lids,     &import_procs,
00928                                   &num_exported,    &export_gids,
00929                                   &export_lids,     &export_procs );
00930   if (status != ZOLTAN_OK) {
00931     throw std::runtime_error("Zoltan_Balance() returned error code " + status);
00932   }
00933 
00934   //: Initialize destination processor IDs (dest_proc_ids)
00935   reset_dest_proc_data();
00936 
00937   int actual_exported = 0;
00938   if ( new_decomp && ( num_exported != -1 ) ) {
00939     const unsigned parallel_rank = parallel_machine_rank(comm_);
00940     /* New Decomposition was generated */
00941     for (int j=0; j < num_exported; ++j ) {
00942 
00943       //: Get exported region, local, global, and processor ids
00944       //const unsigned rid = export_gids[ j*::num_gid_entries() ];  // Region ID variable
00945       //if (!rid) throw runtime_error ("Region ID variable should be non-zero.");
00946       const unsigned lid = export_lids[ j*::num_lid_entries() ];  // Local  ID variable
00947       const unsigned pid = export_procs[ j ];                     // Exported Processor ID (i.e., MPI "rank" )
00948 
00949       if (parallel_rank != pid) {
00950         ++actual_exported;
00951         set_destination_proc(lid, pid);
00952       }
00953     }
00954   }
00955 
00956   RebalancingNeeded = 0 ;
00957   if (new_decomp) {
00958     int rebalneeded=0;
00959     stk::all_reduce_sum(comm_, &actual_exported, &rebalneeded, 1);
00960     if (rebalneeded)  RebalancingNeeded = 1;
00961   }
00962 
00966   if ( ZOLTAN_OK !=
00967        Zoltan_LB_Free_Data( &import_gids, &import_lids, &import_procs,
00968                             &export_gids, &export_lids, &export_procs )) {
00969     throw runtime_error (" FATAL ERROR in Zoltan_LB_Free_Data.");
00970   }
00971 
00972 }
00973 
00974 void Zoltan::convert_names_and_values(const Parameters &from, Parameters &to)
00975 {
00976   /* First time through, fill the conversion tables. */
00977   if (!Name_Conversion) {
00978     Name_Conversion = new Parameters;
00979     fill_name_conversion (*Name_Conversion);
00980   }
00981   if (!Value_Conversion) {
00982     Value_Conversion = new Parameters;
00983     fill_value_conversion(*Value_Conversion);
00984   }
00985 
00986   // NOTE: "LOAD BALANCING METHOD" default
00987   // is also hard coded in fill_default_values();
00988   std::string algorithm;
00989   const std::string keyname("LOAD BALANCING METHOD");
00990   if( from.isParameter(keyname) )
00991     algorithm = from.get<std::string>(keyname);
00992   else
00993     algorithm = "0";
00994 
00995   const Parameters & General = Name_Conversion->sublist("General");
00996   const Parameters & Algorithm = Name_Conversion->sublist(algorithm);
00997 
00998   Parameters::ConstIterator
00999     from_iter  = from.begin(),
01000     from_end   = from.end();
01001 
01002   /* Iterate over all of the input parameters to find proper
01003      Zoltan names and Zoltan parameter values. */
01004   for (; from_iter != from_end; ++from_iter) {
01005     std::string
01006       from_name = from.name(from_iter),
01007       to_name;
01008 
01009     /* Check to see if this is a general parameter name
01010        and if not, check if it is an algorithm specific name. */
01011 
01012     if (General.isParameter(from_name)) {
01013       to_name = General.get<std::string>(from_name);
01014     } else {
01015       to_name = Algorithm.get<std::string>(from_name);
01016     }
01017 
01018     /* Now convert the parameter value to the correct form.
01019        Only a couple of parameters have parameter conversion.
01020        The ones converted are nested in Value_Conversion.
01021     */
01022     std::string to_value = Teuchos::getValue<string>(from.entry(from_iter));
01023     //if (Value_Conversion->isParameter(from_name)) to_value = Value_Conversion->get<std::string>(to_value);
01024     if (Value_Conversion->isParameter(from_name)) to_value = Value_Conversion->sublist(from_name).get<std::string>(to_value);
01025     if (!to_name.empty()) to.set(to_name, to_value);
01026   }
01027 }
01028 
01029 void Zoltan::merge_default_values(const Parameters &from,
01030                                         Parameters &to)
01031 {
01032   Parameters default_values;
01033   fill_default_values(default_values);
01034   to.setParameters(default_values);
01035   to.setParameters(from);
01036 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends