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.

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

Ensure consistent ghost cells for terrain fields #1353

Merged
merged 10 commits into from
Nov 18, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
9 changes: 3 additions & 6 deletions amr-wind/physics/TerrainDrag.H
Original file line number Diff line number Diff line change
Expand Up @@ -22,13 +22,11 @@ public:

~TerrainDrag() override = default;

void
initialize_fields(int /*level*/, const amrex::Geometry& /*geom*/) override
{}
void initialize_fields(int level, const amrex::Geometry& geom) override;

void pre_init_actions() override;
void pre_init_actions() override {}

void post_init_actions() override;
void post_init_actions() override {}

void post_regrid_actions() override {}

Expand All @@ -42,7 +40,6 @@ private:
CFDSim& m_sim;
const FieldRepo& m_repo;
const amrex::AmrCore& m_mesh;
Field& m_velocity;
//! Blanking Field for Terrain or Buildings
IntField& m_terrain_blank;
//! Terrain Drag Force Term
Expand Down
178 changes: 83 additions & 95 deletions amr-wind/physics/TerrainDrag.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,6 @@ TerrainDrag::TerrainDrag(CFDSim& sim)
: m_sim(sim)
, m_repo(sim.repo())
, m_mesh(sim.mesh())
, m_velocity(sim.repo().get_field("velocity"))
, m_terrain_blank(sim.repo().declare_int_field("terrain_blank", 1, 1, 1))
, m_terrain_drag(sim.repo().declare_int_field("terrain_drag", 1, 1, 1))
, m_terrainz0(sim.repo().declare_field("terrainz0", 1, 1, 1))
Expand All @@ -44,104 +43,93 @@ TerrainDrag::TerrainDrag(CFDSim& sim)
m_sim.io_manager().register_output_int_var("terrain_blank");
m_sim.io_manager().register_io_var("terrainz0");
m_sim.io_manager().register_io_var("terrain_height");

m_terrain_blank.setVal(0.0);
m_terrain_drag.setVal(0.0);
m_terrainz0.set_default_fillpatch_bc(m_sim.time());
m_terrain_height.set_default_fillpatch_bc(m_sim.time());
}

