diff --git a/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp b/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp index 63be791f75..9aeb3dfe40 100644 --- a/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp +++ b/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp @@ -2,7 +2,7 @@ #include "amr-wind/CFDSim.H" #include "amr-wind/core/FieldUtils.H" #include "amr-wind/wind_energy/ABL.H" - +#include "amr-wind/utilities/linear_interpolation.H" #include "AMReX_ParmParse.H" namespace amr_wind::pde::icns { @@ -73,28 +73,17 @@ void ABLMeanBoussinesq::operator()( m_gravity[0], m_gravity[1], m_gravity[2]}; // Mean temperature profile used to compute background forcing term - // - // Assumes that the temperature profile is at the cell-centers of the level - // 0 grid. For finer meshes, it will extrapolate beyond the Level0 - // cell-centers for the lo/hi cells. - // + const int idir = m_axis; - const int nh_max = static_cast(m_theta_ht.size()) - 2; - const int lp1 = lev + 1; const amrex::Real* theights = m_theta_ht.data(); const amrex::Real* tvals = m_theta_vals.data(); + const amrex::Real* theights_end = m_theta_ht.end(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { - amrex::Real temp = T0; amrex::IntVect iv(i, j, k); const amrex::Real ht = problo[idir] + (iv[idir] + 0.5) * dx[idir]; - - const int il = amrex::min(k / lp1, nh_max); - const int ir = il + 1; - temp = tvals[il] + - ((tvals[ir] - tvals[il]) / (theights[ir] - theights[il])) * - (ht - theights[il]); - + const amrex::Real temp = + amr_wind::interp::linear(theights, theights_end, tvals, ht); const amrex::Real fac = beta * (temp - T0); src_term(i, j, k, 0) += gravity[0] * fac; src_term(i, j, k, 1) += gravity[1] * fac; diff --git a/amr-wind/equation_systems/icns/source_terms/ABLMesoForcingMom.cpp b/amr-wind/equation_systems/icns/source_terms/ABLMesoForcingMom.cpp index c851270f2a..239fcedf7b 100644 --- a/amr-wind/equation_systems/icns/source_terms/ABLMesoForcingMom.cpp +++ b/amr-wind/equation_systems/icns/source_terms/ABLMesoForcingMom.cpp @@ -338,9 +338,8 @@ void ABLMesoForcingMom::operator()( // averaged velocities (non tendency) const amrex::Real* vheights_begin = (m_tendency) ? m_meso_ht.data() : m_vavg_ht.data(); - const amrex::Real* vheights_end = (m_tendency) - ? m_meso_ht.data() + m_meso_ht.size() - : m_vavg_ht.data() + m_vavg_ht.size(); + const amrex::Real* vheights_end = + (m_tendency) ? m_meso_ht.end() : m_vavg_ht.end(); const amrex::Real* u_error_val = m_error_meso_avg_U.data(); const amrex::Real* v_error_val = m_error_meso_avg_V.data(); const int idir = (int)m_axis; diff --git a/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp b/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp index 94e144f046..abebb49b0a 100644 --- a/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp +++ b/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp @@ -121,10 +121,9 @@ void BodyForce::operator()( const auto& nph_time = 0.5 * (m_time.current_time() + m_time.new_time()); const auto& problo = m_mesh.Geom(lev).ProbLoArray(); const auto& dx = m_mesh.Geom(lev).CellSizeArray(); - const int lp1 = lev + 1; - const int nh_max = (int)m_prof_x.size() - 2; const amrex::Real* force_ht = m_ht.data(); + const amrex::Real* force_ht_end = m_ht.end(); const amrex::Real* force_x = m_prof_x.data(); const amrex::Real* force_y = m_prof_y.data(); @@ -134,19 +133,10 @@ void BodyForce::operator()( bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { amrex::IntVect iv(i, j, k); const amrex::Real ht = problo[2] + (iv[2] + 0.5) * dx[2]; - const int il = amrex::min(k / lp1, nh_max); - const int ir = il + 1; - amrex::Real fx; - amrex::Real fy; - - fx = force_x[il] + ((force_x[ir] - force_x[il]) / - (force_ht[ir] - force_ht[il])) * - (ht - force_ht[il]); - - fy = force_y[il] + ((force_y[ir] - force_y[il]) / - (force_ht[ir] - force_ht[il])) * - (ht - force_ht[il]); - + const amrex::Real fx = amr_wind::interp::linear( + force_ht, force_ht_end, force_x, ht); + const amrex::Real fy = amr_wind::interp::linear( + force_ht, force_ht_end, force_y, ht); src_term(i, j, k, 0) += fx; src_term(i, j, k, 1) += fy; }); diff --git a/amr-wind/equation_systems/icns/source_terms/HurricaneForcing.cpp b/amr-wind/equation_systems/icns/source_terms/HurricaneForcing.cpp index 0df4a05378..888ba0f10c 100644 --- a/amr-wind/equation_systems/icns/source_terms/HurricaneForcing.cpp +++ b/amr-wind/equation_systems/icns/source_terms/HurricaneForcing.cpp @@ -4,6 +4,7 @@ #include "amr-wind/core/vs/vstraits.H" #include "amr-wind/wind_energy/ABL.H" #include "amr-wind/core/FieldUtils.H" +#include "amr-wind/utilities/linear_interpolation.H" #include "AMReX_ParmParse.H" #include "AMReX_Gpu.H" @@ -61,32 +62,20 @@ void HurricaneForcing::operator()( const amrex::Real f = m_coriolis_factor; // Mean velocity profile used to compute background hurricane forcing term - // - // Assumes that the velocity profile is at the cell-centers of the finest - // level grid. For finer meshes, it will extrapolate beyond the Level0 - // cell-centers for the lo/hi cells. const int idir = m_axis; - const int nh_max = static_cast(m_vel_ht.size()) - 2; - const int lp1 = lev + 1; const amrex::Real* heights = m_vel_ht.data(); + const amrex::Real* heights_end = m_vel_ht.end(); const amrex::Real* vals = m_vel_vals.data(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { amrex::IntVect iv(i, j, k); const amrex::Real ht = problo[idir] + (iv[idir] + 0.5) * dx[idir]; - const int il = amrex::min(k / lp1, nh_max); - const int ir = il + 1; - const amrex::Real umean = - vals[3 * il + 0] + ((vals[3 * ir + 0] - vals[3 * il + 0]) / - (heights[ir] - heights[il])) * - (ht - heights[il]); + amr_wind::interp::linear(heights, heights_end, vals, ht, 3, 0); const amrex::Real vmean = - vals[3 * il + 1] + ((vals[3 * ir + 1] - vals[3 * il + 1]) / - (heights[ir] - heights[il])) * - (ht - heights[il]); + amr_wind::interp::linear(heights, heights_end, vals, ht, 3, 1); // Gradient velocities varying with height const amrex::Real V_z = V * (Vzh - ht) / Vzh; diff --git a/amr-wind/equation_systems/temperature/source_terms/ABLMesoForcingTemp.cpp b/amr-wind/equation_systems/temperature/source_terms/ABLMesoForcingTemp.cpp index 2f9543fc54..f9c0930410 100644 --- a/amr-wind/equation_systems/temperature/source_terms/ABLMesoForcingTemp.cpp +++ b/amr-wind/equation_systems/temperature/source_terms/ABLMesoForcingTemp.cpp @@ -302,8 +302,7 @@ void ABLMesoForcingTemp::operator()( const amrex::Real* theights_begin = (m_tendency) ? m_meso_ht.data() : m_theta_ht.data(); const amrex::Real* theights_end = - (m_tendency) ? m_meso_ht.data() + m_meso_ht.size() - : m_theta_ht.data() + m_theta_ht.size(); + (m_tendency) ? m_meso_ht.end() : m_theta_ht.end(); const amrex::Real* theta_error_val = m_error_meso_avg_theta.data(); const int idir = (int)m_axis; diff --git a/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp b/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp index ea198c5d76..ff02c41a7c 100644 --- a/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp +++ b/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp @@ -1,6 +1,7 @@ #include "amr-wind/equation_systems/temperature/source_terms/BodyForce.H" #include "amr-wind/CFDSim.H" #include "amr-wind/utilities/trig_ops.H" +#include "amr-wind/utilities/linear_interpolation.H" #include "AMReX_ParmParse.H" #include "AMReX_Gpu.H" @@ -72,10 +73,9 @@ void BodyForce::operator()( const auto& problo = m_mesh.Geom(lev).ProbLoArray(); const auto& dx = m_mesh.Geom(lev).CellSizeArray(); - const int lp1 = lev + 1; - const int nh_max = (int)m_prof_theta.size() - 2; const amrex::Real* force_ht = m_ht.data(); + const amrex::Real* force_ht_end = m_ht.end(); const amrex::Real* force_theta = m_prof_theta.data(); if (m_type == "height_varying" || m_type == "height-varying") { @@ -84,14 +84,8 @@ void BodyForce::operator()( bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { amrex::IntVect iv(i, j, k); const amrex::Real ht = problo[2] + (iv[2] + 0.5) * dx[2]; - const int il = amrex::min(k / lp1, nh_max); - const int ir = il + 1; - amrex::Real ftheta; - - ftheta = - force_theta[il] + ((force_theta[ir] - force_theta[il]) / - (force_ht[ir] - force_ht[il])) * - (ht - force_ht[il]); + const amrex::Real ftheta = amr_wind::interp::linear( + force_ht, force_ht_end, force_theta, ht); src_term(i, j, k, 0) += ftheta; }); diff --git a/amr-wind/equation_systems/temperature/source_terms/HurricaneTempForcing.cpp b/amr-wind/equation_systems/temperature/source_terms/HurricaneTempForcing.cpp index f41ef06974..d74b2eaadd 100644 --- a/amr-wind/equation_systems/temperature/source_terms/HurricaneTempForcing.cpp +++ b/amr-wind/equation_systems/temperature/source_terms/HurricaneTempForcing.cpp @@ -4,6 +4,7 @@ #include "amr-wind/core/vs/vstraits.H" #include "amr-wind/wind_energy/ABL.H" #include "amr-wind/core/FieldUtils.H" +#include "amr-wind/utilities/linear_interpolation.H" #include "AMReX_ParmParse.H" #include "AMReX_Gpu.H" @@ -44,24 +45,16 @@ void HurricaneTempForcing::operator()( const amrex::Real dTzh = m_dTzh; // Mean velocity profile used to compute background hurricane forcing term - // - // Assumes that the velocity profile is at the cell-centers of the finest - // level grid. For finer meshes, it will extrapolate beyond the Level0 - // cell-centers for the lo/hi cells. const int idir = m_axis; - const int nh_max = static_cast(m_vel_ht.size()) - 2; - const int lp1 = lev + 1; const amrex::Real* heights = m_vel_ht.data(); + const amrex::Real* heights_end = m_vel_ht.end(); const amrex::Real* vals = m_vel_vals.data(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { amrex::IntVect iv(i, j, k); const amrex::Real ht = problo[idir] + (iv[idir] + 0.5) * dx[idir]; - const int il = amrex::min(k / lp1, nh_max); - const int ir = il + 1; - /*const amrex::Real umean = vals[3 * il + 0] + ((vals[3 * ir + 0] - vals[3 * il + 0]) / (heights[ir] - heights[il])) * @@ -69,9 +62,7 @@ void HurricaneTempForcing::operator()( */ const amrex::Real dTdR_z = dTdR * (dTzh - ht) / dTzh; const amrex::Real vmean = - vals[3 * il + 1] + ((vals[3 * ir + 1] - vals[3 * il + 1]) / - (heights[ir] - heights[il])) * - (ht - heights[il]); + amr_wind::interp::linear(heights, heights_end, vals, ht, 3, 1); src_term(i, j, k) -= vmean * dTdR_z; }); diff --git a/amr-wind/turbulence/LES/AMD.cpp b/amr-wind/turbulence/LES/AMD.cpp index 968471ed81..0cef614881 100644 --- a/amr-wind/turbulence/LES/AMD.cpp +++ b/amr-wind/turbulence/LES/AMD.cpp @@ -86,8 +86,7 @@ void AMD::update_turbulent_viscosity( tpa_deriv_d.begin()); const amrex::Real* p_tpa_coord_begin = tpa_coord_d.data(); - const amrex::Real* p_tpa_coord_end = - tpa_coord_d.data() + tpa_coord_d.size(); + const amrex::Real* p_tpa_coord_end = tpa_coord_d.end(); const amrex::Real* p_tpa_deriv = tpa_deriv_d.data(); const int nlevels = repo.num_active_levels(); for (int lev = 0; lev < nlevels; ++lev) {