Skip to content

Commit

Permalink
Enable multilevel inflow/outflow BC with native format (Exawind#970)
Browse files Browse the repository at this point in the history
  • Loading branch information
marchdf committed Feb 9, 2024
1 parent bf1c3cd commit db9add5
Show file tree
Hide file tree
Showing 7 changed files with 289 additions and 73 deletions.
1 change: 1 addition & 0 deletions amr-wind/wind_energy/ABLBoundaryPlane.H
Original file line number Diff line number Diff line change
Expand Up @@ -168,6 +168,7 @@ private:
const int,
const Field*);
#endif
int boundary_native_file_levels();

std::string m_title{"ABL boundary planes"};

Expand Down
177 changes: 104 additions & 73 deletions amr-wind/wind_energy/ABLBoundaryPlane.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -543,49 +543,52 @@ void ABLBoundaryPlane::write_file()
const std::string chkname =
m_filename + amrex::Concatenate("/bndry_output", t_step);

amrex::Print() << "Writing abl boundary checkpoint file " << chkname
amrex::Print() << "Writing ABL boundary checkpoint file " << chkname
<< " at time " << time << std::endl;

const int nlevels = m_repo.num_active_levels();
const std::string level_prefix = "Level_";
amrex::PreBuildDirectorHierarchy(chkname, level_prefix, 1, true);
amrex::PreBuildDirectorHierarchy(chkname, level_prefix, nlevels, true);
for (int lev = 0; lev < nlevels; ++lev) {
for (auto* fld : m_fields) {

// for now only output level 0
const int lev = 0;
auto& field = *fld;

for (auto* fld : m_fields) {
const auto& geom = field.repo().mesh().Geom();

auto& field = *fld;
// note: by using the entire domain box we end up using 1
// processor to hold all boundaries
amrex::Box domain = geom[lev].Domain();
amrex::BoxArray ba(domain);
amrex::DistributionMapping dm{ba};

const auto& geom = field.repo().mesh().Geom();
amrex::BndryRegister bndry(
ba, dm, m_in_rad, m_out_rad, m_extent_rad,
field.num_comp());

// note: by using the entire domain box we end up using 1 processor
// to hold all boundaries
amrex::Box domain = geom[lev].Domain();
amrex::BoxArray ba(domain);
amrex::DistributionMapping dm{ba};
bndry.setVal(1.0e13);

amrex::BndryRegister bndry(
ba, dm, m_in_rad, m_out_rad, m_extent_rad, field.num_comp());
bndry.copyFrom(
field(lev), 0, 0, 0, field.num_comp(),
geom[lev].periodicity());

bndry.copyFrom(
field(lev), 0, 0, 0, field.num_comp(), geom[lev].periodicity());
std::string filename = amrex::MultiFabFileFullPrefix(
lev, chkname, level_prefix, field.name());

std::string filename = amrex::MultiFabFileFullPrefix(
lev, chkname, level_prefix, field.name());
// print individual faces
for (amrex::OrientationIter oit; oit != nullptr; ++oit) {
auto ori = oit();
const std::string plane = m_plane_names[ori];

// print individual faces
for (amrex::OrientationIter oit; oit != nullptr; ++oit) {
auto ori = oit();
const std::string plane = m_plane_names[ori];
if (std::find(m_planes.begin(), m_planes.end(), plane) ==
m_planes.end()) {
continue;
}

if (std::find(m_planes.begin(), m_planes.end(), plane) ==
m_planes.end()) {
continue;
std::string facename =
amrex::Concatenate(filename + '_', ori, 1);
bndry[ori].write(facename);
}

std::string facename =
amrex::Concatenate(filename + '_', ori, 1);
bndry[ori].write(facename);
}
}
}
Expand Down Expand Up @@ -727,25 +730,31 @@ void ABLBoundaryPlane::read_header()
nc += fld->num_comp();
}

// TODO: need to generalize to lev > 0 somehow
const int lev = 0;
for (amrex::OrientationIter oit; oit != nullptr; ++oit) {
auto ori = oit();

// TODO: would be safer and less storage to not allocate all of
// these but we do not use m_planes for input and need to detect
// mass inflow from field bcs same for define level data below
if (std::all_of(
m_fields.begin(), m_fields.end(), [ori](const auto* fld) {
return fld->bc_type()[ori] != BC::mass_inflow;
})) {
continue;
}

m_in_data.define_plane(ori);

const amrex::Box& minBox = m_mesh.boxArray(lev).minimalBox();
const int nlevels = boundary_native_file_levels();
for (int lev = 0; lev < nlevels; ++lev) {

amrex::IntVect plo(minBox.loVect());
amrex::IntVect phi(minBox.hiVect());
const int normal = ori.coordDir();
plo[normal] = ori.isHigh() ? minBox.hiVect()[normal] + 1 : -1;
phi[normal] = ori.isHigh() ? minBox.hiVect()[normal] + 1 : -1;
const amrex::Box pbx(plo, phi);
m_in_data.define_level_data(ori, pbx, nc);
const amrex::Box& minBox = m_mesh.boxArray(lev).minimalBox();

amrex::IntVect plo(minBox.loVect());
amrex::IntVect phi(minBox.hiVect());
const int normal = ori.coordDir();
plo[normal] = ori.isHigh() ? minBox.hiVect()[normal] + 1 : -1;
phi[normal] = ori.isHigh() ? minBox.hiVect()[normal] + 1 : -1;
const amrex::Box pbx(plo, phi);
m_in_data.define_level_data(ori, pbx, nc);
}
}
}
}
Expand Down Expand Up @@ -807,47 +816,51 @@ void ABLBoundaryPlane::read_file()

