diff --git a/src/driver/driver.cpp b/src/driver/driver.cpp index fe761ded..0935b380 100644 --- a/src/driver/driver.cpp +++ b/src/driver/driver.cpp @@ -518,6 +518,20 @@ Real Driver::UpdateWallClock() { void Driver::InitBoundaryValuesAndPrimitives(Mesh *pm) { // Note: with MPI, sends on ALL MBs must be complete before receives execute + // Initialize Z4c + z4c::Z4c *pz4c = pm->pmb_pack->pz4c; + if (pz4c != nullptr) { + (void) pz4c->RestrictU(this, 0); + (void) pz4c->InitRecv(this, -1); // stage < 0 suppresses InitFluxRecv + (void) pz4c->SendU(this, 0); + (void) pz4c->ClearSend(this, -1); + (void) pz4c->ClearRecv(this, -1); + (void) pz4c->RecvU(this, 0); + (void) pz4c->Z4cBoundaryRHS(this, 0); + (void) pz4c->ApplyPhysicalBCs(this, 0); + (void) pz4c->Prolongate(this, 0); + } + // Initialize HYDRO: ghost zones and primitive variables (everywhere) // includes communications for shearing box boundaries hydro::Hydro *phydro = pm->pmb_pack->phydro; @@ -563,7 +577,10 @@ void Driver::InitBoundaryValuesAndPrimitives(Mesh *pm) { if (pdyngr == nullptr) { (void) pmhd->ConToPrim(this, 0); } else { - pdyngr->ConToPrim(this, 0); + if (pz4c != nullptr) { + (void) pz4c->Z4cToADM_(this, 0); + } + (void) pdyngr->ConToPrim(this, 0); } } @@ -581,19 +598,5 @@ void Driver::InitBoundaryValuesAndPrimitives(Mesh *pm) { (void) prad->Prolongate(this, 0); } - // Initialize Z4c - z4c::Z4c *pz4c = pm->pmb_pack->pz4c; - if (pz4c != nullptr) { - (void) pz4c->RestrictU(this, 0); - (void) pz4c->InitRecv(this, -1); // stage < 0 suppresses InitFluxRecv - (void) pz4c->SendU(this, 0); - (void) pz4c->ClearSend(this, -1); - (void) pz4c->ClearRecv(this, -1); - (void) pz4c->RecvU(this, 0); - (void) pz4c->Z4cBoundaryRHS(this, 0); - (void) pz4c->ApplyPhysicalBCs(this, 0); - (void) pz4c->Prolongate(this, 0); - } - return; } diff --git a/src/mesh/mesh_refinement.cpp b/src/mesh/mesh_refinement.cpp index d9796cca..1a644398 100644 --- a/src/mesh/mesh_refinement.cpp +++ b/src/mesh/mesh_refinement.cpp @@ -631,7 +631,7 @@ void MeshRefinement::RedistAndRefineMeshBlocks(ParameterInput *pin, int nnew, in RefineFC(new_to_old, pmhd->b0, pmhd->coarse_b0); } if (pz4c != nullptr) { - RefineCC(new_to_old, pz4c->u0, pz4c->coarse_u0); + RefineCC(new_to_old, pz4c->u0, pz4c->coarse_u0, true); } } @@ -1027,11 +1027,10 @@ void MeshRefinement::CopyForRefinementFC(DvceFaceFld4D &b,DvceFaceFld4D &n2o, DvceArray5D &a, - DvceArray5D &ca) { + DvceArray5D &ca, bool is_z4c) { int nvar = a.extent_int(1); // TODO(@user): 2nd index from L of in array must be NVAR auto &new_nmb = new_nmb_eachrank[global_variable::my_rank]; MeshBlockPack* pmbp = pmy_mesh->pmb_pack; - bool not_z4c = (pmbp->pz4c == nullptr)? true : false; auto &indcs = pmy_mesh->mb_indcs; auto &cis = indcs.cis, &cie = indcs.cie; auto &cjs = indcs.cjs, &cje = indcs.cje; @@ -1074,7 +1073,7 @@ void MeshRefinement::RefineCC(DualArray1D &n2o, DvceArray5D &a, int fk = 2*k - cks; // correct when cks=ks // call inlined prolongation operator for CC variables - if (not_z4c) { + if (!is_z4c) { ProlongCC(m,v,k,j,i,fk,fj,fi,multi_d,three_d,ca,a); } else { switch (indcs.ng) { @@ -1184,7 +1183,6 @@ void MeshRefinement::RestrictCC(DvceArray5D &u, DvceArray5D &cu, int nvar = u.extent_int(1); // TODO(@user): 2nd index from L of in array must be NVAR MeshBlockPack* pmbp = pmy_mesh->pmb_pack; - bool not_z4c = (pmbp->pz4c == nullptr)? true : false; auto &indcs = pmy_mesh->mb_indcs; auto &cis = indcs.cis, &cie = indcs.cie; auto &cjs = indcs.cjs, &cje = indcs.cje; @@ -1219,7 +1217,7 @@ void MeshRefinement::RestrictCC(DvceArray5D &u, DvceArray5D &cu, int finei = 2*i - cis; // correct when cis=is int finej = 2*j - cjs; // correct when cjs=js int finek = 2*k - cks; // correct when cks=ks - if (not_z4c) { + if (!is_z4c) { cu(m,n,k,j,i) = 0.125*(u(m,n,finek ,finej ,finei) + u(m,n,finek ,finej ,finei+1) + u(m,n,finek ,finej+1,finei) + u(m,n,finek ,finej+1,finei+1) diff --git a/src/mesh/mesh_refinement.hpp b/src/mesh/mesh_refinement.hpp index 3618f674..e8070ece 100644 --- a/src/mesh/mesh_refinement.hpp +++ b/src/mesh/mesh_refinement.hpp @@ -107,7 +107,8 @@ class MeshRefinement { void CopyForRefinementCC(DvceArray5D &a, DvceArray5D &ca); void CopyForRefinementFC(DvceFaceFld4D &b, DvceFaceFld4D &cb); - void RefineCC(DualArray1D &n2o, DvceArray5D &a, DvceArray5D &ca); + void RefineCC(DualArray1D &n2o, DvceArray5D &a, DvceArray5D &ca, + bool is_z4c=false); void RefineFC(DualArray1D &n2o, DvceFaceFld4D &b, DvceFaceFld4D &cb); void RestrictCC(DvceArray5D &a, DvceArray5D &ca, bool is_z4c=false); diff --git a/src/z4c/z4c_tasks.cpp b/src/z4c/z4c_tasks.cpp index c9af84d1..813dc08b 100644 --- a/src/z4c/z4c_tasks.cpp +++ b/src/z4c/z4c_tasks.cpp @@ -421,7 +421,7 @@ TaskStatus Z4c::RestrictWeyl(Driver *pdrive, int stage) { float time_32 = static_cast(pmy_pack->pmesh->time); if ((last_output_time==time_32) && (stage == pdrive->nexp_stages)) { if (pmy_pack->pmesh->multilevel) { - pmy_pack->pmesh->pmr->RestrictCC(u_weyl, coarse_u_weyl); + pmy_pack->pmesh->pmr->RestrictCC(u_weyl, coarse_u_weyl, true); } } return TaskStatus::complete;