Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Improvements to wave generation - ocean waves #1264

Merged
merged 8 commits into from
Oct 2, 2024
+81 −92
Merged

Improvements to wave generation - ocean waves #1264

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
fd63acb
put loops where they should be already
mbkuhn Sep 25, 2024
6e1f244
set up rampf properly for generation region
mbkuhn Sep 25, 2024
2de6c1a
duplicate changes for outprofile
mbkuhn Sep 25, 2024
e334a32
adding limitation to velocity forcing, rephrasing things
mbkuhn Sep 27, 2024
ce064e5
tweak for when liquid is added
mbkuhn Sep 27, 2024
e874c36
duplicate for outlet region
mbkuhn Sep 27, 2024
ba2ebf2
remove false default, update docs
mbkuhn Oct 2, 2024
9517ae3
correcting/updating some stuff in the reg tests while we're at it
mbkuhn Oct 2, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion amr-wind/ocean_waves/relaxation_zones/RelaxationZones.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -34,7 +34,7 @@ struct RelaxZonesBaseData
bool has_beach{true};
bool has_outprofile{false};

amrex::Real ramp_period{2.0};
amrex::Real ramp_period;

bool regrid_occurred{false};
};
Expand Down
133 changes: 72 additions & 61 deletions amr-wind/ocean_waves/relaxation_zones/relaxation_zones_ops.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ void read_inputs(

wdata.has_ramp = pp.contains("timeramp_period");
if (wdata.has_ramp) {
pp.query("timeramp_period", wdata.ramp_period);
pp.get("timeramp_period", wdata.ramp_period);
}

amrex::ParmParse pp_multiphase("MultiPhase");
Expand All @@ -59,6 +59,7 @@ void apply_relaxation_zones(CFDSim& sim, const RelaxZonesBaseData& wdata)
const auto& mphase = sim.physics_manager().get<MultiPhase>();
const amrex::Real rho1 = mphase.rho1();
const amrex::Real rho2 = mphase.rho2();
constexpr amrex::Real vof_tiny = 1e-12;

for (int lev = 0; lev < nlevels; ++lev) {
auto& ls = m_ow_levelset(lev);
Expand Down Expand Up @@ -119,40 +120,47 @@ void apply_relaxation_zones(CFDSim& sim, const RelaxZonesBaseData& wdata)
if (x <= problo[0] + gen_length) {
const amrex::Real Gamma =
utils::gamma_generate(x - problo[0], gen_length);
const amrex::Real vf =
(1. - Gamma) * target_volfrac(i, j, k) * rampf +
// Get bounded new vof, incorporate with increment
amrex::Real new_vof =
(1. - Gamma) * target_volfrac(i, j, k) +
Gamma * volfrac(i, j, k);
volfrac(i, j, k) = (vf > 1. - 1.e-10) ? 1.0 : vf;
// Force liquid velocity, update according to mom.
new_vof = (new_vof > 1. - vof_tiny)
? 1.0
: (new_vof < vof_tiny ? 0.0 : new_vof);
const amrex::Real dvf = new_vof - volfrac(i, j, k);
volfrac(i, j, k) += rampf * dvf;
// Force liquid velocity only if target vof present
const amrex::Real fvel_liq =
(target_volfrac(i, j, k) > vof_tiny) ? 1.0 : 0.0;
amrex::Real rho_ = rho1 * volfrac(i, j, k) +
rho2 * (1.0 - volfrac(i, j, k));
vel(i, j, k, 0) =
(rho1 * volfrac(i, j, k) *
(rampf * (1. - Gamma) *
target_vel(i, j, k, 0) +
Gamma * vel(i, j, k, 0)) +
rho2 * (1. - volfrac(i, j, k)) * vel(i, j, k, 0)) /
rho_;
vel(i, j, k, 1) =
(rho1 * volfrac(i, j, k) *
(rampf * (1. - Gamma) *
target_vel(i, j, k, 1) +
Gamma * vel(i, j, k, 1)) +
rho2 * (1. - volfrac(i, j, k)) * vel(i, j, k, 1)) /
rho_;
vel(i, j, k, 2) =
(rho1 * volfrac(i, j, k) *
(rampf * (1. - Gamma) *
target_vel(i, j, k, 2) +
Gamma * vel(i, j, k, 2)) +
rho2 * (1. - volfrac(i, j, k)) * vel(i, j, k, 2)) /
rho_;
for (int n = 0; n < vel.ncomp; ++n) {
// Get updated liquid velocity
amrex::Real vel_liq = vel(i, j, k, n);
const amrex::Real dvel_liq =
((1. - Gamma) * target_vel(i, j, k, n) +
Gamma * vel_liq) -
vel_liq;
vel_liq += rampf * fvel_liq * dvel_liq;
// If liquid was added, that liquid has target_vel
amrex::Real integrated_vel_liq =
volfrac(i, j, k) * vel_liq;
integrated_vel_liq +=
rampf * fvel_liq * amrex::max(0.0, dvf) *
(target_vel(i, j, k, n) - vel(i, j, k, n));
// Update overall velocity using momentum
vel(i, j, k, n) = (rho1 * integrated_vel_liq +
rho2 * (1. - volfrac(i, j, k)) *
vel(i, j, k, n)) /
rho_;
}
}
// Numerical beach (sponge layer)
// Outlet region
if (x + beach_length >= probhi[0]) {
const amrex::Real Gamma = utils::gamma_absorb(
x - (probhi[0] - beach_length), beach_length,
beach_length_factor);
// Numerical beach (sponge layer)
if (has_beach) {
volfrac(i, j, k) =
(1.0 - Gamma) *
Expand All @@ -162,45 +170,48 @@ void apply_relaxation_zones(CFDSim& sim, const RelaxZonesBaseData& wdata)
amrex::Real rho_ = rho1 * volfrac(i, j, k) +
rho2 * (1.0 - volfrac(i, j, k));
// Target solution in liquid is vel = 0
vel(i, j, k, 0) = (rho1 * volfrac(i, j, k) * Gamma +
rho2 * (1. - volfrac(i, j, k))) *
vel(i, j, k, 0) / rho_;
vel(i, j, k, 1) = (rho1 * volfrac(i, j, k) * Gamma +
rho2 * (1. - volfrac(i, j, k))) *
vel(i, j, k, 1) / rho_;
vel(i, j, k, 2) = (rho1 * volfrac(i, j, k) * Gamma +
rho2 * (1. - volfrac(i, j, k))) *
vel(i, j, k, 2) / rho_;
for (int n = 0; n < vel.ncomp; ++n) {
vel(i, j, k, n) =
(rho1 * volfrac(i, j, k) * Gamma +
rho2 * (1. - volfrac(i, j, k))) *
vel(i, j, k, n) / rho_;
}
}
// Forcing to wave profile instead
if (has_outprofile) {
const amrex::Real vf =
(1. - Gamma) * target_volfrac(i, j, k) * rampf +
// Same steps as in wave generation region
amrex::Real new_vof =
(1. - Gamma) * target_volfrac(i, j, k) +
Gamma * volfrac(i, j, k);
volfrac(i, j, k) = (vf > 1. - 1.e-10) ? 1.0 : vf;
// Force liquid velocity, update according to mom.
new_vof =
(new_vof > 1. - vof_tiny)
? 1.0
: (new_vof < vof_tiny ? 0.0 : new_vof);
const amrex::Real dvf = new_vof - volfrac(i, j, k);
volfrac(i, j, k) += dvf;
const amrex::Real fvel_liq =
(target_volfrac(i, j, k) > vof_tiny) ? 1.0
: 0.0;
amrex::Real rho_ = rho1 * volfrac(i, j, k) +
rho2 * (1.0 - volfrac(i, j, k));
vel(i, j, k, 0) = (rho1 * volfrac(i, j, k) *
(rampf * (1. - Gamma) *
target_vel(i, j, k, 0) +
Gamma * vel(i, j, k, 0)) +
rho2 * (1. - volfrac(i, j, k)) *
vel(i, j, k, 0)) /
rho_;
vel(i, j, k, 1) = (rho1 * volfrac(i, j, k) *
(rampf * (1. - Gamma) *
target_vel(i, j, k, 1) +
Gamma * vel(i, j, k, 1)) +
rho2 * (1. - volfrac(i, j, k)) *
vel(i, j, k, 1)) /
rho_;
vel(i, j, k, 2) = (rho1 * volfrac(i, j, k) *
(rampf * (1. - Gamma) *
target_vel(i, j, k, 2) +
Gamma * vel(i, j, k, 2)) +
rho2 * (1. - volfrac(i, j, k)) *
vel(i, j, k, 2)) /
rho_;
for (int n = 0; n < vel.ncomp; ++n) {
amrex::Real vel_liq = vel(i, j, k, n);
const amrex::Real dvel_liq =
((1. - Gamma) * target_vel(i, j, k, n) +
Gamma * vel_liq) -
vel_liq;
vel_liq += rampf * fvel_liq * dvel_liq;
amrex::Real integrated_vel_liq =
volfrac(i, j, k) * vel_liq;
integrated_vel_liq +=
rampf * fvel_liq * amrex::max(0.0, dvf) *
(target_vel(i, j, k, n) - vel(i, j, k, n));
vel(i, j, k, n) =
(rho1 * integrated_vel_liq +
rho2 * (1. - volfrac(i, j, k)) *
vel(i, j, k, n)) /
rho_;
}
}
}

Expand Down
2 changes: 1 addition & 1 deletion docs/sphinx/user/inputs_ocean_waves.rst
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -76,7 +76,7 @@ This section is for setting up wave forcing and relaxation zones.

.. input_param:: OceanWaves.label.timeramp_period

**type:** Real, optional, default = 2.0
**type:** Real, optional

An initial ramp-up period for the wave forcing. Without specifying a period, the wave
forcing will begin at full strength.
Expand Down
16 changes: 3 additions & 13 deletions test/test_files/ow_linear/ow_linear.inp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,9 +20,6 @@ time.checkpoint_interval = -1 # Steps between checkpoint files
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# PHYSICS #
#.......................................#
incflo.use_godunov = 1
incflo.diffusion_type = 2
incflo.godunov_type="weno_z"
transport.model = TwoPhaseTransport
transport.viscosity_fluid1=0.0
transport.viscosity_fluid2=0.0
Expand All @@ -39,8 +36,6 @@ OceanWaves.Wave1.relax_zone_gen_length=4.0
OceanWaves.Wave1.numerical_beach_length=8.0
MultiPhase.density_fluid1=1000.
MultiPhase.density_fluid2=1.
MultiPhase.interface_smoothing=0
MultiPhase.interface_smoothing_frequency=1
ICNS.source_terms = GravityForcing
ICNS.use_perturb_pressure = true
ICNS.reconstruct_true_pressure = true
Expand All @@ -57,15 +52,10 @@ geometry.prob_lo = 0.0 0. -1 # Lo corner coordinates
geometry.prob_hi = 30. 1. 1 # Hi corner coordinates
geometry.is_periodic = 0 1 0 # Periodicity x y z (0/1)

xlo.type = "pressure_inflow"
xlo.type = "wave_generation"
xhi.type = "pressure_outflow"
xlo.vof_type = "zero_gradient"
xhi.vof_type = "zero_gradient"

zlo.type = "slip_wall"
zhi.type = "pressure_outflow"
zlo.vof_type = "zero_gradient"
zhi.vof_type = "zero_gradient"

incflo.verbose=1
zhi.type = "slip_wall"

incflo.verbose=1
16 changes: 3 additions & 13 deletions test/test_files/ow_stokes/ow_stokes.inp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,9 +20,6 @@ time.checkpoint_interval = -1 # Steps between checkpoint files
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# PHYSICS #
#.......................................#
incflo.use_godunov = 1
incflo.diffusion_type = 2
incflo.godunov_type="weno_z"
transport.model = TwoPhaseTransport
transport.viscosity_fluid1=0.0
transport.viscosity_fluid2=0.0
Expand All @@ -41,8 +38,6 @@ OceanWaves.Wave1.numerical_beach_length=2.0
OceanWaves.Wave1.numerical_beach_length_factor=2.0
MultiPhase.density_fluid1=1000.
MultiPhase.density_fluid2=1.
MultiPhase.interface_smoothing=0
MultiPhase.interface_smoothing_frequency=1
ICNS.source_terms = GravityForcing
ICNS.use_perturb_pressure = true
MultiPhase.verbose=1
Expand All @@ -58,15 +53,10 @@ geometry.prob_lo = 0.0 0. -1 # Lo corner coordinates
geometry.prob_hi = 24. 1. 1 # Hi corner coordinates
geometry.is_periodic = 0 1 0 # Periodicity x y z (0/1)

xlo.type = "pressure_inflow"
xlo.type = "wave_generation"
xhi.type = "pressure_outflow"
xlo.vof_type = "zero_gradient"
xhi.vof_type = "zero_gradient"

zlo.type = "slip_wall"
zhi.type = "pressure_outflow"
zlo.vof_type = "zero_gradient"
zhi.vof_type = "zero_gradient"

incflo.verbose=1
zhi.type = "slip_wall"

incflo.verbose=1
4 changes: 1 addition & 3 deletions test/test_files/ow_w2a/ow_w2a.inp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@ time.cfl = 0.95 # CFL factor
#.......................................#
time.plot_interval = 10 # Steps between plot files
time.checkpoint_interval = -1 # Steps between checkpoint files
io.outputs = density velocity p vof ow_levelset ow_velocity w2a_levelset w2a_velocity
io.outputs = density velocity p vof ow_levelset ow_velocity

OceanWaves.label = W2A1
OceanWaves.W2A1.type = W2AWaves
Expand All @@ -29,8 +29,6 @@ OceanWaves.W2A1.number_interp_above_surface = 5
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# PHYSICS #
#.......................................#
incflo.use_godunov = 1
incflo.godunov_type = weno_z
incflo.mflux_type = minmod
transport.model = TwoPhaseTransport
transport.viscosity_fluid1=1.0e-3
Expand Down