const std::string level_prefix = "Level_";

const int lev = 0;
for (auto* fld : m_fields) {
const int nlevels = boundary_native_file_levels();
for (int lev = 0; lev < nlevels; ++lev) {
for (auto* fld : m_fields) {

auto& field = *fld;
const auto& geom = field.repo().mesh().Geom();
auto& field = *fld;
const auto& geom = field.repo().mesh().Geom();

amrex::Box domain = geom[lev].Domain();
amrex::BoxArray ba(domain);
amrex::DistributionMapping dm{ba};
amrex::Box domain = geom[lev].Domain();
amrex::BoxArray ba(domain);
amrex::DistributionMapping dm{ba};

amrex::BndryRegister bndry1(
ba, dm, m_in_rad, m_out_rad, m_extent_rad, field.num_comp());
amrex::BndryRegister bndry2(
ba, dm, m_in_rad, m_out_rad, m_extent_rad, field.num_comp());
amrex::BndryRegister bndry1(
ba, dm, m_in_rad, m_out_rad, m_extent_rad,
field.num_comp());
amrex::BndryRegister bndry2(
ba, dm, m_in_rad, m_out_rad, m_extent_rad,
field.num_comp());

bndry1.setVal(1.0e13);
bndry2.setVal(1.0e13);
bndry1.setVal(1.0e13);
bndry2.setVal(1.0e13);

std::string filename1 = amrex::MultiFabFileFullPrefix(
lev, chkname1, level_prefix, field.name());
std::string filename2 = amrex::MultiFabFileFullPrefix(
lev, chkname2, level_prefix, field.name());
std::string filename1 = amrex::MultiFabFileFullPrefix(
lev, chkname1, level_prefix, field.name());
std::string filename2 = amrex::MultiFabFileFullPrefix(
lev, chkname2, level_prefix, field.name());

for (amrex::OrientationIter oit; oit != nullptr; ++oit) {
auto ori = oit();
for (amrex::OrientationIter oit; oit != nullptr; ++oit) {
auto ori = oit();

if ((!m_in_data.is_populated(ori)) ||
(field.bc_type()[ori] != BC::mass_inflow)) {
continue;
}
if ((!m_in_data.is_populated(ori)) ||
(field.bc_type()[ori] != BC::mass_inflow)) {
continue;
}

std::string facename1 =
amrex::Concatenate(filename1 + '_', ori, 1);
std::string facename2 =
amrex::Concatenate(filename2 + '_', ori, 1);
std::string facename1 =
amrex::Concatenate(filename1 + '_', ori, 1);
std::string facename2 =
amrex::Concatenate(filename2 + '_', ori, 1);

bndry1[ori].read(facename1);
bndry2[ori].read(facename2);
bndry1[ori].read(facename1);
bndry2[ori].read(facename2);

m_in_data.read_data_native(
oit, bndry1, bndry2, lev, fld, time, m_in_times);
m_in_data.read_data_native(
oit, bndry1, bndry2, lev, fld, time, m_in_times);
}
}
}
}
Expand Down Expand Up @@ -1052,6 +1065,24 @@ void ABLBoundaryPlane::impl_buffer_field(
}
#endif

// Count the number of levels defined by the native boundary files
int ABLBoundaryPlane::boundary_native_file_levels()
{
int nlevels = 0;
const std::string chkname =
m_filename + amrex::Concatenate("/bndry_output", m_in_timesteps[0]);
const std::string level_prefix = "Level_";
for (int lev = 0; lev < m_repo.num_active_levels(); ++lev) {
const std::string levname = amrex::LevelFullPath(lev, chkname);
if (amrex::FileExists(levname)) {
nlevels = lev + 1;
} else {
break;
}
}
return nlevels;
}

//! True if box intersects the boundary
bool ABLBoundaryPlane::box_intersects_boundary(
const amrex::Box& bx, const int lev, const amrex::Orientation ori) const
Expand Down
2 changes: 2 additions & 0 deletions test/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -217,6 +217,7 @@ add_test_re(ib_sphere_Re_100)
add_test_re(vortex_ring_collision)
add_test_re(fat_cored_vortex_ring)
add_test_re(abl_bndry_output_native)
add_test_re(abl_bndry_output_amr_native)
add_test_re(vortex_patch_scalar_vel)
add_test_re(zalesak_disk_scalar_vel)
add_test_re(rain_drop)
Expand Down Expand Up @@ -284,6 +285,7 @@ endif()
add_test_red(abl_bndry_input_native abl_bndry_output_native)
add_test_red(abl_godunov_restart abl_godunov)
add_test_red(abl_bndry_input_amr_native abl_bndry_output_native)
add_test_red(abl_bndry_input_amr_native_mlbc abl_bndry_output_amr_native)

if(AMR_WIND_ENABLE_NETCDF)
add_test_red(abl_bndry_input abl_bndry_output)
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,88 @@
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# SIMULATION STOP #
#.......................................#
time.stop_time = 22000.0 # Max (simulated) time to evolve
time.max_step = 10 # Max number of time steps
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# TIME STEP COMPUTATION #
#.......................................#
time.fixed_dt = 0.4 # Use this constant dt if > 0
time.cfl = 0.95 # CFL factor
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# INPUT AND OUTPUT #
#.......................................#
io.restart_file = "../abl_bndry_output_amr_native/chk00005"
time.plot_interval = 10 # Steps between plot files
time.checkpoint_interval = -1 # Steps between checkpoint files
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# PHYSICS #
#.......................................#
incflo.gravity = 0. 0. -9.81 # Gravitational force (3D)
incflo.density = 1.0 # Reference density
incflo.use_godunov = 1
incflo.diffusion_type = 2
transport.viscosity = 1.0e-5
transport.laminar_prandtl = 0.7
transport.turbulent_prandtl = 0.3333
turbulence.model = Smagorinsky
Smagorinsky_coeffs.Cs = 0.135
incflo.physics = ABL
ICNS.source_terms = CoriolisForcing GeostrophicForcing
BoussinesqBuoyancy.reference_temperature = 290.0
ABL.reference_temperature = 290.0
CoriolisForcing.east_vector = 1.0 0.0 0.0
CoriolisForcing.north_vector = 0.0 1.0 0.0
CoriolisForcing.latitude = 90.0
CoriolisForcing.rotational_time_period = 125663.706143592
GeostrophicForcing.geostrophic_wind = 10.0 0.0 0.0
incflo.velocity = 10.0 0.0 0.0
ABL.temperature_heights = 0.0 2000.0
ABL.temperature_values = 290.0 290.0
ABL.perturb_temperature = false
ABL.cutoff_height = 50.0
ABL.perturb_velocity = true
ABL.perturb_ref_height = 50.0
ABL.Uperiods = 4.0
ABL.Vperiods = 4.0
ABL.deltaU = 1.0
ABL.deltaV = 1.0
ABL.kappa = .41
ABL.surface_roughness_z0 = 0.01
ABL.bndry_file = "../abl_bndry_output_amr_native/bndry_files"
ABL.bndry_io_mode = 1
ABL.bndry_var_names = velocity temperature
ABL.bndry_output_format = native
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# ADAPTIVE MESH REFINEMENT #
#.......................................#
amr.n_cell = 48 48 48 # Grid cells at coarsest AMRlevel
amr.max_level = 1 # Max AMR level in hierarchy

tagging.labels = static
tagging.static.type = "CartBoxRefinement"
tagging.static.static_refinement_def = "static_box.txt"

#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# GEOMETRY #
#.......................................#
geometry.prob_lo = 0. 0. 0. # Lo corner coordinates
geometry.prob_hi = 1000. 1000. 1000. # Hi corner coordinates
geometry.is_periodic = 0 0 0 # Periodicity x y z (0/1)
incflo.delp = 0. 0. 0. # Prescribed (cyclic) pressure gradient
# Boundary conditions
xlo.type = "mass_inflow"
xlo.density = 1.0
xlo.temperature = 0.0
xhi.type = "pressure_outflow"
ylo.type = "mass_inflow"
ylo.density = 1.0
ylo.temperature = 0.0
yhi.type = "pressure_outflow"
zlo.type = "wall_model"
zhi.type = "slip_wall"
zhi.temperature_type = "fixed_gradient"
zhi.temperature = 0.0
#¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨¨#
# VERBOSITY #
#.......................................#
incflo.verbose = 0 # incflo_level
Original file line number Diff line number Diff line change
@@ -0,0 +1,5 @@
2
1
300.0 300.0 40.0 700.0 700.0 400.0
1
450.0 450.0 100.0 550.0 550.0 200.0
Loading

0 comments on commit db9add5

Please sign in to comment.