| 467 | } |
| 468 | |
| 469 | void |
| 470 | Gravity::gravity_sync (int crse_level, int fine_level, const Vector<MultiFab*>& drho, const Vector<MultiFab*>& dphi) |
| 471 | { |
| 472 | BL_PROFILE("Gravity::gravity_sync()"); |
| 473 | |
| 474 | // There is no need to do a synchronization if |
| 475 | // we didn't solve on the fine levels. |
| 476 | |
| 477 | if (fine_level > gravity::max_solve_level) { |
| 478 | return; |
| 479 | } else { |
| 480 | fine_level = amrex::min(fine_level, gravity::max_solve_level); |
| 481 | } |
| 482 | |
| 483 | BL_ASSERT(parent->finestLevel()>crse_level); |
| 484 | if (gravity::verbose > 1 && ParallelDescriptor::IOProcessor()) { |
| 485 | std::cout << " ... gravity_sync at crse_level " << crse_level << '\n'; |
| 486 | std::cout << " ... up to finest_level " << fine_level << '\n'; |
| 487 | } |
| 488 | |
| 489 | const Geometry& crse_geom = parent->Geom(crse_level); |
| 490 | const Box& crse_domain = crse_geom.Domain(); |
| 491 | |
| 492 | int nlevs = fine_level - crse_level + 1; |
| 493 | |
| 494 | // Construct delta(phi) and delta(grad_phi). delta(phi) |
| 495 | // needs a ghost zone for holding the boundary condition |
| 496 | // in the same way that phi does. |
| 497 | |
| 498 | Vector<std::unique_ptr<MultiFab> > delta_phi(nlevs); |
| 499 | |
| 500 | for (int lev = crse_level; lev <= fine_level; ++lev) { |
| 501 | delta_phi[lev - crse_level] = std::make_unique<MultiFab>(grids[lev], dmap[lev], 1, 1); |
| 502 | delta_phi[lev - crse_level]->setVal(0.0); |
| 503 | } |
| 504 | |
| 505 | Vector< Vector<std::unique_ptr<MultiFab> > > ec_gdPhi(nlevs); |
| 506 | |
| 507 | for (int lev = crse_level; lev <= fine_level; ++lev) { |
| 508 | ec_gdPhi[lev - crse_level].resize(AMREX_SPACEDIM); |
| 509 | |
| 510 | const DistributionMapping& dm = LevelData[lev]->DistributionMap(); |
| 511 | for (int n = 0; n < AMREX_SPACEDIM; ++n) { |
| 512 | ec_gdPhi[lev - crse_level][n] = std::make_unique<MultiFab>(LevelData[lev]->getEdgeBoxArray(n), dm, 1, 0); |
| 513 | ec_gdPhi[lev - crse_level][n]->setVal(0.0); |
| 514 | } |
| 515 | } |
| 516 | |
| 517 | // Construct a container for the right-hand-side (4 * pi * G * drho + dphi). |
| 518 | // dphi appears in the construction of the boundary conditions because it |
| 519 | // indirectly represents a change in mass on the domain (the mass motion that |
| 520 | // occurs on the fine grid, whose gravitational effects are now indirectly |
| 521 | // being propagated to the coarse grid). |
| 522 | |
| 523 | // We will temporarily leave the RHS divided by (4 * pi * G) because that |
| 524 | // is the form expected by the boundary condition routine. |
| 525 | |
| 526 | Vector<std::unique_ptr<MultiFab> > g_rhs(nlevs); |