From 05e7c8e22994d39efe5a13562c420af06c9951bd Mon Sep 17 00:00:00 2001 From: Prakash Mohan Date: Fri, 8 Nov 2024 13:50:29 -0700 Subject: [PATCH 1/6] use interp::linear instead of hand rolling the interpolation --- .../icns/source_terms/ABLMeanBoussinesq.cpp | 21 +++++-------------- 1 file changed, 5 insertions(+), 16 deletions(-) diff --git a/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp b/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp index 63be791f75..b1190c9524 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.data() + m_theta_ht.size(); 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; From db0d91218babe93ee71dcee53aaf8e3b00863fc5 Mon Sep 17 00:00:00 2001 From: Prakash Mohan Date: Tue, 12 Nov 2024 14:44:19 -0700 Subject: [PATCH 2/6] formatting changes --- .../icns/source_terms/ABLMeanBoussinesq.cpp | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp b/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp index b1190c9524..581235e441 100644 --- a/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp +++ b/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp @@ -73,7 +73,7 @@ void ABLMeanBoussinesq::operator()( m_gravity[0], m_gravity[1], m_gravity[2]}; // Mean temperature profile used to compute background forcing term - + const int idir = m_axis; const amrex::Real* theights = m_theta_ht.data(); const amrex::Real* tvals = m_theta_vals.data(); @@ -82,8 +82,8 @@ void ABLMeanBoussinesq::operator()( 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 amrex::Real temp = amr_wind::interp::linear( - theights, theights_end, tvals, ht); + 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; From 470ee05ddd48afcff1aa0e64b770b1467cbb42b5 Mon Sep 17 00:00:00 2001 From: prakash Date: Wed, 13 Nov 2024 08:14:53 -0700 Subject: [PATCH 3/6] Replace manual interp with utils function --- .../icns/source_terms/BodyForce.cpp | 20 +++++-------------- .../icns/source_terms/HurricaneForcing.cpp | 19 ++++-------------- .../temperature/source_terms/BodyForce.cpp | 14 ++++--------- .../source_terms/HurricaneTempForcing.cpp | 15 +++----------- 4 files changed, 16 insertions(+), 52 deletions(-) diff --git a/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp b/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp index 94e144f046..f015abb378 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.data() + m_ht.size(); 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..f6d1049fd6 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.data() + m_vel_ht.size(); 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/BodyForce.cpp b/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp index ea198c5d76..9731c27027 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.data() + m_ht.size(); 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..8053495a65 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.data() + m_vel_ht.size(); 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; }); From dd4bcadf2adb1d0b9ef674173a46cd1cc8cfe3ae Mon Sep 17 00:00:00 2001 From: prakash Date: Thu, 14 Nov 2024 12:57:03 -0700 Subject: [PATCH 4/6] replace data() + size() with end() --- .../equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp | 2 +- amr-wind/equation_systems/icns/source_terms/BodyForce.cpp | 2 +- .../equation_systems/icns/source_terms/HurricaneForcing.cpp | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp b/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp index 581235e441..9aeb3dfe40 100644 --- a/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp +++ b/amr-wind/equation_systems/icns/source_terms/ABLMeanBoussinesq.cpp @@ -77,7 +77,7 @@ void ABLMeanBoussinesq::operator()( const int idir = m_axis; const amrex::Real* theights = m_theta_ht.data(); const amrex::Real* tvals = m_theta_vals.data(); - const amrex::Real* theights_end = m_theta_ht.data() + m_theta_ht.size(); + const amrex::Real* theights_end = m_theta_ht.end(); amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept { amrex::IntVect iv(i, j, k); diff --git a/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp b/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp index f015abb378..abebb49b0a 100644 --- a/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp +++ b/amr-wind/equation_systems/icns/source_terms/BodyForce.cpp @@ -123,7 +123,7 @@ void BodyForce::operator()( const auto& dx = m_mesh.Geom(lev).CellSizeArray(); const amrex::Real* force_ht = m_ht.data(); - const amrex::Real* force_ht_end = m_ht.data() + m_ht.size(); + 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(); diff --git a/amr-wind/equation_systems/icns/source_terms/HurricaneForcing.cpp b/amr-wind/equation_systems/icns/source_terms/HurricaneForcing.cpp index f6d1049fd6..888ba0f10c 100644 --- a/amr-wind/equation_systems/icns/source_terms/HurricaneForcing.cpp +++ b/amr-wind/equation_systems/icns/source_terms/HurricaneForcing.cpp @@ -65,7 +65,7 @@ void HurricaneForcing::operator()( const int idir = m_axis; const amrex::Real* heights = m_vel_ht.data(); - const amrex::Real* heights_end = m_vel_ht.data() + m_vel_ht.size(); + 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 { From d31dc1c25da5dd622b651afa1617372855f3407c Mon Sep 17 00:00:00 2001 From: prakash Date: Thu, 14 Nov 2024 13:54:32 -0700 Subject: [PATCH 5/6] more places of data + size vs end --- .../equation_systems/temperature/source_terms/BodyForce.cpp | 2 +- .../temperature/source_terms/HurricaneTempForcing.cpp | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp b/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp index 9731c27027..ff02c41a7c 100644 --- a/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp +++ b/amr-wind/equation_systems/temperature/source_terms/BodyForce.cpp @@ -75,7 +75,7 @@ void BodyForce::operator()( const auto& dx = m_mesh.Geom(lev).CellSizeArray(); const amrex::Real* force_ht = m_ht.data(); - const amrex::Real* force_ht_end = m_ht.data() + m_ht.size(); + 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") { diff --git a/amr-wind/equation_systems/temperature/source_terms/HurricaneTempForcing.cpp b/amr-wind/equation_systems/temperature/source_terms/HurricaneTempForcing.cpp index 8053495a65..d74b2eaadd 100644 --- a/amr-wind/equation_systems/temperature/source_terms/HurricaneTempForcing.cpp +++ b/amr-wind/equation_systems/temperature/source_terms/HurricaneTempForcing.cpp @@ -48,7 +48,7 @@ void HurricaneTempForcing::operator()( const int idir = m_axis; const amrex::Real* heights = m_vel_ht.data(); - const amrex::Real* heights_end = m_vel_ht.data() + m_vel_ht.size(); + 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 { From 4e2b21e76290ca8c37884543154bc362cd8736ca Mon Sep 17 00:00:00 2001 From: prakash Date: Thu, 14 Nov 2024 14:06:26 -0700 Subject: [PATCH 6/6] Even more places where data() + size() is used instead of end() --- .../equation_systems/icns/source_terms/ABLMesoForcingMom.cpp | 5 ++--- .../temperature/source_terms/ABLMesoForcingTemp.cpp | 3 +-- amr-wind/turbulence/LES/AMD.cpp | 3 +-- 3 files changed, 4 insertions(+), 7 deletions(-) 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/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/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) {