From 3fb5bf288642e3843e325d6c0db0f8aee9e28ea3 Mon Sep 17 00:00:00 2001 From: knarf Date: Wed, 25 Aug 2010 18:29:48 -0400 Subject: implement unusedpoints_mask --- Carpet/CarpetLib/src/dh.cc | 526 ++++++++++++++++++++++++++++++++++++--------- 1 file changed, 420 insertions(+), 106 deletions(-) (limited to 'Carpet/CarpetLib/src/dh.cc') diff --git a/Carpet/CarpetLib/src/dh.cc b/Carpet/CarpetLib/src/dh.cc index 9679b3fd7..31c60876f 100644 --- a/Carpet/CarpetLib/src/dh.cc +++ b/Carpet/CarpetLib/src/dh.cc @@ -28,6 +28,27 @@ list dh::alldh; +vect,dim> +dh::fast_dboxes:: +fast_ref_refl_sendrecv; + +void +dh::fast_dboxes:: +init_fast_ref_refl_sendrecv () +{ + static bool initialised = false; + if (initialised) return; + initialised = true; + fast_ref_refl_sendrecv[0][0] = & fast_dboxes::fast_ref_refl_sendrecv_0_0; + fast_ref_refl_sendrecv[0][1] = & fast_dboxes::fast_ref_refl_sendrecv_0_1; + fast_ref_refl_sendrecv[1][0] = & fast_dboxes::fast_ref_refl_sendrecv_1_0; + fast_ref_refl_sendrecv[1][1] = & fast_dboxes::fast_ref_refl_sendrecv_1_1; + fast_ref_refl_sendrecv[2][0] = & fast_dboxes::fast_ref_refl_sendrecv_2_0; + fast_ref_refl_sendrecv[2][1] = & fast_dboxes::fast_ref_refl_sendrecv_2_1; +} + + + // Constructors dh:: dh (gh & h_, @@ -37,6 +58,7 @@ dh (gh & h_, ghost_widths(ghost_widths_), buffer_widths(buffer_widths_), prolongation_orders_space(prolongation_orders_space_) { + fast_dboxes::init_fast_ref_refl_sendrecv(); size_t const maxreflevels = h.reffacts.size(); assert (ghost_widths.size() >= maxreflevels); assert (buffer_widths.size() >= maxreflevels); @@ -923,8 +945,8 @@ regrid (bool const do_init) full_dboxes const & box = full_level.AT(c); full_dboxes const & obox0 = full_boxes.AT(ml).AT(orl).AT(0); - // Refinement restriction may fill all active points, and - // must use all active points + // Refinement restriction may fill all coarse active + // points, and must use all fine active points ibset needrecv (box.active.contracted_for (obox0.interior)); @@ -932,7 +954,7 @@ regrid (bool const do_init) full_dboxes & obox = full_boxes.AT(ml).AT(orl).AT(cc); ibset const contracted_active - (box.active.contracted_for (obox0.interior)); + (box.active.contracted_for (obox.interior)); ibset const ovlp = obox.active & contracted_active; for (ibset::const_iterator @@ -975,6 +997,9 @@ regrid (bool const do_init) static Carpet::Timer timer_mask ("CarpetLib::dh::regrid::mask"); timer_mask.start(); + // Declare this here to save it for 'unused-mask' + ibset all_refined; + if (rl > 0) { int const orl = rl - 1; full_cboxes const& full_olevel = full_boxes.AT(ml).AT(orl); @@ -985,101 +1010,150 @@ regrid (bool const do_init) ivect const reffact = h.reffacts.AT(rl) / h.reffacts.AT(orl); // This works only when the refinement factor is 2 - if (all (reffact == 2)) { - - ibbox const& base = domain_exterior; - ibbox const& obase = h.baseextent(ml,orl); - - // Calculate the union of all coarse regions - ibset const allointr (full_olevel, & full_dboxes::interior); - - // Project to current level - ivect const rf(reffact); - // TODO: Why expand by rf-1? This expansion shouldn't - // matter, since the coarse grid is larger than the fine - // grid anyway, except at the outer boundary - // ibset const parent (allointr.expand(rf-1,rf-1).contracted_for(base)); - ibset const parent (allointr.expanded_for(base)); - - // Subtract the active region - ibset const notrefined = parent - allactive; - - // Enlarge this set - vect enlarged; - for (int d=0; d enlarged; + for (int d=0; d all_boundaries; + for (int d=0; d all_boundaries; + // Set prolongation information for current level for (int d=0; d 0) { + int const orl = rl - 1; + full_cboxes const& full_olevel = full_boxes.AT(ml).AT(orl); + // Local boxes are not communicated or post-processed, and + // thus can be modified even for coarser levels + local_cboxes& local_olevel = local_boxes.AT(ml).AT(orl); + + // This works only when the refinement factor is 2 + ivect const reffact = h.reffacts.AT(rl) / h.reffacts.AT(orl); + if (all (reffact == 2)) { + // use the already computed 'all_refined' to get region from where + // no information will be used later (overwritten) + // First: get the region which will get restricted + ibset restricted_region = all_refined.contracted_for(h.baseextent(ml,orl)); + // This is too big - during MoL-substeps information within this + // region will be used to update points outside -> need to + // shrink it by some points + // The way we shrink it is to invert it, expand that, and invert + // again. To invert we define an enclosing box and subtract it from that. + i2vect to_shrink = buffer_widths[orl] + ghost_widths[orl]; + ibbox enclosing = restricted_region.container().expand(ivect(1)+to_shrink); + ibset unused_region = enclosing - (enclosing - restricted_region).expand(to_shrink); + // Now we have the interesting region in 'unused_region' and need to store + // the union of this and the local regions + for (int lc = 0; lc < h.local_components(orl); ++ lc) { + int const c = h.get_component(orl, lc); full_dboxes const& obox = full_olevel.AT(c); - local_dboxes & local_box = local_level.AT(lc); // Local boxes are not communicated or post-processed, and // thus can be modified even for coarser levels local_dboxes & local_obox = local_olevel.AT(lc); - - // Set restriction information for next coarser level - local_obox.restricted_region = - all_refined.contracted_for(obox.exterior) & obox.owned; - - // Set prolongation information for current level - for (int d=0; d 0) { int const orl = rl - 1; - full_cboxes const& full_olevel = full_boxes.AT(ml).AT(orl); + light_cboxes & light_olevel = light_boxes.AT(ml).AT(orl); local_cboxes & local_olevel = local_boxes.AT(ml).AT(orl); + full_cboxes const& full_olevel = full_boxes.AT(ml).AT(orl); + fast_dboxes & fast_olevel = fast_boxes.AT(ml).AT(orl); // This works only with cell centred grids if (h.refcent == cell_centered) { @@ -1099,6 +1175,9 @@ regrid (bool const do_init) ivect const izero = ivect (0); ivect const ione = ivect (1); + assert (all (h.reffacts.AT(rl) % h.reffacts.AT(orl) == 0)); + ivect const reffact = h.reffacts.AT(rl) / h.reffacts.AT(orl); + vect,dim> all_fine_boundary; for (int dir=0; dir 0 @@ -1452,7 +1706,25 @@ regrid (bool const do_init) timer_bcast_comm_ref_rest.start(); broadcast_schedule (fast_level_otherprocs, fast_olevel, & fast_dboxes::fast_ref_rest_sendrecv); - timer_bcast_comm_ref_rest.stop(); + timer_bcast_comm_ref_rest.stop(); + } + + if (rl > 0) { + int const orl = rl - 1; + fast_dboxes & fast_olevel = fast_boxes.AT(ml).AT(orl); + static Carpet::Timer timer_bcast_comm_ref_refl + ("CarpetLib::dh::regrid::bcast_comm::ref_refl"); + timer_bcast_comm_ref_refl.start(); + for (int dir = 0; dir < dim; ++ dir) { + for (int face = 0; face < 2; ++ face) { + srpvect fast_dboxes::* const fast_ref_refl_sendrecv = + fast_dboxes::fast_ref_refl_sendrecv[dir][face]; + + broadcast_schedule (fast_level_otherprocs, fast_olevel, + fast_ref_refl_sendrecv); + } + } + timer_bcast_comm_ref_refl.stop(); } // TODO: Maybe broadcast old2new schedule only if do_init is @@ -1499,25 +1771,7 @@ regrid (bool const do_init) cout << box; } // for c - for (int c = 0; c < h.components(rl); ++ c) { - light_dboxes const & box = light_boxes.AT(ml).AT(rl).AT(c); - cout << eol; - cout << "ml=" << ml << " rl=" << rl << " c=" << c << eol; - cout << box; - } // for c - - for (int lc = 0; lc < h.local_components(rl); ++ lc) { - int const c = h.get_component (rl, lc); - local_dboxes const & box = local_boxes.AT(ml).AT(rl).AT(lc); - cout << eol; - cout << "ml=" << ml << " rl=" << rl << " lc=" << lc << " c=" << c << eol; - cout << box; - } // for lc - - fast_dboxes const & fast_box = fast_boxes.AT(ml).AT(rl); - cout << eol; - cout << "ml=" << ml << " rl=" << rl << eol; - cout << fast_box; + // light, local, and fast boxes are output later } // if output_bboxes @@ -1563,6 +1817,46 @@ regrid (bool const do_init) // Output: if (output_bboxes or there_was_an_error) { + for (int ml = 0; ml < h.mglevels(); ++ ml) { + for (int rl = 0; rl < h.reflevels(); ++ rl) { + + cout << eol; + cout << "ml=" << ml << " rl=" << rl << eol; + cout << "baseextent=" << h.baseextent(ml,rl) << eol; + + for (int c = 0; c < h.components(rl); ++ c) { + cout << eol; + cout << "ml=" << ml << " rl=" << rl << " c=" << c << eol; + cout << "extent=" << h.extent(ml,rl,c) << eol; + cout << "outer_boundaries=" << h.outer_boundaries(rl,c) << eol; + cout << "processor=" << h.outer_boundaries(rl,c) << eol; + } // for c + + // full boxes have already been output (and deallocated) + + for (int c = 0; c < h.components(rl); ++ c) { + light_dboxes const & box = light_boxes.AT(ml).AT(rl).AT(c); + cout << eol; + cout << "ml=" << ml << " rl=" << rl << " c=" << c << eol; + cout << box; + } // for c + + for (int lc = 0; lc < h.local_components(rl); ++ lc) { + int const c = h.get_component (rl, lc); + local_dboxes const & box = local_boxes.AT(ml).AT(rl).AT(lc); + cout << eol; + cout << "ml=" << ml << " rl=" << rl << " lc=" << lc << " c=" << c << eol; + cout << box; + } // for lc + + fast_dboxes const & fast_box = fast_boxes.AT(ml).AT(rl); + cout << eol; + cout << "ml=" << ml << " rl=" << rl << eol; + cout << fast_box; + + } // for rl + } // for ml + cout << eol; cout << "memoryof(gh)=" << memoryof(h) << eol; cout << "memoryof(dh)=" << memoryof(*this) << eol; @@ -1823,6 +2117,12 @@ mpi_datatype (dh::fast_dboxes const &) ENTRY (dh::srpvect, fast_ref_bnd_prol_sendrecv), ENTRY (dh::srpvect, fast_old2new_sync_sendrecv), ENTRY (dh::srpvect, fast_old2new_ref_prol_sendrecv), + ENTRY (dh::srpvect, fast_ref_refl_sendrecv_0_0), + ENTRY (dh::srpvect, fast_ref_refl_sendrecv_0_1), + ENTRY (dh::srpvect, fast_ref_refl_sendrecv_1_0), + ENTRY (dh::srpvect, fast_ref_refl_sendrecv_1_1), + ENTRY (dh::srpvect, fast_ref_refl_sendrecv_2_0), + ENTRY (dh::srpvect, fast_ref_refl_sendrecv_2_1), {1, sizeof s, MPI_UB, "MPI_UB", "MPI_UB"} }; #undef ENTRY @@ -1960,6 +2260,13 @@ memory () memoryof (fast_ref_rest_sendrecv) + memoryof (fast_sync_sendrecv) + memoryof (fast_ref_bnd_prol_sendrecv) + + memoryof (fast_ref_refl_sendrecv_0_0) + + memoryof (fast_ref_refl_sendrecv_0_1) + + memoryof (fast_ref_refl_sendrecv_1_0) + + memoryof (fast_ref_refl_sendrecv_1_1) + + memoryof (fast_ref_refl_sendrecv_2_0) + + memoryof (fast_ref_refl_sendrecv_2_1) + + memoryof (do_init) + memoryof (fast_old2new_sync_sendrecv) + memoryof (fast_old2new_ref_prol_sendrecv); } @@ -2233,6 +2540,13 @@ output (ostream & os) << " fast_ref_rest_sendrecv: " << fast_ref_rest_sendrecv << eol << " fast_sync_sendrecv: " << fast_sync_sendrecv << eol << " fast_ref_bnd_prol_sendrecv: " << fast_ref_bnd_prol_sendrecv << eol + << " fast_ref_refl_sendrecv_0_0: " << fast_ref_refl_sendrecv_0_0 << eol + << " fast_ref_refl_sendrecv_0_1: " << fast_ref_refl_sendrecv_0_1 << eol + << " fast_ref_refl_sendrecv_1_0: " << fast_ref_refl_sendrecv_1_0 << eol + << " fast_ref_refl_sendrecv_1_1: " << fast_ref_refl_sendrecv_1_1 << eol + << " fast_ref_refl_sendrecv_2_0: " << fast_ref_refl_sendrecv_2_0 << eol + << " fast_ref_refl_sendrecv_2_1: " << fast_ref_refl_sendrecv_2_1 << eol + << " do_init: " << do_init << eol << " fast_old2new_sync_sendrecv: " << fast_old2new_sync_sendrecv << eol << " fast_old2new_ref_prol_sendrecv: " << fast_old2new_ref_prol_sendrecv << eol << "}" << eol; -- cgit v1.2.3