AFEPack
|
00001 00053 #ifndef __MPI_ULoadBalance_h__ 00054 #define __MPI_ULoadBalance_h__ 00055 00056 #ifdef __MPI_LoadBalance_h__ 00057 #error "MPI_LoadBalance.h 和 MPI_ULoadBalance.h 是冲突的,不能一起用。" 00058 #endif 00059 00060 #include <sys/types.h> 00061 #include <sys/stat.h> 00062 #include <unistd.h> 00063 #include <cstdio> 00064 00065 #include <boost/program_options.hpp> 00066 #include <boost/archive/binary_iarchive.hpp> 00067 #include <boost/archive/binary_oarchive.hpp> 00068 #include <boost/archive/text_iarchive.hpp> 00069 #include <boost/archive/text_oarchive.hpp> 00070 00071 #include <AFEPack/Serialization.h> 00072 #include "MPI.h" 00073 #include "MPI_UGeometry_archive.h" 00074 #include "MPI_Migration.h" 00075 00076 namespace MPI { 00077 00078 template <class FOREST> 00079 class HLoadBalance { 00080 public: 00081 typedef FOREST forest_t; 00082 enum {dim = FOREST::dim, dow = FOREST::dow}; 00083 typedef typename forest_t::matcher_t matcher_t; 00084 00085 struct rank_map_val { 00086 int old_rank; 00087 int priority; 00092 }; 00093 typedef std::map<int,rank_map_val> rank_map_t; 00094 00095 template <class GEO> bool 00096 is_refined_on_new_rank(const GEO& geo, 00097 int new_rank) const { 00098 bool result = false; 00099 if (geo.isRefined() && 00100 is_on_this_new_rank(*geo.child[0], new_rank)) { 00101 result = true; 00102 } 00103 return result; 00104 } 00105 00106 private: 00107 // 在计算每个几何体上的负载是使用的性质 00108 property_id_t<u_int> _pid_loading; 00109 00110 // 标记一个共享几何体是否已经存储了的性质 00111 property_id_t<> _pid_is_saved; 00112 public: 00113 template <class GEO> void 00114 set_is_saved(const GEO& geo) const { 00115 if (geo.get_property(_pid_is_saved) == NULL) { 00116 geo.new_property(_pid_is_saved); 00117 } 00118 } 00119 00120 private: 00125 property_id_t<rank_map_t> _pid_rank_map; 00126 00127 // 部分需要拆分和合并的几何体的全局标号 00128 property_id_t<unsigned long> _pid_global_idx; 00129 00130 // 元素为(全局标号,对象指针)的对 00131 std::map<unsigned long, void *> _global_pointer_to_merge; 00132 std::map<unsigned long, std::list<void *> > _global_pointer_to_merge_later; 00133 00134 // (进程,列表(全局标号,维数,对象指针)组)的映射 00135 typedef std::map<unsigned long, std::pair<int, void *> > tuple_map_t; 00136 std::map<int, tuple_map_t> _global_pointer_to_share; 00137 00138 typedef BirdView<forest_t> birdview_t; 00139 typedef BirdViewSet<forest_t> birdview_set_t; 00140 typedef HLoadBalance<forest_t> this_t; 00141 00142 forest_t * _forest; 00143 std::vector<u_int> _cut_point; 00144 std::vector<int> _new_rank; 00145 00146 typename forest_t::container_t _nre; 00147 00148 public: 00149 HLoadBalance(forest_t& forest) : _forest(&forest) {} 00150 ~HLoadBalance() {} 00151 00152 forest_t& forest() const { return *_forest; } 00153 int rank() const { return _forest->rank(); } 00154 00155 private: 00156 template <class GEO> rank_map_t * 00157 new_rank_map(const GEO& geo) const { 00158 return geo.new_property(_pid_rank_map); 00159 } 00160 template <class GEO> unsigned long * 00161 new_global_idx(const GEO& geo) const { 00162 return geo.new_property(_pid_global_idx); 00163 } 00164 00165 public: 00166 void clear() { 00167 free_property_id(_pid_loading); 00168 free_property_id(_pid_rank_map); 00169 free_property_id(_pid_global_idx); 00170 _global_pointer_to_merge.clear(); 00171 _global_pointer_to_merge_later.clear(); 00172 _global_pointer_to_share.clear(); 00173 _cut_point.clear(); 00174 _new_rank.clear(); 00175 } 00176 00177 template <class GEO> rank_map_t * 00178 get_rank_map(const GEO& geo) const { 00179 return geo.get_property(_pid_rank_map); 00180 } 00181 template <class GEO> unsigned long * 00182 get_global_idx(const GEO& geo) const { 00183 return geo.get_property(_pid_global_idx); 00184 } 00185 00187 00190 template <class GEO> void 00191 print_geo_rank_map(const GEO& geo) const { 00192 std::cerr << "[" << _forest->rank() << ":" << geo.dimension; 00193 rank_map_t * p_map = get_rank_map(geo); 00194 if (p_map != NULL) print_rank_map(p_map); 00195 std::cerr << "]"; 00196 } 00197 00198 void print_rank_map(const rank_map_t * p_map) const { 00199 typename rank_map_t::const_iterator 00200 the_pair = p_map->begin(), 00201 end_pair = p_map->end(); 00202 for (;the_pair != end_pair;++ the_pair) { 00203 std::cerr << "(" << the_pair->first << ":" 00204 << the_pair->second.old_rank << ":" 00205 << the_pair->second.priority << ")"; 00206 } 00207 } 00209 00213 template <class GEO> bool 00214 is_on_this_new_rank(const GEO& geo, 00215 int new_rank) const { 00216 rank_map_t * p_map = this->get_rank_map(geo); 00217 if (p_map == NULL) return false; 00218 typename rank_map_t::iterator 00219 the_pair = p_map->find(new_rank); 00220 return (the_pair != p_map->end()); 00221 } 00222 00231 template <class GEO> bool 00232 is_save_on_this_rank(const GEO& geo, 00233 int new_rank) const { 00234 rank_map_t * p_map = this->get_rank_map(geo); 00235 assert (p_map != NULL); 00236 typename rank_map_t::iterator 00237 the_pair = p_map->find(new_rank); 00239 if (the_pair == p_map->end()) return false; 00240 00242 return (the_pair->second.old_rank == _forest->rank()); 00243 } 00244 00256 template <class GEO> 00257 void merge_global_pointer(int type, 00258 unsigned long global_idx, 00259 GEO *& p_geo) { 00260 std::map<unsigned long, void *>::iterator 00261 the_pair = _global_pointer_to_merge.find(global_idx); 00262 if (type == 2) { 00263 if (the_pair == _global_pointer_to_merge.end()) { 00264 _global_pointer_to_merge.insert(std::pair<unsigned long,void*>(global_idx, p_geo)); 00265 00269 std::map<unsigned long, std::list<void *> >::iterator 00270 the_entry = _global_pointer_to_merge_later.find(global_idx); 00271 if (the_entry != _global_pointer_to_merge_later.end()) { 00272 std::list<void *>::iterator 00273 the_ptr = the_entry->second.begin(), 00274 end_ptr = the_entry->second.end(); 00275 for (;the_ptr != end_ptr;++ the_ptr) { 00276 GEO ** geo_ptr_ptr = (GEO **)(*the_ptr); 00277 (*geo_ptr_ptr) = p_geo; 00278 } 00279 _global_pointer_to_merge_later.erase(the_entry); 00280 } 00281 } else { 00282 assert (the_pair->second == p_geo); 00283 } 00284 } else { 00285 assert (type == 4); 00286 if (the_pair != _global_pointer_to_merge.end()) { 00287 p_geo = (GEO *)(the_pair->second); 00288 } else { 00289 _global_pointer_to_merge_later[global_idx].push_back((void *)(&p_geo)); 00290 } 00291 } 00292 } 00293 00302 template <class GEO> 00303 void share_global_pointer(int rank, 00304 unsigned long global_idx, 00305 GEO *& p_geo) { 00306 _global_pointer_to_share[rank][global_idx] = 00307 std::pair<int,void *>(p_geo->dimension, (void *)(p_geo)); 00308 } 00309 00310 private: 00314 void share_global_pointer() { 00315 std::cerr << "Exchanging shared global pointers ..." << std::endl; 00316 00317 typedef std::map<int, tuple_map_t> map_t; 00318 00319 int n = 0; 00320 std::list<int> target; 00321 std::list<BinaryBuffer<> > out_buf, in_buf; 00322 typename map_t::iterator 00323 the_pair = _global_pointer_to_share.begin(), 00324 end_pair = _global_pointer_to_share.end(); 00325 for (;the_pair != end_pair;++ the_pair, ++ n) { 00326 target.push_back(the_pair->first); 00327 00328 out_buf.push_back(BinaryBuffer<>()); 00329 AFEPack::ostream<> os(out_buf.back()); 00330 00331 in_buf.push_back(BinaryBuffer<>()); 00332 00333 u_int n_item = the_pair->second.size(); 00334 os << n_item; 00335 #if 0 00336 std::cerr << "sending " << the_pair->second.size() 00337 << " items from " << _forest->rank() 00338 << " to " << the_pair->first << std::endl; 00339 #endif 00340 typename tuple_map_t::iterator 00341 the_tuple = the_pair->second.begin(), 00342 end_tuple = the_pair->second.end(); 00343 for (;the_tuple != end_tuple;++ the_tuple) { 00344 const unsigned long& global_idx = the_tuple->first; 00345 int& dimension = the_tuple->second.first; 00346 void *& p_obj = the_tuple->second.second; 00347 os << global_idx << dimension << p_obj; 00348 } 00349 } 00350 00351 MPI_Barrier(_forest->communicator()); 00352 sendrecv_data(_forest->communicator(), n, out_buf.begin(), 00353 in_buf.begin(), target.begin()); 00354 00355 typename std::list<int>::iterator 00356 the_rank = target.begin(), 00357 end_rank = target.end(); 00358 typename std::list<BinaryBuffer<> >::iterator 00359 the_buf = in_buf.begin(); 00360 for (;the_rank != end_rank;++ the_rank, ++ the_buf) { 00361 AFEPack::istream<> is(*the_buf); 00362 u_int n_item; 00363 is >> n_item; 00364 #if 0 00365 std::cerr << "received " << n_item 00366 << " items from " << *the_rank 00367 << " to " << _forest->rank() << std::endl; 00368 #endif 00369 for (u_int i = 0;i < n_item;++ i) { 00370 unsigned long global_idx; 00371 int dimension; 00372 void * remote_obj; 00373 is >> global_idx >> dimension >> remote_obj; 00374 std::pair<int,void *>& the_pair = 00375 _global_pointer_to_share[*the_rank][global_idx]; 00376 assert ((dimension == the_pair.first) && 00377 (dimension >= 0) && 00378 (dimension <= 3)); 00379 00380 if (dimension == 0) { 00381 #define GDIM 0 00382 #define GEO HGeometry<GDIM,dow> 00383 #define OBJ Shared_object<GEO> 00384 GEO * p_geo = (GEO *)(the_pair.second); 00385 OBJ * p_info = _forest->get_shared_info(*p_geo); 00386 if (p_info == NULL) p_info = _forest->new_shared_info(*p_geo); 00387 p_info->add_clone(*the_rank, (GEO *)remote_obj); 00388 #undef OBJ 00389 #undef GEO 00390 #undef GDIM 00391 } else if (dimension == 1) { 00392 #define GDIM 1 00393 #define GEO HGeometry<GDIM,dow> 00394 #define OBJ Shared_object<GEO> 00395 GEO * p_geo = (GEO *)(the_pair.second); 00396 OBJ * p_info = _forest->get_shared_info(*p_geo); 00397 if (p_info == NULL) p_info = _forest->new_shared_info(*p_geo); 00398 p_info->add_clone(*the_rank, (GEO *)remote_obj); 00399 #undef OBJ 00400 #undef GEO 00401 #undef GDIM 00402 } else if (dimension == 2) { 00403 #define GDIM 2 00404 #define GEO HGeometry<GDIM,dow> 00405 #define OBJ Shared_object<GEO> 00406 GEO * p_geo = (GEO *)(the_pair.second); 00407 OBJ * p_info = _forest->get_shared_info(*p_geo); 00408 if (p_info == NULL) p_info = _forest->new_shared_info(*p_geo); 00409 p_info->add_clone(*the_rank, (GEO *)remote_obj); 00410 #undef OBJ 00411 #undef GEO 00412 #undef GDIM 00413 } else if (dimension == 3) { 00414 #define GDIM 3 00415 #define GEO HGeometry<GDIM,dow> 00416 #define OBJ Shared_object<GEO> 00417 GEO * p_geo = (GEO *)(the_pair.second); 00418 OBJ * p_info = _forest->get_shared_info(*p_geo); 00419 if (p_info == NULL) p_info = _forest->new_shared_info(*p_geo); 00420 p_info->add_clone(*the_rank, (GEO *)remote_obj); 00421 #undef OBJ 00422 #undef GEO 00423 #undef GDIM 00424 } 00425 } 00426 } 00427 00428 #if 1 00429 00432 the_pair = _global_pointer_to_share.begin(); 00433 for (;the_pair != end_pair;++ the_pair) { 00434 typename tuple_map_t::iterator 00435 the_tuple = the_pair->second.begin(), 00436 end_tuple = the_pair->second.end(); 00437 for (;the_tuple != end_tuple;++ the_tuple) { 00438 int& dimension = the_tuple->second.first; 00439 void *& p_obj = the_tuple->second.second; 00440 if (dimension == 0) { 00441 #define GDIM 0 00442 #define GEO HGeometry<GDIM,dow> 00443 GEO * p_geo = (GEO *)p_obj; 00444 assert (_forest->get_shared_info(*p_geo) != NULL); 00445 #undef GEO 00446 #undef GDIM 00447 } else if (dimension == 1) { 00448 #define GDIM 1 00449 #define GEO HGeometry<GDIM,dow> 00450 GEO * p_geo = (GEO *)p_obj; 00451 assert (_forest->get_shared_info(*p_geo) != NULL); 00452 #undef GEO 00453 #undef GDIM 00454 } else if (dimension == 2) { 00455 #define GDIM 2 00456 #define GEO HGeometry<GDIM,dow> 00457 GEO * p_geo = (GEO *)p_obj; 00458 assert (_forest->get_shared_info(*p_geo) != NULL); 00459 #undef GEO 00460 #undef GDIM 00461 } else if (dimension == 3) { 00462 #define GDIM 3 00463 #define GEO HGeometry<GDIM,dow> 00464 GEO * p_geo = (GEO *)p_obj; 00465 assert (_forest->get_shared_info(*p_geo) != NULL); 00466 #undef GEO 00467 #undef GDIM 00468 } 00469 } 00470 } 00471 #endif 00472 } 00473 00475 public: 00481 template <class LOADER> 00482 void config(birdview_t& ir_mesh, 00483 LOADER& loader, 00484 u_int (LOADER::*value)(const GeometryBM&)) { 00485 if (&(ir_mesh.getForest()) != _forest) { 00486 std::cerr << "The mesh should be on the same forest of the HLoadBalance." 00487 << std::endl; 00488 abort(); 00489 } 00490 00491 RegularMesh<dim,dow>& mesh = ir_mesh.regularMesh(); 00492 new_property_id(_pid_loading); 00493 u_int n_ele = mesh.n_geometry(dim); 00494 for (u_int i = 0;i < n_ele;++ i) { 00495 HGeometry<dim,dow> * p_geo = mesh.template h_geometry<dim>(i); 00496 u_int * p_loading = p_geo->get_property(_pid_loading); 00497 if (p_loading == NULL) { 00498 p_loading = p_geo->new_property(_pid_loading); 00499 } 00500 (*p_loading) = (loader.*value)(mesh.geometry(dim,i)); 00501 } 00502 } 00503 00504 void config(birdview_t& ir_mesh, 00505 u_int (*value)(const GeometryBM&) = &_default_element_loading) { 00506 if (&(ir_mesh.getForest()) != _forest) { 00507 std::cerr << "The mesh should be on the same forest of the HLoadBalance." 00508 << std::endl; 00509 abort(); 00510 } 00511 00512 RegularMesh<dim,dow>& mesh = ir_mesh.regularMesh(); 00513 new_property_id(_pid_loading); 00514 u_int n_ele = mesh.n_geometry(dim); 00515 for (u_int i = 0;i < n_ele;++ i) { 00516 const HGeometry<dim,dow> * p_geo = mesh.template h_geometry<dim>(i); 00517 u_int * p_loading = p_geo->get_property(_pid_loading); 00518 if (p_loading == NULL) { 00519 p_loading = p_geo->new_property(_pid_loading); 00520 } 00521 (*p_loading) = (*value)(mesh.geometry(dim,i)); 00522 } 00523 } 00524 00525 static u_int _default_element_loading(const GeometryBM&) { return 1; } 00527 00529 public: 00537 void partition(u_int n_new_rank = 0, 00538 double tol_percent = 0.05, 00539 bool is_renumerate_root = false, 00540 void (*f)(double, double, double, 00541 double&,double&,double&) = NULL) { 00542 int rank = _forest->rank(); 00543 int n_rank = _forest->n_rank(); 00544 00545 u_int loading = this->lump_loading(); 00546 u_int total_loading, partial_loading; 00547 MPI_Comm comm = _forest->communicator(); 00548 MPI_Scan(&loading, &partial_loading, 1, MPI_UNSIGNED, MPI_SUM, comm); 00549 00550 if (rank == n_rank - 1) { 00551 total_loading = partial_loading; 00552 } 00553 if (n_new_rank == 0) n_new_rank = n_rank; 00554 00555 MPI_Bcast(&total_loading, 1, MPI_UNSIGNED, n_rank - 1, comm); 00556 double mean_loading = n_new_rank; 00557 mean_loading = total_loading/mean_loading; 00558 00562 this->refine_root_element(mean_loading*tol_percent, is_renumerate_root, f); 00563 00564 _cut_point.clear(); 00565 _new_rank.clear(); 00566 00567 u_int idx = 0; 00568 _cut_point.push_back(idx); 00569 00570 loading = partial_loading - loading; 00571 u_int current_rank = (u_int)(loading/mean_loading); 00572 _new_rank.push_back(current_rank); 00573 00574 typename forest_t::RootIterator 00575 the_ele = _forest->beginRootElement(), 00576 end_ele = _forest->endRootElement(); 00577 for (;the_ele != end_ele;++ the_ele) { 00578 idx += 1; 00579 loading += this->get_loading(*the_ele); 00580 u_int the_rank = (u_int)(loading/mean_loading); 00581 the_rank = std::min(the_rank, n_new_rank - 1); 00582 if (the_rank > current_rank) { 00583 _cut_point.push_back(idx); 00584 _new_rank.push_back(++ current_rank); 00585 } 00586 } 00587 if (_cut_point.back() != idx) { 00588 _cut_point.push_back(idx); 00589 } else { 00590 _new_rank.pop_back(); 00591 } 00592 00593 free_property_id(_pid_loading); 00594 } 00595 00596 private: 00601 u_int lump_loading() { 00602 u_int loading = 0; 00603 typename forest_t::RootIterator 00604 the_ele = _forest->beginRootElement(), 00605 end_ele = _forest->endRootElement(); 00606 for (;the_ele != end_ele;++ the_ele) { 00607 loading += this->get_loading(*the_ele); 00608 } 00609 return loading; 00610 } 00611 00612 private: 00613 00618 template <class GEO> u_int 00619 get_loading(const GEO& geo) const { 00620 if (&geo == NULL) return 0; 00621 u_int * p_loading = geo.get_property(_pid_loading); 00622 if (p_loading == NULL) { 00623 p_loading = geo.new_property(_pid_loading); 00624 (*p_loading) = 0; 00625 for (u_int i = 0;i < geo.n_child;++ i) { 00626 (*p_loading) += get_loading(*geo.child[i]); 00627 } 00628 } 00629 return (*p_loading); 00630 } 00632 00634 00644 void refine_root_element(double tol, bool is_renum_root, 00645 void (*f)(double,double,double, 00646 double&,double&,double&)) { 00647 _nre.clear(); 00648 00649 typename forest_t::RootIterator 00650 the_ele = _forest->beginRootElement(), 00651 end_ele = _forest->endRootElement(); 00652 for (;the_ele != end_ele;++ the_ele) { 00653 this->refine_root_element(tol, *the_ele, _nre); 00654 } 00655 _forest->rootElement().swap(_nre); 00656 if (is_renum_root) { 00657 _forest->renumerateRootElement(f); 00658 } 00659 } 00660 00661 void refine_root_element(double tol, 00662 HGeometry<dim,dow>& geo, 00663 typename forest_t::container_t& nre) const { 00664 u_int * p_loading = geo.get_property(_pid_loading); 00665 u_int loading = 0; 00666 if (p_loading != NULL) loading = *p_loading; 00667 if ((loading > tol) && 00668 (geo.isRefined()) && 00669 (geo.child[0]->get_property(_pid_loading) != NULL)) { 00670 for (u_int i = 0;i < geo.n_child;++ i) { 00671 this->refine_root_element(tol, *geo.child[i], nre); 00672 } 00673 } else { 00674 nre.push_back(&geo); 00675 } 00676 } 00677 00682 void coarse_root_element() { 00683 property_id_t<bool> pid_is_found; 00684 new_property_id(pid_is_found); 00685 00686 _nre.clear(); 00687 typename forest_t::RootIterator 00688 the_ele = _forest->beginRootElement(), 00689 end_ele = _forest->endRootElement(); 00690 for (;the_ele != end_ele;++ the_ele) { 00691 this->coarse_root_element(pid_is_found, *the_ele, _nre); 00692 } 00693 _forest->rootElement().swap(_nre); 00694 } 00695 00696 void coarse_root_element(const property_id_t<bool>& pid, 00697 HGeometry<dim,dow>& geo, 00698 typename forest_t::container_t& nre) const { 00699 if (geo.get_property(pid) != NULL) return; 00700 if ((geo.parent != NULL) && 00701 (_forest->get_shared_info(*geo.parent) == NULL)) { 00702 this->coarse_root_element(pid, *geo.parent, nre); 00703 } else { 00704 nre.push_back(&geo); 00705 geo.new_property(pid); 00706 } 00707 } 00709 00710 00712 private: 00718 void set_new_rank() { 00719 new_property_id(_pid_rank_map); 00720 00721 set_new_rank_down(); 00722 set_new_rank_up(); 00723 00724 if (_forest->n_rank() <= 1) return; 00725 00729 Shared_type_filter::only<0> type_filter; 00730 MPI_Comm comm = _forest->communicator(); 00731 #define SYNC_DATA(D) \ 00732 if (dim >= D) { \ 00733 sync_data(comm, _forest->template get_shared_list<D>(), *this, \ 00734 &this_t::template pack_set_new_rank<D>, \ 00735 &this_t::template unpack_set_new_rank<D>, \ 00736 type_filter); \ 00737 } 00738 SYNC_DATA(0); 00739 SYNC_DATA(1); 00740 SYNC_DATA(2); 00741 SYNC_DATA(3); 00742 #undef SYNC_DATA 00743 00744 //check_new_rank(); /// 调试时用 00745 } 00746 00747 public: 00748 template <int GDIM> void 00749 pack_set_new_rank(HGeometry<GDIM,dow> * geo, 00750 int remote_rank, 00751 AFEPack::ostream<>& os) { 00752 assert ((_forest->get_shared_info(*geo) != NULL)); 00753 00754 rank_map_t * p_map = get_rank_map(*geo); 00755 if (p_map == NULL) { 00756 p_map = new_rank_map(*geo); 00757 } 00758 00759 std::size_t n = p_map->size(); 00760 os << n; 00761 typename rank_map_t::iterator 00762 the_pair = p_map->begin(), 00763 end_pair = p_map->end(); 00764 for (;the_pair != end_pair;++ the_pair) { 00765 const int& new_rank = the_pair->first; 00766 const int& old_rank = the_pair->second.old_rank; 00767 int& priority = the_pair->second.priority; 00768 if (is_refined_on_new_rank(*geo, new_rank)) { 00769 priority = 2; 00770 } else if (_forest->is_dummy(*geo)) { 00771 priority = 0; 00772 } else { 00773 priority = 1; 00774 } 00775 00776 os << new_rank << old_rank << priority; 00777 } 00778 } 00779 00780 template <int GDIM> void 00781 unpack_set_new_rank(HGeometry<GDIM,dow> * geo, 00782 int remote_rank, 00783 AFEPack::istream<>& is) { 00784 assert ((_forest->get_shared_info(*geo) != NULL)); 00785 00794 rank_map_t * p_map = get_rank_map(*geo); 00795 assert (p_map != NULL); 00796 00797 int new_rank; 00798 std::size_t n; 00799 is >> n; 00800 for (u_int i = 0;i < n;++ i) { 00801 rank_map_val remote_rank_map; 00802 is >> new_rank 00803 >> remote_rank_map.old_rank 00804 >> remote_rank_map.priority; 00805 typename rank_map_t::iterator 00806 the_pair = p_map->find(new_rank); 00807 if (the_pair == p_map->end()) { 00808 (*p_map)[new_rank] = remote_rank_map; 00809 } else { 00810 if (the_pair->second.priority == remote_rank_map.priority) { 00812 if (the_pair->second.old_rank > remote_rank_map.old_rank) { 00813 the_pair->second.old_rank = remote_rank_map.old_rank; 00814 } 00815 } else if (remote_rank_map.priority > the_pair->second.priority) { 00816 the_pair->second = remote_rank_map; 00817 } 00818 } 00819 } 00820 } 00821 00822 private: 00827 void check_new_rank() { 00828 Shared_type_filter::only<0> type_filter; 00829 MPI_Comm comm = _forest->communicator(); 00830 00831 #define SYNC_DATA(D) \ 00832 if (dim >= D) { \ 00833 sync_data(comm, _forest->template get_shared_list<D>(), *this, \ 00834 &this_t::template pack_check_new_rank<D>, \ 00835 &this_t::template unpack_check_new_rank<D>, \ 00836 type_filter); \ 00837 } 00838 SYNC_DATA(0); 00839 SYNC_DATA(1); 00840 SYNC_DATA(2); 00841 SYNC_DATA(3); 00842 #undef SYNC_DATA 00843 } 00844 00845 public: 00846 template <int GDIM> void 00847 pack_check_new_rank(HGeometry<GDIM,dow> * geo, 00848 int remote_rank, 00849 AFEPack::ostream<>& os) { 00850 rank_map_t * p_map = get_rank_map(*geo); 00851 assert (p_map != NULL); 00852 00853 os << p_map->size(); 00854 typename rank_map_t::iterator 00855 the_pair = p_map->begin(), 00856 end_pair = p_map->end(); 00857 for (;the_pair != end_pair;++ the_pair) { 00858 const int& new_rank = the_pair->first; 00859 const int& old_rank = the_pair->second.old_rank; 00860 int& priority = the_pair->second.priority; 00861 os << new_rank << old_rank << priority; 00862 } 00863 } 00864 00865 template <int GDIM> void 00866 unpack_check_new_rank(HGeometry<GDIM,dow> * geo, 00867 int remote_rank, 00868 AFEPack::istream<>& is) { 00869 assert ((_forest->get_shared_info(*geo) != NULL)); 00870 00871 rank_map_t * p_map = get_rank_map(*geo); 00872 assert (p_map != NULL); 00873 00874 std::size_t n; 00875 is >> n; 00876 00877 rank_map_t remote_rank_map; 00878 for (u_int i = 0;i < n;++ i) { 00879 int new_rank; 00880 rank_map_val remote_rank_map_val; 00881 is >> new_rank 00882 >> remote_rank_map_val.old_rank 00883 >> remote_rank_map_val.priority; 00884 remote_rank_map[new_rank] = remote_rank_map_val; 00885 } 00886 00887 typename rank_map_t::iterator 00888 the_remote_rank_map = remote_rank_map.begin(), 00889 end_remote_rank_map = remote_rank_map.end(); 00890 for (;the_remote_rank_map != end_remote_rank_map;++ the_remote_rank_map) { 00891 const int& remote_new_rank = the_remote_rank_map->first; 00892 const int& remote_old_rank = the_remote_rank_map->second.old_rank; 00893 const int& remote_priority = the_remote_rank_map->second.priority; 00894 00895 typename rank_map_t::iterator 00896 the_pair = p_map->find(remote_new_rank); 00897 int errcode = 1; 00898 if (the_pair == p_map->end()) { 00899 errcode *= 2; 00900 } else { 00901 if (the_pair->second.old_rank != remote_old_rank) { 00902 errcode *= 3; 00903 } 00904 if (the_pair->second.priority != remote_priority) { 00905 errcode *= 5; 00906 } 00907 } 00908 if (errcode != 1) { 00909 std::cerr << " [errcode:" << errcode << "]"; 00910 print_geo_rank_map(*geo); 00911 std::cerr << "-[" << remote_rank 00912 << ":" << GDIM; 00913 print_rank_map(&remote_rank_map); 00914 std::cerr << "] "; 00915 MPI_Abort(_forest->communicator(), errcode); 00916 } 00917 } 00918 } 00919 00920 private: 00925 template <class GEO> void 00926 geometry_set_new_rank_down(const GEO& geo, int new_rank) const { 00927 rank_map_t * p_map = get_rank_map(geo); 00928 if (p_map == NULL) { 00929 p_map = new_rank_map(geo); 00930 } else if (p_map->find(new_rank) != p_map->end()) { 00931 return; 00932 } 00933 00934 (*p_map)[new_rank].old_rank = _forest->rank(); 00935 00937 for (u_int i = 0;i < geo.n_vertex;++ i) { 00938 geometry_set_new_rank_down(*geo.vertex[i], new_rank); 00939 } 00940 for (u_int i = 0;i < geo.n_boundary;++ i) { 00941 geometry_set_new_rank_down(*geo.boundary[i], new_rank); 00942 } 00943 if (geo.isRefined()) { 00944 for (u_int i = 0;i < geo.n_child;++ i) { 00945 geometry_set_new_rank_down(*geo.child[i], new_rank); 00946 } 00947 } 00948 } 00949 00953 void set_new_rank_down() { 00954 typename forest_t::RootIterator 00955 the_ele = _forest->beginRootElement(); 00956 u_int n_part = _new_rank.size(); 00957 for (u_int i = 0;i < n_part;++ i) { 00958 for (u_int j = _cut_point[i];j < _cut_point[i + 1];++ j, ++ the_ele) { 00959 geometry_set_new_rank_down(*the_ele, _new_rank[i]); 00960 } 00961 } 00962 assert (the_ele == _forest->endRootElement()); 00963 } 00964 00971 template <class GEO> void 00972 geometry_set_new_rank_up(const GEO& geo, int new_rank) const { 00973 rank_map_t * p_map = get_rank_map(geo); 00974 if (p_map == NULL) p_map = new_rank_map(geo); 00975 00976 (*p_map)[new_rank].old_rank = _forest->rank(); 00977 00979 for (u_int i = 0;i < geo.n_vertex;++ i) { 00980 geometry_set_new_rank_up(*geo.vertex[i], new_rank); 00981 } 00982 for (u_int i = 0;i < geo.n_boundary;++ i) { 00983 geometry_set_new_rank_up(*geo.boundary[i], new_rank); 00984 } 00985 GEO* const& p_parent = geo.parent; 00986 if (p_parent != NULL) { 00987 if (! is_on_this_new_rank(*p_parent, new_rank)) { 00988 geometry_set_new_rank_up(*p_parent, new_rank); 00989 } 00990 00991 for (u_int i = 0;i < p_parent->n_child;++ i) { 00992 GEO* const& p_sibling = p_parent->child[i]; 00993 if (p_sibling == &geo) continue; 00994 if (! is_on_this_new_rank(*p_sibling, new_rank)) { 00995 geometry_set_new_rank_up(*p_sibling, new_rank); 00996 } 00997 } 00998 } 00999 } 01000 01004 void set_new_rank_up() { 01005 typename forest_t::RootIterator 01006 the_ele = _forest->beginRootElement(); 01007 u_int n_part = _new_rank.size(); 01008 for (u_int i = 0;i < n_part;++ i) { 01009 for (u_int j = _cut_point[i];j < _cut_point[i + 1];++ j, ++ the_ele) { 01010 geometry_set_new_rank_up(*the_ele, _new_rank[i]); 01011 } 01012 } 01013 assert (the_ele == _forest->endRootElement()); 01014 } 01016 01018 private: 01034 void global_index() { 01038 new_property_id(_pid_global_idx); 01039 unsigned long idx = 0; 01040 typename forest_t::RootIterator 01041 the_ele = _forest->beginRootElement(), 01042 end_ele = _forest->endRootElement(); 01043 for (;the_ele != end_ele;++ the_ele) { 01044 HGeometry<dim,dow> * p_geo = &(*the_ele); 01045 do { 01046 if (p_geo->parent != NULL) { 01047 p_geo = p_geo->parent; 01048 } else break; 01049 } while (true); 01050 geometry_global_index(*p_geo, idx); 01051 } 01052 free_property_id(_pid_global_idx); 01053 01057 MPI_Comm comm = _forest->communicator(); 01058 unsigned long global_idx = idx; 01059 MPI_Scan(&idx, &global_idx, 1, MPI_UNSIGNED_LONG, MPI_SUM, comm); 01060 idx = global_idx - idx; 01061 01065 new_property_id(_pid_global_idx); 01066 the_ele = _forest->beginRootElement(); 01067 for (;the_ele != end_ele;++ the_ele) { 01068 HGeometry<dim,dow> * p_geo = &(*the_ele); 01069 do { 01070 if (p_geo->parent != NULL) { 01071 p_geo = p_geo->parent; 01072 } else break; 01073 } while (true); 01074 geometry_global_index(*p_geo, idx); 01075 } 01076 01077 if (_forest->n_rank() <= 1) return; 01082 Shared_type_filter::only<0> type_filter; 01083 #define SYNC_DATA(D) \ 01084 if (dim >= D) { \ 01085 sync_data(comm, _forest->template get_shared_list<D>(), *this, \ 01086 &this_t::template pack_global_index<D>, \ 01087 &this_t::template unpack_global_index<D>, \ 01088 type_filter); \ 01089 } 01090 SYNC_DATA(0); 01091 SYNC_DATA(1); 01092 SYNC_DATA(2); 01093 SYNC_DATA(3); 01094 #undef SYNC_DATA 01095 } 01096 01097 public: 01098 template <int GDIM> void 01099 pack_global_index(HGeometry<GDIM,dow> * geo, 01100 int remote_rank, 01101 AFEPack::ostream<>& os) { 01102 unsigned long * p_idx = get_global_idx(*geo); 01103 assert (p_idx != NULL); 01104 01105 os << *p_idx; 01106 } 01107 01113 template <int GDIM> void 01114 unpack_global_index(HGeometry<GDIM,dow> * geo, 01115 int remote_rank, 01116 AFEPack::istream<>& is) { 01117 assert ((_forest->get_shared_info(*geo) != NULL)); 01118 01119 unsigned long * p_idx = get_global_idx(*geo); 01120 assert (p_idx != NULL); 01121 01122 unsigned long remote_idx; 01123 is >> remote_idx; 01124 if (*p_idx > remote_idx) { 01125 *p_idx = remote_idx; 01126 } 01127 } 01128 01129 private: 01133 template <class GEO> void 01134 geometry_global_index(GEO& geo, 01135 unsigned long& idx) const { 01136 unsigned long * p_idx = get_global_idx(geo); 01137 if (p_idx != NULL) return; 01138 01139 rank_map_t * p_map = get_rank_map(geo); 01140 if (_forest->get_shared_info(geo) != NULL || 01141 (p_map != NULL && p_map->size() > 1)) { 01142 p_idx = new_global_idx(geo); 01143 *p_idx = idx ++; 01144 } 01145 01147 for (u_int i = 0;i < geo.n_vertex;++ i) { 01148 geometry_global_index(*geo.vertex[i], idx); 01149 } 01150 for (u_int i = 0;i < geo.n_boundary;++ i) { 01151 geometry_global_index(*geo.boundary[i], idx); 01152 } 01153 if (geo.isRefined()) { 01154 for (u_int i = 0;i < geo.n_child;++ i) { 01155 geometry_global_index(*geo.child[i], idx); 01156 } 01157 } 01158 } 01160 01162 public: 01163 void write_config_file(const std::string& dirname) { 01164 char filename[1024]; 01165 sprintf(filename, "%s/.config", dirname.c_str()); 01166 std::ofstream os(filename); 01167 os << _forest->n_rank() << "\t# number of old rank" << std::endl; 01168 os.close(); 01169 Migration::save_config(dirname); 01170 } 01171 01172 private: 01173 void export_birdview(birdview_t& ir_mesh, 01174 Migration::data_id_t data_id, 01175 u_int idx) { 01176 typename birdview_t::ActiveIterator 01177 the_ele = ir_mesh.beginActiveElement(), 01178 end_ele = ir_mesh.endActiveElement(); 01179 for (;the_ele != end_ele;++ the_ele) { 01180 HGeometry<dim,dow> * p_geo = the_ele->h_element; 01181 AFEPack::ostream<> os(p_geo->buffer[data_id]); 01182 os << idx; 01183 } 01184 } 01185 01186 std::vector<u_int> _n_orphans; 01187 std::vector<std::list<void *> > _orphans_ptr; 01188 01189 template <class GEO> void 01190 geometry_collect_orphans(GEO& geo, 01191 const property_id_t<>& pid, 01192 int new_rank) { 01193 if (geo.get_property(pid) != NULL) return; 01194 geo.new_property(pid); 01195 01196 if ((geo.get_property(_pid_is_saved) == NULL) && 01197 is_save_on_this_rank(geo, new_rank)) { 01198 _n_orphans[GEO::dim] += 1; 01199 _orphans_ptr[GEO::dim].push_back((void *)(&geo)); 01200 } 01201 01202 for (u_int i = 0;i < geo.n_vertex;++ i) { 01203 geometry_collect_orphans(*geo.vertex[i], pid, new_rank); 01204 } 01205 for (u_int i = 0;i < geo.n_boundary;++ i) { 01206 geometry_collect_orphans(*geo.boundary[i], pid, new_rank); 01207 } 01208 if (geo.isRefined()) { 01209 for (u_int i = 0;i < geo.n_child;++ i) { 01210 geometry_collect_orphans(*geo.child[i], pid, new_rank); 01211 } 01212 } 01213 } 01214 01215 void collect_orphans(int new_rank) { 01216 _n_orphans.clear(); 01217 _n_orphans.resize(dim + 1, 0); 01218 01219 _orphans_ptr.clear(); 01220 _orphans_ptr.resize(dim + 1); 01221 01222 property_id_t<> pid_is_collected; 01223 new_property_id(pid_is_collected); 01224 01225 typename forest_t::RootIterator 01226 the_ele = _forest->beginRootElement(), 01227 end_ele = _forest->endRootElement(); 01228 for (;the_ele != end_ele;++ the_ele) { 01229 HGeometry<dim,dow> * p_geo = &(*the_ele); 01230 do { 01231 if (p_geo->parent != NULL) { 01232 p_geo = p_geo->parent; 01233 } else break; 01234 } while (true); 01235 geometry_collect_orphans(*p_geo, pid_is_collected, new_rank); 01236 } 01237 } 01238 01239 template <int D> void 01240 save_dim_orphans(boost::archive::HGeometry_oarchive<this_t>& oa) { 01241 typedef HGeometry<D,dow> GEO; 01242 u_int n_orphans = _n_orphans[D]; 01243 #if 1 01244 std::cerr << D << "D(" << n_orphans << ") "; 01245 #endif 01246 assert (n_orphans == _orphans_ptr[D].size()); 01247 oa & boost::serialization::make_binary_object(&n_orphans, sizeof(u_int)); 01248 typename std::list<void *>::iterator 01249 the_ptr = _orphans_ptr[D].begin(), 01250 end_ptr = _orphans_ptr[D].end(); 01251 for (;the_ptr != end_ptr;++ the_ptr) { 01252 GEO * p_geo = (GEO *)(*the_ptr); 01253 oa & p_geo; 01254 } 01255 } 01256 01257 void save_orphans(boost::archive::HGeometry_oarchive<this_t>& oa) { 01258 #if 1 01259 std::cerr << "saving orphans: "; 01260 #endif 01261 collect_orphans(oa.new_rank()); 01262 if (dim >= 0) this->template save_dim_orphans<0>(oa); 01263 if (dim >= 1) this->template save_dim_orphans<1>(oa); 01264 if (dim >= 2) this->template save_dim_orphans<2>(oa); 01265 if (dim >= 3) this->template save_dim_orphans<3>(oa); 01266 #if 1 01267 std::cerr << std::endl; 01268 #endif 01269 } 01270 01271 public: 01272 void save_data(const std::string& dirname, 01273 birdview_set_t& bvs) { 01274 set_new_rank(); 01275 global_index(); 01276 01277 int dummy; 01278 char command[1024], *filename = command; 01279 01280 Migration::data_id_t id = Migration::name_to_id("__internal_birdviewflag__"); 01281 if (! Migration::is_valid(id)) { 01282 id = Migration::register_data_name("__internal_birdviewflag__"); 01283 } 01284 typename birdview_set_t::iterator 01285 the_bv = bvs.begin(), end_bv = bvs.end(); 01286 for (u_int idx = 0;the_bv != end_bv;++ the_bv) { 01287 export_birdview(*the_bv, id, idx ++); 01288 } 01289 Migration::sync_data_buffer(*_forest); 01290 01291 typedef boost::archive::HGeometry_oarchive<this_t> archive_t; 01292 typename forest_t::RootIterator 01293 the_ele = _forest->beginRootElement(); 01294 u_int n_part = _new_rank.size(); 01295 for (u_int i = 0;i < n_part;++ i) { 01296 sprintf(command, "mkdir -p %s && mkdir -p %s/%d", 01297 dirname.c_str(), dirname.c_str(), _new_rank[i]); 01298 dummy = system(command); 01299 sprintf(filename, "%s/%d/%d.dat", dirname.c_str(), 01300 _new_rank[i], _forest->rank()); 01301 std::ofstream os(filename, std::ios::binary); 01302 01303 { 01304 archive_t oa(*this, _new_rank[i], os); 01305 new_property_id(_pid_is_saved); 01306 01307 u_int n_ele = _cut_point[i + 1] - _cut_point[i]; 01308 oa & boost::serialization::make_binary_object(&n_ele, sizeof(u_int)); 01309 std::cerr << "saving data in " << filename 01310 << ", # macro element = " << n_ele << std::endl; 01311 01312 for (u_int j = _cut_point[i];j < _cut_point[i + 1];++ j) { 01313 HGeometry<dim,dow> * p_ele = *(the_ele ++); 01314 oa & p_ele; 01315 } 01316 01317 save_orphans(oa); 01318 free_property_id(_pid_is_saved); 01319 } 01320 01321 os.close(); 01322 } 01323 01324 if (_forest->rank() == 0) write_config_file(dirname); 01325 _forest->rootElement().swap(_nre); 01326 _nre.clear(); 01327 MPI_Barrier(_forest->communicator()); 01328 } 01330 01332 public: 01333 int load_config_file(const std::string& dirname) { 01334 Migration::load_config(dirname); 01335 01336 char filename[1024]; 01337 sprintf(filename, "%s/.config", dirname.c_str()); 01338 filtering_istream is; 01339 Migration::ensured_open_filtering_stream(filename, is); 01340 01341 int old_rank; 01342 is >> old_rank; 01343 return old_rank; 01344 } 01345 01346 private: 01347 void reconstruct_birdview(HElement<dim,dow>& ele, 01348 Migration::data_id_t data_id, 01349 u_int idx) { 01350 typedef Migration::buffer_t buffer_t; 01351 buffer_t& hbuf = ele.h_element->buffer; 01352 typename buffer_t::iterator it = hbuf.find(data_id); 01353 bool is_active = true; 01354 ele.value = 0; 01355 if (it == hbuf.end()) { 01356 is_active = false; 01357 } else { 01358 BinaryBuffer<>& buf = it->second; 01359 int n = buf.size()/sizeof(u_int), i; 01360 u_int idx1; 01361 AFEPack::istream<> is(buf); 01362 for (i = 0;i < n;++ i) { 01363 is >> idx1; 01364 if (idx1 == idx) break; 01365 } 01366 if (i == n) is_active = false; 01367 } 01368 if (! is_active) { 01369 assert (ele.h_element->isRefined()); 01370 ele.value = 1; 01371 ele.refine(); 01372 for (int i = 0;i < ele.n_child;++ i) { 01373 reconstruct_birdview(*ele.child[i], data_id, idx); 01374 } 01375 } 01376 } 01377 01378 void reconstruct_birdview(birdview_t& ir_mesh, 01379 Migration::data_id_t data_id, 01380 u_int idx) { 01381 typename birdview_t::RootIterator 01382 the_ele = ir_mesh.beginRootElement(), 01383 end_ele = ir_mesh.endRootElement(); 01384 for (;the_ele != end_ele;++ the_ele) { 01385 reconstruct_birdview(*the_ele, data_id, idx); 01386 } 01387 } 01388 01393 template <class GEO> 01394 void set_dummy_flag(GEO& geo) const { 01395 if ((_forest->get_shared_info(geo) != NULL) && 01396 (! _forest->is_dummy(geo))) { 01397 _forest->dummy(geo); 01398 } 01399 01400 for (u_int i = 0;i < geo.n_boundary;++ i) { 01401 this->set_dummy_flag(*geo.boundary[i]); 01402 } 01403 if (geo.isRefined()) { 01404 for (u_int i = 0;i < geo.n_child;++ i) { 01405 this->set_dummy_flag(*geo.child[i]); 01406 } 01407 } 01408 } 01409 01413 template <class GEO> 01414 void clear_dummy_flag(GEO& geo) const { 01415 if (_forest->is_dummy(geo)) { 01416 _forest->undummy(geo); 01417 } 01418 01419 for (u_int i = 0;i < geo.n_boundary;++ i) { 01420 this->clear_dummy_flag(*geo.boundary[i]); 01421 } 01422 if (geo.isRefined()) { 01423 for (u_int i = 0;i < geo.n_child;++ i) { 01424 this->clear_dummy_flag(*geo.child[i]); 01425 } 01426 } 01427 } 01428 01439 void set_dummy_flag() { 01440 typename forest_t::RootIterator 01441 the_ele = _forest->beginRootElement(), 01442 end_ele = _forest->endRootElement(); 01443 for (;the_ele != end_ele;++ the_ele) { 01444 const HGeometry<dim,dow> * p_geo = &(*the_ele); 01445 do { 01446 if (p_geo->parent != NULL) { 01447 p_geo = p_geo->parent; 01448 } else break; 01449 } while (true); 01450 set_dummy_flag(*p_geo); 01451 } 01452 01453 the_ele = _forest->beginRootElement(); 01454 for (;the_ele != end_ele;++ the_ele) { 01455 clear_dummy_flag(*the_ele); 01456 } 01457 } 01458 01459 private: 01460 template <int D> 01461 void load_dim_orphans(boost::archive::HGeometry_iarchive<this_t>& ia) const { 01462 typedef HGeometry<D,dow> GEO; 01463 u_int n_orphans; 01464 ia & boost::serialization::make_binary_object(&n_orphans, sizeof(u_int)); 01465 #if 1 01466 std::cerr << D << "D(" << n_orphans << "), "; 01467 #endif 01468 for (u_int i = 0;i < n_orphans;++ i) { 01469 GEO * p_geo; 01470 ia & p_geo; 01471 } 01472 } 01473 01474 void load_orphans(boost::archive::HGeometry_iarchive<this_t>& ia) const { 01475 #if 1 01476 std::cerr << "loading orphans: "; 01477 #endif 01478 if (dim >= 0) this->template load_dim_orphans<0>(ia); 01479 if (dim >= 1) this->template load_dim_orphans<1>(ia); 01480 if (dim >= 2) this->template load_dim_orphans<2>(ia); 01481 if (dim >= 3) this->template load_dim_orphans<3>(ia); 01482 #if 1 01483 std::cerr << std::endl; 01484 #endif 01485 } 01486 01488 template <class GEO> void 01489 geometry_check_integrity(const GEO& geo) const { 01493 if (geo.parent != NULL) { 01494 u_int i = 0; 01495 for (;i < geo.parent->n_child;++ i) { 01496 if (&geo == geo.parent->child[i]) break; 01497 } 01498 assert ((i < geo.parent->n_child)); 01499 } 01500 01504 for (u_int i = 0;i < geo.n_boundary;++ i) { 01505 assert ((geo.boundary[i] != NULL)); 01506 geometry_check_integrity(*geo.boundary[i]); 01507 } 01508 if (geo.isRefined()) { 01513 for (u_int i = 0;i < geo.n_child;++ i) { 01514 assert ((geo.child[i] != NULL)); 01515 assert ((&geo == geo.child[i]->parent)); 01516 geometry_check_integrity(*geo.child[i]); 01517 } 01518 } else { 01522 for (u_int i = 0;i < geo.n_child;++ i) { 01523 assert (geo.child[i] == NULL); 01524 } 01525 } 01526 } 01527 01528 public: 01529 void check_integrity() const { 01530 typename forest_t::RootIterator 01531 the_ele = _forest->beginRootElement(), 01532 end_ele = _forest->endRootElement(); 01533 for (;the_ele != end_ele;++ the_ele) { 01534 HGeometry<dim,dow> * p_geo = &(*the_ele); 01535 do { 01536 if (p_geo->parent != NULL) { 01537 p_geo = p_geo->parent; 01538 } else break; 01539 } while (true); 01540 geometry_check_integrity(*p_geo); 01541 } 01542 } 01543 01544 public: 01545 void load_data(const std::string& dirname, 01546 birdview_set_t& bvs, 01547 bool is_has_orphans = false) { 01548 clear(); 01549 01550 if (rank() == 6) 01551 std::cerr << "is_has_orphans=" << is_has_orphans << "\n"; 01552 01553 int n_old_rank = load_config_file(dirname); 01554 01556 MPI_Barrier(_forest->communicator()); 01557 01559 char filename[1024]; 01560 01561 typedef boost::archive::HGeometry_iarchive<this_t> archive_t; 01562 01563 int rank = 0; 01564 do { 01565 std::ifstream is; 01566 do { 01567 sprintf(filename, "%s/%d/%d.dat", dirname.c_str(), 01568 _forest->rank(), rank ++); 01569 if (rank > n_old_rank) { 01570 if (! _global_pointer_to_merge_later.empty()) { 01571 std::cerr << "Rank " << _forest->rank() << ": "; 01572 std::map<unsigned long, std::list<void *> >::iterator 01573 the_entry = _global_pointer_to_merge_later.begin(), 01574 end_entry = _global_pointer_to_merge_later.end(); 01575 for (;the_entry != end_entry;++ the_entry) { 01576 std::cerr << "(" << the_entry->first 01577 << "->" << the_entry->second.size() 01578 << ")"; 01579 } 01580 std::cerr << std::endl; 01581 assert (false); 01582 } 01583 01584 this->share_global_pointer(); 01585 this->coarse_root_element(); 01586 this->set_dummy_flag(); 01587 _forest->set_is_used_up(); 01588 01589 //_forest->verify_shared_object(); 01590 //this->check_integrity(); /// 调试时候使用 01591 01593 Migration::data_id_t id = Migration::name_to_id("__internal_birdviewflag__"); 01594 typename birdview_set_t::iterator 01595 the_bv = bvs.begin(), end_bv = bvs.end(); 01596 for (u_int idx = 0;the_bv != end_bv;++ the_bv) { 01597 the_bv->reinit(*_forest); 01598 reconstruct_birdview(*the_bv, id, idx ++); 01599 } 01600 return; 01601 } 01602 is.open(filename, std::ios::binary); 01603 if (is.good()) break; 01604 } while (1); 01605 01606 { 01607 u_int n_ele; 01608 std::cerr << "loading data from " << filename; 01609 archive_t ia(*this, _forest->rank(), is); 01610 ia & boost::serialization::make_binary_object(&n_ele, sizeof(u_int)); 01611 std::cerr << ", # macro element = " << n_ele << std::endl; 01612 01613 for (u_int i = 0;i < n_ele;++ i) { 01614 HGeometry<dim,dow> * p_ele; 01615 ia & p_ele; 01616 _forest->rootElement().push_back(p_ele); 01617 } 01618 01619 if (is_has_orphans) load_orphans(ia); 01620 } 01621 01622 is.close(); 01623 } while (1); 01624 01625 } 01626 01627 }; 01628 01629 template <class FOREST> 01630 void load_forest(const std::string& dirname, 01631 FOREST& forest, 01632 bool is_has_orphans = false) { 01633 HLoadBalance<FOREST> hlb(forest); 01634 BirdViewSet<FOREST> bvs; 01635 hlb.load_data(dirname, bvs, is_has_orphans); 01636 } 01637 01638 template <class FOREST> 01639 void load_mesh(const std::string& dirname, 01640 FOREST& forest, 01641 BirdView<FOREST>& mesh, 01642 bool is_has_orphans = false) { 01643 HLoadBalance<FOREST> hlb(forest); 01644 BirdViewSet<FOREST> bvs; 01645 bvs.add(mesh); 01646 hlb.load_data(dirname, bvs, is_has_orphans); 01647 } 01648 01649 template <class FOREST> 01650 void load_mesh(const std::string& dirname, 01651 bool is_has_orphans, 01652 FOREST& forest, 01653 u_int n_mesh, ...) { 01654 HLoadBalance<FOREST> hlb(forest); 01655 typedef BirdView<FOREST> * mesh_ptr_t; 01656 01657 BirdViewSet<FOREST> bvs; 01658 va_list ap; 01659 va_start(ap, n_mesh); 01660 for (u_int i = 0;i < n_mesh;++ i) { 01661 mesh_ptr_t p_mesh = va_arg(ap, mesh_ptr_t); 01662 bvs.add(*p_mesh); 01663 } 01664 va_end(ap); 01665 01666 hlb.load_data(dirname, bvs, is_has_orphans); 01667 } 01668 01669 template <class FOREST> 01670 void load_mesh_set(const std::string& dirname, 01671 FOREST& forest, 01672 BirdViewSet<FOREST>& bvs, 01673 bool is_has_orphans = false) { 01674 HLoadBalance<FOREST> hlb(forest); 01675 hlb.load_data(dirname, bvs, is_has_orphans); 01676 } 01677 01685 void lb_collect_local_data_dir(MPI_Comm comm, 01686 const std::string& src_dir, 01687 const std::string& dst_dir); 01693 void lb_sync_local_data_dir(MPI_Comm comm, 01694 const std::string& src_dir, 01695 const std::string& dst_dir); 01696 } 01697 01698 #endif // __MPI_ULoadBalance_h__ 01699