From 9b6e0895db1df02e3cc43a98354a45afc725773e Mon Sep 17 00:00:00 2001 From: Michael Kuhn Date: Fri, 20 Jun 2025 16:01:39 -0600 Subject: [PATCH 1/3] make w2a initialization more accurate, especially for low resolutions --- .../relaxation_zones/waves2amr_ops.H | 68 ++++++++++++++----- 1 file changed, 52 insertions(+), 16 deletions(-) diff --git a/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H b/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H index ad1a5581d2..1dc8157159 100644 --- a/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H +++ b/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H @@ -422,9 +422,8 @@ struct ReadInputsOp } // Set up plan for FFTW if (wdata.is_ocean) { - wdata.plan_vector.emplace_back( - modes_hosgrid::plan_ifftw( - wdata.n0, wdata.n1, wdata.c_eta_mptr, plan_f)); + wdata.plan_vector.emplace_back(modes_hosgrid::plan_ifftw( + wdata.n0, wdata.n1, wdata.c_eta_mptr, plan_f)); } else { modes_hosgrid::plan_ifftw_nwt( wdata.n0, wdata.n1, wdata.plan_vector, wdata.r_eta_mptr, @@ -589,10 +588,8 @@ struct InitDataOp void // cppcheck-suppress constParameterReference operator()( - W2AWaves::DataType& data, - int level, - const amrex::Geometry& geom, - bool multiphase_mode) + W2AWaves::DataType & + data, int level, const amrex::Geometry & geom, bool multiphase_mode) { #ifdef AMR_WIND_USE_W2A @@ -652,6 +649,16 @@ struct InitDataOp const bool has_beach = wdata.has_beach && multiphase_mode; const bool init_wave_field = wdata.init_wave_field || !multiphase_mode; + amrex::MultiFab phi_tmp_fab( + velocity(level).boxArray(), velocity(level).DistributionMap(), 1, + 3); + amrex::MultiFab vel_tmp_fab( + velocity(level).boxArray(), velocity(level).DistributionMap(), + AMREX_SPACEDIM, 3); + + const auto& phi_tmp = phi_tmp_fab.arrays(); + const auto& vel_tmp = vel_tmp_fab.arrays(); + amrex::ParallelFor( velocity(level), amrex::IntVect(3), [=] AMREX_GPU_DEVICE(int nbx, int i, int j, int k) noexcept { @@ -673,17 +680,13 @@ struct InitDataOp x, problo[0], gen_length, probhi[0], beach_length, wave_sol, bulk, outlet); - const amrex::Real phi = local_profile[3] - z; - const amrex::Real cell_length_2D = - std::sqrt(dx[0] * dx[0] + dx[2] * dx[2]); - if (phi + cell_length_2D >= 0) { - vel[nbx](i, j, k, 0) = local_profile[0]; - vel[nbx](i, j, k, 1) = local_profile[1]; - vel[nbx](i, j, k, 2) = local_profile[2]; - } + phi_tmp[nbx](i, j, k) = local_profile[3] - z; + vel_tmp[nbx](i, j, k, 0) = local_profile[0]; + vel_tmp[nbx](i, j, k, 1) = local_profile[1]; + vel_tmp[nbx](i, j, k, 2) = local_profile[2]; if (multiphase_mode) { - phi_arrs[nbx](i, j, k) = phi; + phi_arrs[nbx](i, j, k) = phi_tmp[nbx](i, j, k); } // Default w2a values matter for where no updates happen @@ -693,6 +696,39 @@ struct InitDataOp w2a_vel[nbx](i, j, k, 2) = quiescent[0]; }); amrex::Gpu::streamSynchronize(); + + amrex::ParallelFor( + velocity(level), amrex::IntVect(2), + [=] AMREX_GPU_DEVICE(int nbx, int i, int j, int k) noexcept { + const amrex::Real eps = 2. * std::cbrt(dx[0] * dx[1] * dx[2]); + const amrex::Real vof_local = + multiphase::levelset_to_vof(i, j, k, eps, phi_tmp[nbx]); + + if (vof_local > 1. - constants::TIGHT_TOL) { + // Entire cell is liquid + vel[nbx](i, j, k, 0) = vel_tmp[nbx](i, j, k, 0); + vel[nbx](i, j, k, 1) = vel_tmp[nbx](i, j, k, 1); + vel[nbx](i, j, k, 2) = vel_tmp[nbx](i, j, k, 2); + } else if (vof_local > 1e-3) { + // Velocity of liquid becomes velocity of entire cell + // Interpolate using cell-centered values to liquid centroid + // Consider only z-direction to find centroid + const amrex::Real z_cell = problo[2] + (k + 0.5) * dx[2]; + const amrex::Real z_cell_below = z_cell - dx[2]; + const amrex::Real z_cell_bottom = z_cell - 0.5 * dx[2]; + const amrex::Real z_liq_center = + z_cell_bottom + 0.5 * vof_local * dx[2]; + for (int n = 0; n < AMREX_SPACEDIM; ++n) { + vel[nbx](i, j, k, n) = + vel_tmp[nbx](i, j, k - 1, n) + + (vel_tmp[nbx](i, j, k, n) - + vel_tmp[nbx](i, j, k - 1, n)) * + ((z_liq_center - z_cell_below) / + (z_cell - z_cell_below)); + } + } + }); + amrex::Gpu::streamSynchronize(); #else amrex::ignore_unused(data, level, geom, multiphase_mode); #endif From c080e9466026d340dc3a33d84eeeff5ec72cd881 Mon Sep 17 00:00:00 2001 From: Michael Kuhn Date: Mon, 23 Jun 2025 15:56:10 -0600 Subject: [PATCH 2/3] formatting --- amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H b/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H index 1dc8157159..7fae582d04 100644 --- a/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H +++ b/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H @@ -422,8 +422,9 @@ struct ReadInputsOp } // Set up plan for FFTW if (wdata.is_ocean) { - wdata.plan_vector.emplace_back(modes_hosgrid::plan_ifftw( - wdata.n0, wdata.n1, wdata.c_eta_mptr, plan_f)); + wdata.plan_vector.emplace_back( + modes_hosgrid::plan_ifftw( + wdata.n0, wdata.n1, wdata.c_eta_mptr, plan_f)); } else { modes_hosgrid::plan_ifftw_nwt( wdata.n0, wdata.n1, wdata.plan_vector, wdata.r_eta_mptr, @@ -588,8 +589,10 @@ struct InitDataOp void // cppcheck-suppress constParameterReference operator()( - W2AWaves::DataType & - data, int level, const amrex::Geometry & geom, bool multiphase_mode) + W2AWaves::DataType& data, + int level, + const amrex::Geometry& geom, + bool multiphase_mode) { #ifdef AMR_WIND_USE_W2A From 4058e2e5facb0e3faa844532d678506835df56c3 Mon Sep 17 00:00:00 2001 From: Michael Kuhn Date: Mon, 21 Jul 2025 10:20:51 -0600 Subject: [PATCH 3/3] avoid dividing by 0, just in case --- amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H b/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H index 7fae582d04..4bcda13078 100644 --- a/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H +++ b/amr-wind/ocean_waves/relaxation_zones/waves2amr_ops.H @@ -727,7 +727,7 @@ struct InitDataOp (vel_tmp[nbx](i, j, k, n) - vel_tmp[nbx](i, j, k - 1, n)) * ((z_liq_center - z_cell_below) / - (z_cell - z_cell_below)); + (z_cell - z_cell_below + constants::EPS)); } } });