From 759e7bb4b8ce296d98b9d4aec07c36298900afd9 Mon Sep 17 00:00:00 2001 From: Ann Almgren Date: Wed, 16 Oct 2024 10:18:52 -0700 Subject: [PATCH 1/2] fix for real_bcs, remove some more fb calls --- .../ERF_BoundaryConditions_realbdy.cpp | 2 + Source/BoundaryConditions/ERF_FillPatch.cpp | 66 +++++++++---------- Source/BoundaryConditions/ERF_PhysBCFunct.H | 12 ++-- Source/BoundaryConditions/ERF_PhysBCFunct.cpp | 27 +++++--- 4 files changed, 59 insertions(+), 48 deletions(-) diff --git a/Source/BoundaryConditions/ERF_BoundaryConditions_realbdy.cpp b/Source/BoundaryConditions/ERF_BoundaryConditions_realbdy.cpp index 289349d90..5d6be8b7d 100644 --- a/Source/BoundaryConditions/ERF_BoundaryConditions_realbdy.cpp +++ b/Source/BoundaryConditions/ERF_BoundaryConditions_realbdy.cpp @@ -67,6 +67,8 @@ ERF::fill_from_realbdy (const Vector& mfs, { MultiFab& mf = *mfs[var_idx]; + mf.FillBoundary(geom[lev].periodicity()); + // // Note that "domain" is mapped onto the type of box the data is in // diff --git a/Source/BoundaryConditions/ERF_FillPatch.cpp b/Source/BoundaryConditions/ERF_FillPatch.cpp index 4761a43c3..d8e86b166 100644 --- a/Source/BoundaryConditions/ERF_FillPatch.cpp +++ b/Source/BoundaryConditions/ERF_FillPatch.cpp @@ -111,8 +111,7 @@ ERF::FillPatch (int lev, Real time, // Impose physical bc's on coarse data (note time and 0 are not used) // Note that we call FillBoundary inside the physbcs call // We should not need to call this on old data since that would have been filled before the timestep started - // (*physbcs_cons[lev-1])(vars_old[lev-1][Vars::cons],0,mf_c.nComp(),ngvect_cons,time,BCVars::cons_bc); - (*physbcs_cons[lev-1])(vars_new[lev-1][Vars::cons],0,mf_c.nComp(),ngvect_cons,time,BCVars::cons_bc); + (*physbcs_cons[lev-1])(vars_new[lev-1][Vars::cons],0,mf_c.nComp(),ngvect_cons,time,BCVars::cons_bc,true); // Call FillPatchTwoLevels which ASSUMES that all ghost cells have already been filled FillPatchTwoLevels(mf_c, ngvect_cons, IntVect(0,0,0), @@ -136,8 +135,7 @@ ERF::FillPatch (int lev, Real time, // Impose physical bc's on coarse data (note time and 0 are not used) // Note that we call FillBoundary inside the physbcs call // We should not need to call this on old data since that would have been filled before the timestep started - // (*physbcs_u[lev-1])(vars_old[lev-1][Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc); - (*physbcs_u[lev-1])(vars_new[lev-1][Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc); + (*physbcs_u[lev-1])(vars_new[lev-1][Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc,true); fmf = {&vars_old[lev ][Vars::xvel], &vars_new[lev ][Vars::xvel]}; cmf = {&vars_old[lev-1][Vars::xvel], &vars_new[lev-1][Vars::xvel]}; @@ -156,8 +154,7 @@ ERF::FillPatch (int lev, Real time, // Impose physical bc's on coarse data (note time and 0 are not used) // Note that we call FillBoundary inside the physbcs call // We should not need to call this on old data since that would have been filled before the timestep started - // (*physbcs_v[lev-1])(vars_old[lev-1][Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc); - (*physbcs_v[lev-1])(vars_new[lev-1][Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc); + (*physbcs_v[lev-1])(vars_new[lev-1][Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc,true); fmf = {&vars_old[lev ][Vars::yvel], &vars_new[lev ][Vars::yvel]}; cmf = {&vars_old[lev-1][Vars::yvel], &vars_new[lev-1][Vars::yvel]}; @@ -176,14 +173,10 @@ ERF::FillPatch (int lev, Real time, // Impose physical bc's on coarse data (note time and 0 are not used) // Note that we call FillBoundary inside the physbcs call // We should not need to call this on old data since that would have been filled before the timestep started - // (*physbcs_w[lev-1])(vars_old[lev-1][Vars::zvel], - // vars_old[lev-1][Vars::xvel], - // vars_old[lev-1][Vars::yvel], - // ngvect_vels,time,BCVars::zvel_bc); (*physbcs_w[lev-1])(vars_new[lev-1][Vars::zvel], vars_new[lev-1][Vars::xvel], vars_new[lev-1][Vars::yvel], - ngvect_vels,time,BCVars::zvel_bc); + ngvect_vels,time,BCVars::zvel_bc,true); fmf = {&vars_old[lev ][Vars::zvel], &vars_new[lev ][Vars::zvel]}; cmf = {&vars_old[lev-1][Vars::zvel], &vars_new[lev-1][Vars::zvel]}; @@ -203,10 +196,13 @@ ERF::FillPatch (int lev, Real time, int icomp_cons = 0; int ncomp_cons = mfs_vel[Vars::cons]->nComp(); + bool do_fb = true; + #ifdef ERF_USE_NETCDF // We call this here because it is an ERF routine if (use_real_bcs && (lev==0)) { fill_from_realbdy(mfs_vel,time,cons_only,icomp_cons,ncomp_cons,ngvect_cons,ngvect_vels); + do_fb = false; } #endif @@ -214,12 +210,12 @@ ERF::FillPatch (int lev, Real time, // We call these even if init_type == real because these will fill the vertical bcs // Note that we call FillBoundary inside the physbcs call - (*physbcs_cons[lev])(*mfs_vel[Vars::cons],icomp_cons,ncomp_cons,ngvect_cons,time,BCVars::cons_bc); + (*physbcs_cons[lev])(*mfs_vel[Vars::cons],icomp_cons,ncomp_cons,ngvect_cons,time,BCVars::cons_bc, do_fb); if (!cons_only) { - (*physbcs_u[lev])(*mfs_vel[Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc); - (*physbcs_v[lev])(*mfs_vel[Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc); + (*physbcs_u[lev])(*mfs_vel[Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc, do_fb); + (*physbcs_v[lev])(*mfs_vel[Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc, do_fb); (*physbcs_w[lev])(*mfs_vel[Vars::zvel],*mfs_vel[Vars::xvel],*mfs_vel[Vars::yvel], - ngvect_vels,time,BCVars::zvel_bc); + ngvect_vels,time,BCVars::zvel_bc, do_fb); } } @@ -318,8 +314,8 @@ ERF::FillIntermediatePatch (int lev, Real time, Vector ftime = {time,time}; // Impose physical bc's on coarse data (note time and 0 are not used) - (*physbcs_cons[lev-1])(vars_old[lev-1][Vars::cons],0,ncomp_cons,IntVect{ng_cons},time,BCVars::cons_bc); - (*physbcs_cons[lev-1])(vars_new[lev-1][Vars::cons],0,ncomp_cons,IntVect{ng_cons},time,BCVars::cons_bc); + (*physbcs_cons[lev-1])(vars_old[lev-1][Vars::cons],0,ncomp_cons,IntVect{ng_cons},time,BCVars::cons_bc,true); + (*physbcs_cons[lev-1])(vars_new[lev-1][Vars::cons],0,ncomp_cons,IntVect{ng_cons},time,BCVars::cons_bc,true); // Call FillPatchTwoLevels which ASSUMES that all ghost cells have already been filled mapper = &cell_cons_interp; @@ -350,8 +346,8 @@ ERF::FillIntermediatePatch (int lev, Real time, cmf = {&vars_old[lev-1][Vars::xvel], &vars_new[lev-1][Vars::xvel]}; // Impose physical bc's on coarse data (note time and 0 are not used) - (*physbcs_u[lev-1])(vars_old[lev-1][Vars::xvel],0,1,IntVect{ng_vel},time,BCVars::xvel_bc); - (*physbcs_u[lev-1])(vars_new[lev-1][Vars::xvel],0,1,IntVect{ng_vel},time,BCVars::xvel_bc); + (*physbcs_u[lev-1])(vars_old[lev-1][Vars::xvel],0,1,IntVect{ng_vel},time,BCVars::xvel_bc,true); + (*physbcs_u[lev-1])(vars_new[lev-1][Vars::xvel],0,1,IntVect{ng_vel},time,BCVars::xvel_bc,true); // Call FillPatchTwoLevels which ASSUMES that all ghost cells have already been filled FillPatchTwoLevels(mfu, IntVect{ng_vel}, IntVect(0,0,0), @@ -368,8 +364,8 @@ ERF::FillIntermediatePatch (int lev, Real time, cmf = {&vars_old[lev-1][Vars::yvel], &vars_new[lev-1][Vars::yvel]}; // Impose physical bc's on coarse data (note time and 0 are not used) - (*physbcs_v[lev-1])(vars_old[lev-1][Vars::yvel],0,1,IntVect{ng_vel},time,BCVars::yvel_bc); - (*physbcs_v[lev-1])(vars_new[lev-1][Vars::yvel],0,1,IntVect{ng_vel},time,BCVars::yvel_bc); + (*physbcs_v[lev-1])(vars_old[lev-1][Vars::yvel],0,1,IntVect{ng_vel},time,BCVars::yvel_bc,true); + (*physbcs_v[lev-1])(vars_new[lev-1][Vars::yvel],0,1,IntVect{ng_vel},time,BCVars::yvel_bc,true); // Call FillPatchTwoLevels which ASSUMES that all ghost cells have already been filled FillPatchTwoLevels(mfv, IntVect{ng_vel}, IntVect(0,0,0), @@ -389,11 +385,11 @@ ERF::FillIntermediatePatch (int lev, Real time, (*physbcs_w[lev-1])(vars_old[lev-1][Vars::zvel], vars_old[lev-1][Vars::xvel], vars_old[lev-1][Vars::yvel], - IntVect{ng_vel},time,BCVars::zvel_bc); + IntVect{ng_vel},time,BCVars::zvel_bc,true); (*physbcs_w[lev-1])(vars_new[lev-1][Vars::zvel], vars_new[lev-1][Vars::xvel], vars_new[lev-1][Vars::yvel], - IntVect{ng_vel},time,BCVars::zvel_bc); + IntVect{ng_vel},time,BCVars::zvel_bc,true); // Call FillPatchTwoLevels which ASSUMES that all ghost cells have already been filled FillPatchTwoLevels(mfw, IntVect{ng_vel}, IntVect(0,0,0), @@ -410,22 +406,25 @@ ERF::FillIntermediatePatch (int lev, Real time, IntVect ngvect_cons = IntVect(ng_cons,ng_cons,ng_cons); IntVect ngvect_vels = IntVect(ng_vel ,ng_vel ,ng_vel); + bool do_fb = true; + #ifdef ERF_USE_NETCDF // We call this here because it is an ERF routine if (use_real_bcs && (lev==0)) { fill_from_realbdy(mfs_vel,time,cons_only,icomp_cons,ncomp_cons,ngvect_cons, ngvect_vels); + do_fb = false; } #endif if (m_r2d) fill_from_bndryregs(mfs_vel,time); // We call this even if init_type == real because this routine will fill the vertical bcs - (*physbcs_cons[lev])(*mfs_vel[Vars::cons],icomp_cons,ncomp_cons,ngvect_cons,time,BCVars::cons_bc); + (*physbcs_cons[lev])(*mfs_vel[Vars::cons],icomp_cons,ncomp_cons,ngvect_cons,time,BCVars::cons_bc, do_fb); if (!cons_only) { - (*physbcs_u[lev])(*mfs_vel[Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc); - (*physbcs_v[lev])(*mfs_vel[Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc); + (*physbcs_u[lev])(*mfs_vel[Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc, do_fb); + (*physbcs_v[lev])(*mfs_vel[Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc, do_fb); (*physbcs_w[lev])(*mfs_vel[Vars::zvel],*mfs_vel[Vars::xvel],*mfs_vel[Vars::yvel], - ngvect_vels,time,BCVars::zvel_bc); + ngvect_vels,time,BCVars::zvel_bc, do_fb); } // *************************************************************************** @@ -579,21 +578,16 @@ ERF::FillCoarsePatch (int lev, Real time) domain_bcs_type); } - vars_new[lev][Vars::cons].FillBoundary(geom[lev].periodicity()); - vars_new[lev][Vars::xvel].FillBoundary(geom[lev].periodicity()); - vars_new[lev][Vars::yvel].FillBoundary(geom[lev].periodicity()); - vars_new[lev][Vars::zvel].FillBoundary(geom[lev].periodicity()); - // *************************************************************************** // Physical bc's at domain boundary // *************************************************************************** IntVect ngvect_vels = vars_new[lev][Vars::xvel].nGrowVect(); - (*physbcs_cons[lev])(vars_new[lev][Vars::cons],0,ncomp_cons,ngvect_cons,time,BCVars::cons_bc); - ( *physbcs_u[lev])(vars_new[lev][Vars::xvel],0,1 ,ngvect_vels,time,BCVars::xvel_bc); - ( *physbcs_v[lev])(vars_new[lev][Vars::yvel],0,1 ,ngvect_vels,time,BCVars::yvel_bc); + (*physbcs_cons[lev])(vars_new[lev][Vars::cons],0,ncomp_cons,ngvect_cons,time,BCVars::cons_bc,true); + ( *physbcs_u[lev])(vars_new[lev][Vars::xvel],0,1 ,ngvect_vels,time,BCVars::xvel_bc,true); + ( *physbcs_v[lev])(vars_new[lev][Vars::yvel],0,1 ,ngvect_vels,time,BCVars::yvel_bc,true); ( *physbcs_w[lev])(vars_new[lev][Vars::zvel],vars_new[lev][Vars::xvel],vars_new[lev][Vars::yvel], - ngvect_vels,time,BCVars::zvel_bc); + ngvect_vels,time,BCVars::zvel_bc,true); // *************************************************************************** // Since lev > 0 here we don't worry about m_r2d or wrfbdy data diff --git a/Source/BoundaryConditions/ERF_PhysBCFunct.H b/Source/BoundaryConditions/ERF_PhysBCFunct.H index ff09ff93a..d74a9834e 100644 --- a/Source/BoundaryConditions/ERF_PhysBCFunct.H +++ b/Source/BoundaryConditions/ERF_PhysBCFunct.H @@ -52,7 +52,8 @@ public: * @param[in] use_real_bcs if true then we fill boundary conditions for interior locations */ void operator() (amrex::MultiFab& mf, int icomp, int ncomp, - amrex::IntVect const& nghost, const amrex::Real time, int bccomp_cons); + amrex::IntVect const& nghost, const amrex::Real time, int bccomp_cons, + bool do_fb); void impose_lateral_cons_bcs (const amrex::Array4& dest_arr, const amrex::Box& bx, const amrex::Box& domain, @@ -106,7 +107,8 @@ public: * @param[in] use_real_bcs if true then we fill boundary conditions for interior locations */ void operator() (amrex::MultiFab& mf, int icomp, int ncomp, - amrex::IntVect const& nghost, const amrex::Real time, int bccomp); + amrex::IntVect const& nghost, const amrex::Real time, int bccomp, + bool do_fb); void impose_lateral_xvel_bcs (const amrex::Array4& dest_arr, const amrex::Box& bx, const amrex::Box& domain, @@ -164,7 +166,8 @@ public: * @param[in] use_real_bcs if true then we fill boundary conditions for interior locations */ void operator() (amrex::MultiFab& mf, int icomp, int ncomp, - amrex::IntVect const& nghost, const amrex::Real time, int bccomp); + amrex::IntVect const& nghost, const amrex::Real time, int bccomp, + bool do_fb); void impose_lateral_yvel_bcs (const amrex::Array4& dest_arr, const amrex::Box& bx, const amrex::Box& domain, @@ -222,7 +225,8 @@ public: * @param[in] use_real_bcs if true then we fill boundary conditions for interior locations */ void operator() (amrex::MultiFab& mf, amrex::MultiFab& xvel, amrex::MultiFab& yvel, - amrex::IntVect const& nghost, const amrex::Real time, int bccomp); + amrex::IntVect const& nghost, const amrex::Real time, int bccomp, + bool do_fb); void impose_lateral_zvel_bcs (const amrex::Array4& dest_arr, const amrex::Array4& xvel_arr, diff --git a/Source/BoundaryConditions/ERF_PhysBCFunct.cpp b/Source/BoundaryConditions/ERF_PhysBCFunct.cpp index fcfaa91b8..5ab47fb88 100644 --- a/Source/BoundaryConditions/ERF_PhysBCFunct.cpp +++ b/Source/BoundaryConditions/ERF_PhysBCFunct.cpp @@ -14,7 +14,8 @@ using namespace amrex; */ void ERFPhysBCFunct_cons::operator() (MultiFab& mf, int icomp, int ncomp, - IntVect const& nghost, const Real /*time*/, int /*bccomp*/) + IntVect const& nghost, const Real /*time*/, int /*bccomp*/, + bool do_fb) { BL_PROFILE("ERFPhysBCFunct_cons::()"); @@ -50,7 +51,9 @@ void ERFPhysBCFunct_cons::operator() (MultiFab& mf, int icomp, int ncomp, // We fill all of the interior and periodic ghost cells first, so we can fill // those directly inside the lateral and vertical calls. // - mf.FillBoundary(m_geom.periodicity()); + if (do_fb) { + mf.FillBoundary(m_geom.periodicity()); + } #ifdef AMREX_USE_OMP #pragma omp parallel if (Gpu::notInLaunchRegion()) @@ -96,7 +99,8 @@ void ERFPhysBCFunct_cons::operator() (MultiFab& mf, int icomp, int ncomp, } // operator() void ERFPhysBCFunct_u::operator() (MultiFab& mf, int /*icomp*/, int /*ncomp*/, - IntVect const& nghost, const Real time, int bccomp) + IntVect const& nghost, const Real time, int bccomp, + bool do_fb) { BL_PROFILE("ERFPhysBCFunct_u::()"); @@ -130,7 +134,9 @@ void ERFPhysBCFunct_u::operator() (MultiFab& mf, int /*icomp*/, int /*ncomp*/, // We fill all of the interior and periodic ghost cells first, so we can fill // those directly inside the lateral and vertical calls. // - mf.FillBoundary(m_geom.periodicity()); + if (do_fb) { + mf.FillBoundary(m_geom.periodicity()); + } #ifdef AMREX_USE_OMP #pragma omp parallel if (Gpu::notInLaunchRegion()) @@ -179,7 +185,8 @@ void ERFPhysBCFunct_u::operator() (MultiFab& mf, int /*icomp*/, int /*ncomp*/, } // operator() void ERFPhysBCFunct_v::operator() (MultiFab& mf, int /*icomp*/, int /*ncomp*/, - IntVect const& nghost, const Real /*time*/, int bccomp) + IntVect const& nghost, const Real /*time*/, int bccomp, + bool do_fb) { BL_PROFILE("ERFPhysBCFunct_v::()"); @@ -213,7 +220,9 @@ void ERFPhysBCFunct_v::operator() (MultiFab& mf, int /*icomp*/, int /*ncomp*/, // We fill all of the interior and periodic ghost cells first, so we can fill // those directly inside the lateral and vertical calls. // - mf.FillBoundary(m_geom.periodicity()); + if (do_fb) { + mf.FillBoundary(m_geom.periodicity()); + } #ifdef AMREX_USE_OMP #pragma omp parallel if (Gpu::notInLaunchRegion()) @@ -261,7 +270,7 @@ void ERFPhysBCFunct_v::operator() (MultiFab& mf, int /*icomp*/, int /*ncomp*/, void ERFPhysBCFunct_w::operator() (MultiFab& mf, MultiFab& xvel, MultiFab& yvel, IntVect const& nghost, const Real /*time*/, - const int bccomp_w) + const int bccomp_w, bool do_fb) { BL_PROFILE("ERFPhysBCFunct_w::()"); @@ -303,7 +312,9 @@ void ERFPhysBCFunct_w::operator() (MultiFab& mf, MultiFab& xvel, MultiFab& yvel, // We fill all of the interior and periodic ghost cells first, so we can fill // those directly inside the lateral and vertical calls. // - mf.FillBoundary(m_geom.periodicity()); + if (do_fb) { + mf.FillBoundary(m_geom.periodicity()); + } #ifdef AMREX_USE_OMP #pragma omp parallel if (Gpu::notInLaunchRegion()) From 18d47e65c2198caffa3c6ce4ee5bebf96b073722 Mon Sep 17 00:00:00 2001 From: Ann Almgren Date: Wed, 16 Oct 2024 10:23:48 -0700 Subject: [PATCH 2/2] try different https --- Docs/sphinx_doc/theory/WindFarmModels.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Docs/sphinx_doc/theory/WindFarmModels.rst b/Docs/sphinx_doc/theory/WindFarmModels.rst index 5fc8f25fa..3fff9f6dd 100644 --- a/Docs/sphinx_doc/theory/WindFarmModels.rst +++ b/Docs/sphinx_doc/theory/WindFarmModels.rst @@ -311,7 +311,7 @@ An iterative procedure is needed to compute the source terms, and is as follows: Different views of the GAD showing the forces and angles involved: Blade cross section showing the normal (:math:`V_n`) and tangential (:math:`V_t`) components of velocity with the normal (:math:`a_n`) and tangential (:math:`a_t`) induction factors, relative wind velocity :math:`V_r`, blade twist angle :math:`\xi`, angle of relative wind :math:`\psi`, blade pitch angle :math:`\phi`, lift (:math:`L`) and drag (:math:`D`) forces, and normal (:math:`F_n`) and tangential (:math:`F_t`) forces; top view showing the flow direction and inclination angle :math:`\Phi`; and front view showing the actuator disk rotating clockwise. -.. _`Mirocha et. al. 2014`: https://n2t.org/ark:/85065/d7ww7jmh +.. _`Mirocha et. al. 2014`: https://opensky.ucar.edu/islandora/object/articles:13295 .. _`Small Wind Turbines`: https://doi.org/10.1007/978-1-84996-175-2 .. _`turbine specifications`: https://github.com/NREL/openfast-turbine-models/blob/main/IEA-scaled/NREL-2.8-127/NREL-2.82-127_performance.csv .. _`details of the blade geometry`: https://github.com/NREL/openfast-turbine-models/blob/main/IEA-scaled/NREL-2.8-127/20_monolithic_opt2/OpenFAST/NREL-2p8-127_AeroDyn15_blade.dat