void TerrainDrag::post_init_actions()
void TerrainDrag::initialize_fields(int level, const amrex::Geometry& geom)
{
BL_PROFILE("amr-wind::" + this->identifier() + "::post_init_actions");
const auto& geom_vec = m_sim.repo().mesh().Geom();
const int nlevels = m_sim.repo().num_active_levels();
for (int level = 0; level < nlevels; ++level) {
const auto& geom = geom_vec[level];
const auto& dx = geom.CellSizeArray();
const auto& prob_lo = geom.ProbLoArray();
auto& velocity = m_velocity(level);
auto& blanking = m_terrain_blank(level);
auto& terrainz0 = m_terrainz0(level);
auto& terrain_height = m_terrain_height(level);
auto& drag = m_terrain_drag(level);
// copy terrain data to gpu
const auto xterrain_size = m_xterrain.size();
const auto yterrain_size = m_yterrain.size();
const auto zterrain_size = m_zterrain.size();
amrex::Gpu::DeviceVector<amrex::Real> device_xterrain(xterrain_size);
amrex::Gpu::DeviceVector<amrex::Real> device_yterrain(yterrain_size);
amrex::Gpu::DeviceVector<amrex::Real> device_zterrain(zterrain_size);
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_xterrain.begin(), m_xterrain.end(),
device_xterrain.begin());
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_yterrain.begin(), m_yterrain.end(),
device_yterrain.begin());
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_zterrain.begin(), m_zterrain.end(),
device_zterrain.begin());
const auto* xterrain_ptr = device_xterrain.data();
const auto* yterrain_ptr = device_yterrain.data();
const auto* zterrain_ptr = device_zterrain.data();
// Copy Roughness to gpu
const auto xrough_size = m_xrough.size();
const auto yrough_size = m_yrough.size();
const auto z0rough_size = m_z0rough.size();
amrex::Gpu::DeviceVector<amrex::Real> device_xrough(xrough_size);
amrex::Gpu::DeviceVector<amrex::Real> device_yrough(yrough_size);
amrex::Gpu::DeviceVector<amrex::Real> device_z0rough(z0rough_size);
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_xrough.begin(), m_xrough.end(),
device_xrough.begin());
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_yrough.begin(), m_yrough.end(),
device_yrough.begin());
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_z0rough.begin(), m_z0rough.end(),
device_z0rough.begin());
const auto* xrough_ptr = device_xrough.data();
const auto* yrough_ptr = device_yrough.data();
const auto* z0rough_ptr = device_z0rough.data();
for (amrex::MFIter mfi(velocity); mfi.isValid(); ++mfi) {
const auto& vbx = mfi.validbox();
auto levelBlanking = blanking.array(mfi);
auto levelDrag = drag.array(mfi);
auto levelz0 = terrainz0.array(mfi);
auto levelheight = terrain_height.array(mfi);
amrex::ParallelFor(
vbx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
const amrex::Real x = prob_lo[0] + (i + 0.5) * dx[0];
const amrex::Real y = prob_lo[1] + (j + 0.5) * dx[1];
const amrex::Real z = prob_lo[2] + (k + 0.5) * dx[2];
const amrex::Real terrainHt = interp::bilinear(
xterrain_ptr, xterrain_ptr + xterrain_size,
yterrain_ptr, yterrain_ptr + yterrain_size,
zterrain_ptr, x, y);
levelBlanking(i, j, k, 0) =
static_cast<int>(z <= terrainHt);
levelheight(i, j, k, 0) =
std::max(std::abs(z - terrainHt), 0.5 * dx[2]);

amrex::Real roughz0 = 0.1;
if (xrough_size > 0) {
roughz0 = interp::bilinear(
xrough_ptr, xrough_ptr + xrough_size, yrough_ptr,
yrough_ptr + yrough_size, z0rough_ptr, x, y);
}
levelz0(i, j, k, 0) = roughz0;
});
BL_PROFILE("amr-wind::" + this->identifier() + "::initialize_fields");
const auto& dx = geom.CellSizeArray();
const auto& prob_lo = geom.ProbLoArray();
auto& blanking = m_terrain_blank(level);
auto& terrainz0 = m_terrainz0(level);
auto& terrain_height = m_terrain_height(level);
auto& drag = m_terrain_drag(level);
const auto xterrain_size = m_xterrain.size();
const auto yterrain_size = m_yterrain.size();
const auto zterrain_size = m_zterrain.size();
amrex::Gpu::DeviceVector<amrex::Real> device_xterrain(xterrain_size);
amrex::Gpu::DeviceVector<amrex::Real> device_yterrain(yterrain_size);
amrex::Gpu::DeviceVector<amrex::Real> device_zterrain(zterrain_size);
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_xterrain.begin(), m_xterrain.end(),
device_xterrain.begin());
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_yterrain.begin(), m_yterrain.end(),
device_yterrain.begin());
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_zterrain.begin(), m_zterrain.end(),
device_zterrain.begin());
const auto* xterrain_ptr = device_xterrain.data();
const auto* yterrain_ptr = device_yterrain.data();
const auto* zterrain_ptr = device_zterrain.data();
// Copy Roughness to gpu
const auto xrough_size = m_xrough.size();
const auto yrough_size = m_yrough.size();
const auto z0rough_size = m_z0rough.size();
amrex::Gpu::DeviceVector<amrex::Real> device_xrough(xrough_size);
amrex::Gpu::DeviceVector<amrex::Real> device_yrough(yrough_size);
amrex::Gpu::DeviceVector<amrex::Real> device_z0rough(z0rough_size);
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_xrough.begin(), m_xrough.end(),
device_xrough.begin());
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_yrough.begin(), m_yrough.end(),
device_yrough.begin());
amrex::Gpu::copy(
amrex::Gpu::hostToDevice, m_z0rough.begin(), m_z0rough.end(),
device_z0rough.begin());
const auto* xrough_ptr = device_xrough.data();
const auto* yrough_ptr = device_yrough.data();
const auto* z0rough_ptr = device_z0rough.data();
auto levelBlanking = blanking.arrays();
auto levelDrag = drag.arrays();
auto levelz0 = terrainz0.arrays();
auto levelheight = terrain_height.arrays();
amrex::ParallelFor(
blanking, m_terrain_blank.num_grow(),
[=] AMREX_GPU_DEVICE(int nbx, int i, int j, int k) noexcept {
const amrex::Real x = prob_lo[0] + (i + 0.5) * dx[0];
const amrex::Real y = prob_lo[1] + (j + 0.5) * dx[1];
const amrex::Real z = prob_lo[2] + (k + 0.5) * dx[2];
const amrex::Real terrainHt = interp::bilinear(
xterrain_ptr, xterrain_ptr + xterrain_size, yterrain_ptr,
yterrain_ptr + yterrain_size, zterrain_ptr, x, y);
levelBlanking[nbx](i, j, k, 0) = static_cast<int>(z <= terrainHt);
levelheight[nbx](i, j, k, 0) =
std::max(std::abs(z - terrainHt), 0.5 * dx[2]);

amrex::ParallelFor(
vbx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
if ((levelBlanking(i, j, k, 0) == 0) && (k > 0) &&
(levelBlanking(i, j, k - 1, 0) == 1)) {
levelDrag(i, j, k, 0) = 1;
} else {
levelDrag(i, j, k, 0) = 0;
}
});
}
}
amrex::Real roughz0 = 0.1;
if (xrough_size > 0) {
roughz0 = interp::bilinear(
xrough_ptr, xrough_ptr + xrough_size, yrough_ptr,
yrough_ptr + yrough_size, z0rough_ptr, x, y);
}
levelz0[nbx](i, j, k, 0) = roughz0;
});
amrex::ParallelFor(
blanking, [=] AMREX_GPU_DEVICE(int nbx, int i, int j, int k) noexcept {
if ((levelBlanking[nbx](i, j, k, 0) == 0) && (k > 0) &&
(levelBlanking[nbx](i, j, k - 1, 0) == 1)) {
levelDrag[nbx](i, j, k, 0) = 1;
} else {
levelDrag[nbx](i, j, k, 0) = 0;
}
});
}

