diff --git a/amr-wind/CFDSim.H b/amr-wind/CFDSim.H index 12a7299e95..c9ce4ab364 100644 --- a/amr-wind/CFDSim.H +++ b/amr-wind/CFDSim.H @@ -117,6 +117,13 @@ public: bool has_mesh_mapping() const { return m_mesh_mapping; } + void set_during_overset_advance(const bool flag) + { + m_during_overset_advance = flag; + } + + bool is_during_overset_advance() const { return m_during_overset_advance; } + bool is_anelastic() const; private: @@ -145,6 +152,9 @@ private: std::unique_ptr m_helics; bool m_mesh_mapping{false}; + + // State of solver - know if during an overset timestep or not + bool m_during_overset_advance{false}; }; } // namespace amr_wind diff --git a/amr-wind/incflo.cpp b/amr-wind/incflo.cpp index 9139abc7a1..335c4c18ef 100644 --- a/amr-wind/incflo.cpp +++ b/amr-wind/incflo.cpp @@ -267,6 +267,10 @@ void incflo::post_advance_work() if (m_time.write_checkpoint()) { m_sim.io_manager().write_checkpoint_file(); } + + if (m_sim.has_overset()) { + m_sim.set_during_overset_advance(false); + } } /** Perform time-integration for user-defined time or timesteps. @@ -333,6 +337,7 @@ void incflo::do_advance(const int fixed_point_iteration) { if (m_sim.has_overset()) { m_ovst_ops.pre_advance_work(); + m_sim.set_during_overset_advance(true); } if (m_prescribe_vel && fixed_point_iteration == 0) { prescribe_advance(); diff --git a/amr-wind/overset/OversetOps.cpp b/amr-wind/overset/OversetOps.cpp index 0ef0671e24..4702b6bc33 100644 --- a/amr-wind/overset/OversetOps.cpp +++ b/amr-wind/overset/OversetOps.cpp @@ -6,6 +6,8 @@ #include "amr-wind/core/MLMGOptions.H" #include "amr-wind/projection/nodal_projection_ops.H" #include +#include "amr-wind/wind_energy/ABL.H" +#include "amr-wind/wind_energy/ABLBoundaryPlane.H" namespace amr_wind { @@ -93,6 +95,12 @@ void OversetOps::pre_advance_work() (m_gp_copy)->num_grow()); } } + + // Pre advance work for plane was skipped for overset solver, do it here + if (m_sim_ptr->physics_manager().contains("ABL")) { + auto& abl = m_sim_ptr->physics_manager().get(); + abl.bndry_plane().pre_advance_work(); + } } void OversetOps::update_gradp() diff --git a/amr-wind/overset/TiogaInterface.cpp b/amr-wind/overset/TiogaInterface.cpp index 138b1b741e..652414256c 100644 --- a/amr-wind/overset/TiogaInterface.cpp +++ b/amr-wind/overset/TiogaInterface.cpp @@ -192,13 +192,17 @@ void TiogaInterface::register_solution( m_cell_vars = cell_vars; m_node_vars = node_vars; + const amrex::Real time_fillpatch = m_sim.is_during_overset_advance() + ? m_sim.time().new_time() + : m_sim.time().current_time(); + // Move cell variables into scratch field { int icomp = 0; for (const auto& cvar : m_cell_vars) { auto& fld = repo.get_field(cvar); const int ncomp = fld.num_comp(); - fld.fillpatch(m_sim.time().new_time()); + fld.fillpatch(time_fillpatch); field_ops::copy(*m_qcell, fld, 0, icomp, ncomp, num_ghost); icomp += ncomp; } @@ -210,7 +214,7 @@ void TiogaInterface::register_solution( for (const auto& cvar : m_cell_vars) { auto& fld = repo.get_field(cvar); const int ncomp = fld.num_comp(); - fld.fillpatch(m_sim.time().new_time()); + fld.fillpatch(time_fillpatch); // Device to host copy happens here const int nlevels = repo.num_active_levels(); for (int lev = 0; lev < nlevels; ++lev) { @@ -227,7 +231,7 @@ void TiogaInterface::register_solution( auto& fld = repo.get_field(cvar); const int ncomp = fld.num_comp(); - fld.fillpatch(m_sim.time().new_time()); + fld.fillpatch(time_fillpatch); field_ops::copy(*m_qnode, fld, 0, icomp, ncomp, num_ghost); icomp += ncomp; } @@ -239,7 +243,7 @@ void TiogaInterface::register_solution( for (const auto& cvar : m_node_vars) { auto& fld = repo.get_field(cvar); const int ncomp = fld.num_comp(); - fld.fillpatch(m_sim.time().new_time()); + fld.fillpatch(time_fillpatch); // Device to host copy happens here const int nlevels = repo.num_active_levels(); for (int lev = 0; lev < nlevels; ++lev) { @@ -291,6 +295,10 @@ void TiogaInterface::update_solution() { auto& repo = m_sim.repo(); + const amrex::Real time_fillpatch = m_sim.is_during_overset_advance() + ? m_sim.time().new_time() + : m_sim.time().current_time(); + // Update cell variables on device { int icomp = 0; @@ -302,7 +310,7 @@ void TiogaInterface::update_solution() for (int lev = 0; lev < nlevels; ++lev) { htod_memcpy(fld(lev), (*m_qcell_host)(lev), icomp, 0, ncomp); } - fld.fillpatch(m_sim.time().new_time()); + fld.fillpatch(time_fillpatch); icomp += ncomp; } } @@ -318,7 +326,7 @@ void TiogaInterface::update_solution() for (int lev = 0; lev < nlevels; ++lev) { htod_memcpy(fld(lev), (*m_qnode_host)(lev), icomp, 0, ncomp); } - fld.fillpatch(m_sim.time().new_time()); + fld.fillpatch(time_fillpatch); icomp += ncomp; } } diff --git a/amr-wind/wind_energy/ABL.H b/amr-wind/wind_energy/ABL.H index f03fa227f6..68fd28f592 100644 --- a/amr-wind/wind_energy/ABL.H +++ b/amr-wind/wind_energy/ABL.H @@ -115,6 +115,7 @@ public: } const ABLBoundaryPlane& bndry_plane() const { return *m_bndry_plane; } + ABLBoundaryPlane& bndry_plane() { return *m_bndry_plane; } const ABLModulatedPowerLaw& abl_mpl() const { return *m_abl_mpl; } //! Return the ABL statistics calculator diff --git a/amr-wind/wind_energy/ABL.cpp b/amr-wind/wind_energy/ABL.cpp index 7c50d78d4b..8504e7b595 100644 --- a/amr-wind/wind_energy/ABL.cpp +++ b/amr-wind/wind_energy/ABL.cpp @@ -238,7 +238,9 @@ void ABL::pre_advance_work() m_stats->vel_profile_coarse()); } - m_bndry_plane->pre_advance_work(); + if (!m_sim.has_overset()) { + m_bndry_plane->pre_advance_work(); + } m_abl_mpl->pre_advance_work(); }