void TerrainDrag::pre_init_actions()
{
BL_PROFILE("amr-wind::" + this->identifier() + "::pre_init_actions");
}
} // namespace amr_wind::terraindrag
1 change: 1 addition & 0 deletions amr-wind/utilities/sampling/SamplingContainer.H
Original file line number Diff line number Diff line change
Expand Up @@ -95,6 +95,7 @@ public:
for (ParIterType pti(*this, lev); pti.isValid(); ++pti) {
int scomp_curr = scomp;
for (const auto* fld : fields) {
AMREX_ALWAYS_ASSERT(fld->num_grow() > amrex::IntVect{0});
const auto farr = (*fld)(lev).const_array(pti);
interpolate(
pti, farr, lev, fld->field_location(), fld->num_comp(),
Expand Down
2 changes: 1 addition & 1 deletion test/test_files/terrain_box/terrain_box.inp
Original file line number Diff line number Diff line change
Expand Up @@ -56,7 +56,7 @@ incflo.verbose = 0 # incflo_level
incflo.post_processing = sampling
sampling.labels = volume1
sampling.fields = velocity
sampling.int_fields = terrain_blank
sampling.int_fields = terrain_blank terrain_drag
sampling.volume1.type = VolumeSampler
sampling.volume1.hi = 600 600 200
sampling.volume1.lo = 400 400 0.0
Expand Down
6 changes: 5 additions & 1 deletion unit_tests/wind_energy/test_abl_terrain.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -74,7 +74,11 @@ TEST_F(TerrainTest, terrain)
amrex::Vector<std::string> physics{"terrainDrag"};
pp.addarr("physics", physics);
amr_wind::terraindrag::TerrainDrag terrain_drag(sim());
terrain_drag.post_init_actions();
const int nlevels = sim().repo().num_active_levels();
for (int lev = 0; lev < nlevels; ++lev) {
const auto& geom = sim().repo().mesh().Geom(lev);
terrain_drag.initialize_fields(lev, geom);
}
const auto& terrain_blank = sim().repo().get_int_field("terrain_blank");
// Outside Point
const int value_out = utils::field_probe(terrain_blank, 0, 5, 5, 1);
Expand Down