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

Inflow face from file data #1698

Merged
merged 5 commits into from
Jul 19, 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
10 changes: 6 additions & 4 deletions Source/BoundaryConditions/BoundaryConditions_xvel.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -52,6 +52,8 @@ void ERFPhysBCFunct_u::impose_lateral_xvel_bcs (const Array4<Real>& dest_arr,
// First do all ext_dir bcs
if (!is_periodic_in_x)
{
Real* xvel_xlo_ptr = m_u_bc_data[0];
Real* xvel_xhi_ptr = m_u_bc_data[3];
Box bx_xlo(bx); bx_xlo.setBig (0,dom_lo.x-1);
Box bx_xhi(bx); bx_xhi.setSmall(0,dom_hi.x+2);
Box bx_xlo_face(bx); bx_xlo_face.setSmall(0,dom_lo.x ); bx_xlo_face.setBig(0,dom_lo.x );
Expand All @@ -61,7 +63,7 @@ void ERFPhysBCFunct_u::impose_lateral_xvel_bcs (const Array4<Real>& dest_arr,
{
int iflip = dom_lo.x - i;
if (bc_ptr[n].lo(0) == ERFBCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][0];
dest_arr(i,j,k) = (xvel_xlo_ptr) ? xvel_xlo_ptr[k] : l_bc_extdir_vals_d[n][0];
} else if (bc_ptr[n].lo(0) == ERFBCType::foextrap) {
dest_arr(i,j,k) = dest_arr(dom_lo.x,j,k);
} else if (bc_ptr[n].lo(0) == ERFBCType::open) {
Expand All @@ -78,7 +80,7 @@ void ERFPhysBCFunct_u::impose_lateral_xvel_bcs (const Array4<Real>& dest_arr,
bx_xlo_face, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n)
{
if (bc_ptr[n].lo(0) == ERFBCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][0];
dest_arr(i,j,k) = (xvel_xlo_ptr) ? xvel_xlo_ptr[k] : l_bc_extdir_vals_d[n][0];
} else if (bc_ptr[n].lo(0) == ERFBCType::neumann_int) {
dest_arr(i,j,k) = (4.0*dest_arr(dom_lo.x+1,j,k) - dest_arr(dom_lo.x+2,j,k))/3.0;
}
Expand All @@ -89,7 +91,7 @@ void ERFPhysBCFunct_u::impose_lateral_xvel_bcs (const Array4<Real>& dest_arr,
{
int iflip = 2*(dom_hi.x + 1) - i;
if (bc_ptr[n].hi(0) == ERFBCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][3];
dest_arr(i,j,k) = (xvel_xhi_ptr) ? xvel_xhi_ptr[k] : l_bc_extdir_vals_d[n][3];
} else if (bc_ptr[n].hi(0) == ERFBCType::foextrap) {
dest_arr(i,j,k) = dest_arr(dom_hi.x+1,j,k);
} else if (bc_ptr[n].hi(0) == ERFBCType::open) {
Expand All @@ -106,7 +108,7 @@ void ERFPhysBCFunct_u::impose_lateral_xvel_bcs (const Array4<Real>& dest_arr,
bx_xhi_face, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n)
{
if (bc_ptr[n].hi(0) == ERFBCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][3];
dest_arr(i,j,k) = (xvel_xhi_ptr) ? xvel_xhi_ptr[k] : l_bc_extdir_vals_d[n][3];
} else if (bc_ptr[n].hi(0) == ERFBCType::neumann_int) {
dest_arr(i,j,k) = (4.0*dest_arr(dom_hi.x,j,k) - dest_arr(dom_hi.x-1,j,k))/3.0;
}
Expand Down
6 changes: 4 additions & 2 deletions Source/BoundaryConditions/BoundaryConditions_yvel.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -51,13 +51,15 @@ void ERFPhysBCFunct_v::impose_lateral_yvel_bcs (const Array4<Real>& dest_arr,
if (!is_periodic_in_x)
{
// Populate ghost cells on lo-x and hi-x domain boundaries
Real* yvel_xlo_ptr = m_v_bc_data[0];
Real* yvel_xhi_ptr = m_v_bc_data[3];
Box bx_xlo(bx); bx_xlo.setBig (0,dom_lo.x-1);
Box bx_xhi(bx); bx_xhi.setSmall(0,dom_hi.x+1);
ParallelFor(
bx_xlo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int iflip = dom_lo.x - 1- i;
if (bc_ptr[n].lo(0) == ERFBCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][0];
dest_arr(i,j,k) = (yvel_xlo_ptr) ? yvel_xlo_ptr[k] : l_bc_extdir_vals_d[n][0];
} else if (bc_ptr[n].lo(0) == ERFBCType::foextrap) {
dest_arr(i,j,k) = dest_arr(dom_lo.x,j,k);
} else if (bc_ptr[n].lo(0) == ERFBCType::open) {
Expand All @@ -71,7 +73,7 @@ void ERFPhysBCFunct_v::impose_lateral_yvel_bcs (const Array4<Real>& dest_arr,
bx_xhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int iflip = 2*dom_hi.x + 1 - i;
if (bc_ptr[n].hi(0) == ERFBCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][3];
dest_arr(i,j,k) = (yvel_xhi_ptr) ? yvel_xhi_ptr[k] : l_bc_extdir_vals_d[n][3];
} else if (bc_ptr[n].hi(0) == ERFBCType::foextrap) {
dest_arr(i,j,k) = dest_arr(dom_hi.x,j,k);
} else if (bc_ptr[n].hi(0) == ERFBCType::open) {
Expand Down
6 changes: 4 additions & 2 deletions Source/BoundaryConditions/BoundaryConditions_zvel.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -60,13 +60,15 @@ void ERFPhysBCFunct_w::impose_lateral_zvel_bcs (const Array4<Real >& dest_a
// First do all ext_dir bcs
if (!is_periodic_in_x)
{
Real* zvel_xlo_ptr = m_w_bc_data[0];
Real* zvel_xhi_ptr = m_w_bc_data[3];
Box bx_xlo(bx); bx_xlo.setBig (0,dom_lo.x-1);
Box bx_xhi(bx); bx_xhi.setSmall(0,dom_hi.x+1);
ParallelFor(
bx_xlo, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int iflip = dom_lo.x - 1 - i;
if (bc_ptr_w[n].lo(0) == ERFBCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][0];
dest_arr(i,j,k) = (zvel_xlo_ptr) ? zvel_xlo_ptr[k] : l_bc_extdir_vals_d[n][0];
if (l_use_terrain) {
dest_arr(i,j,k) = WFromOmega(i,j,k,dest_arr(i,j,k),xvel_arr,yvel_arr,z_phys_nd,dxInv);
}
Expand All @@ -83,7 +85,7 @@ void ERFPhysBCFunct_w::impose_lateral_zvel_bcs (const Array4<Real >& dest_a
bx_xhi, ncomp, [=] AMREX_GPU_DEVICE (int i, int j, int k, int n) {
int iflip = 2*dom_hi.x + 1 - i;
if (bc_ptr_w[n].hi(0) == ERFBCType::ext_dir) {
dest_arr(i,j,k) = l_bc_extdir_vals_d[n][3];
dest_arr(i,j,k) = (zvel_xhi_ptr) ? zvel_xhi_ptr[k] : l_bc_extdir_vals_d[n][3];
if (l_use_terrain) {
dest_arr(i,j,k) = WFromOmega(i,j,k,dest_arr(i,j,k),xvel_arr,yvel_arr,z_phys_nd,dxInv);
}
Expand Down
8 changes: 4 additions & 4 deletions Source/BoundaryConditions/ERF_FillPatch.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -108,7 +108,7 @@ ERF::FillPatch (int lev, Real time,
FillPatchSingleLevel(*mfs_vel[Vars::zvel], time, fmf, ftime, icomp, icomp, ncomp_w,
geom[lev], *physbcs_w_no_terrain[lev], BCVars::zvel_bc);
(*physbcs_w[lev])(*mfs_vel[Vars::zvel],*mfs_vel[Vars::xvel],*mfs_vel[Vars::yvel],
ngvect_vels,time,BCVars::zvel_bc);
ngvect_vels,time,BCVars::zvel_bc);
} // !cons_only

} else {
Expand Down Expand Up @@ -161,7 +161,7 @@ ERF::FillPatch (int lev, Real time,
*physbcs_w_no_terrain[lev ], BCVars::zvel_bc,
refRatio(lev-1), mapper, domain_bcs_type, BCVars::zvel_bc);
(*physbcs_w[lev])(*mfs_vel[Vars::zvel],*mfs_vel[Vars::xvel],*mfs_vel[Vars::yvel],
ngvect_vels,time,BCVars::zvel_bc);
ngvect_vels,time,BCVars::zvel_bc);
} // !cons_only
} // lev > 0

Expand All @@ -186,7 +186,7 @@ ERF::FillPatch (int lev, Real time,
(*physbcs_u[lev])(*mfs_vel[Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc);
(*physbcs_v[lev])(*mfs_vel[Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc);
(*physbcs_w[lev])(*mfs_vel[Vars::zvel],*mfs_vel[Vars::xvel],*mfs_vel[Vars::yvel],
ngvect_vels,time,BCVars::zvel_bc);
ngvect_vels,time,BCVars::zvel_bc);
}
}

Expand Down Expand Up @@ -410,7 +410,7 @@ ERF::FillIntermediatePatch (int lev, Real time,
(*physbcs_u[lev])(*mfs_vel[Vars::xvel],0,1,ngvect_vels,time,BCVars::xvel_bc);
(*physbcs_v[lev])(*mfs_vel[Vars::yvel],0,1,ngvect_vels,time,BCVars::yvel_bc);
(*physbcs_w[lev])(*mfs_vel[Vars::zvel],*mfs_vel[Vars::xvel],*mfs_vel[Vars::yvel],
ngvect_vels,time,BCVars::zvel_bc);
ngvect_vels,time,BCVars::zvel_bc);
}
// ***************************************************************************

Expand Down
48 changes: 40 additions & 8 deletions Source/BoundaryConditions/ERF_PhysBCFunct.H
Original file line number Diff line number Diff line change
Expand Up @@ -83,15 +83,22 @@ public:
amrex::Array<amrex::Array<amrex::Real,AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> bc_extdir_vals,
amrex::Array<amrex::Array<amrex::Real,AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> bc_neumann_vals,
std::unique_ptr<amrex::MultiFab>& z_phys_nd,
const bool use_real_bcs)
const bool use_real_bcs,
amrex::Vector<amrex::Real*> u_bc_data)
: m_lev(lev), m_geom(geom),
m_domain_bcs_type(domain_bcs_type),
m_domain_bcs_type_d(domain_bcs_type_d),
m_bc_extdir_vals(bc_extdir_vals),
m_bc_neumann_vals(bc_neumann_vals),
m_z_phys_nd(z_phys_nd.get()),
m_use_real_bcs(use_real_bcs)
{}
{
int nval = u_bc_data.size();
m_u_bc_data.resize(nval);
for (int ival(0); ival<nval; ++ival) {
m_u_bc_data[ival] = u_bc_data[ival];
}
}

~ERFPhysBCFunct_u () {}

Expand Down Expand Up @@ -125,6 +132,7 @@ private:
amrex::Array<amrex::Array<amrex::Real, AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> m_bc_neumann_vals;
amrex::MultiFab* m_z_phys_nd;
bool m_use_real_bcs;
amrex::Vector<amrex::Real*> m_u_bc_data;
};

class ERFPhysBCFunct_v
Expand All @@ -136,15 +144,22 @@ public:
amrex::Array<amrex::Array<amrex::Real,AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> bc_extdir_vals,
amrex::Array<amrex::Array<amrex::Real,AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> bc_neumann_vals,
std::unique_ptr<amrex::MultiFab>& z_phys_nd,
const bool use_real_bcs)
const bool use_real_bcs,
amrex::Vector<amrex::Real*> v_bc_data)
: m_lev(lev),
m_geom(geom), m_domain_bcs_type(domain_bcs_type),
m_domain_bcs_type_d(domain_bcs_type_d),
m_bc_extdir_vals(bc_extdir_vals),
m_bc_neumann_vals(bc_neumann_vals),
m_z_phys_nd(z_phys_nd.get()),
m_use_real_bcs(use_real_bcs)
{}
{
int nval = v_bc_data.size();
m_v_bc_data.resize(nval);
for (int ival(0); ival<nval; ++ival) {
m_v_bc_data[ival] = v_bc_data[ival];
}
}

~ERFPhysBCFunct_v () {}

Expand Down Expand Up @@ -179,6 +194,7 @@ private:
amrex::Array<amrex::Array<amrex::Real, AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> m_bc_neumann_vals;
amrex::MultiFab* m_z_phys_nd;
bool m_use_real_bcs;
amrex::Vector<amrex::Real*> m_v_bc_data;
};

class ERFPhysBCFunct_w
Expand All @@ -190,7 +206,8 @@ public:
amrex::Array<amrex::Array<amrex::Real,AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> bc_extdir_vals,
amrex::Array<amrex::Array<amrex::Real,AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> bc_neumann_vals,
const TerrainType& terrain_type, std::unique_ptr<amrex::MultiFab>& z_phys_nd,
const bool use_real_bcs)
const bool use_real_bcs,
amrex::Vector<amrex::Real*> w_bc_data)
: m_lev(lev),
m_geom(geom), m_domain_bcs_type(domain_bcs_type),
m_domain_bcs_type_d(domain_bcs_type_d),
Expand All @@ -199,7 +216,13 @@ public:
m_terrain_type(terrain_type),
m_z_phys_nd(z_phys_nd.get()),
m_use_real_bcs(use_real_bcs)
{}
{
int nval = w_bc_data.size();
m_w_bc_data.resize(nval);
for (int ival(0); ival<nval; ++ival) {
m_w_bc_data[ival] = w_bc_data[ival];
}
}

~ERFPhysBCFunct_w () {}

Expand Down Expand Up @@ -242,6 +265,7 @@ private:
TerrainType m_terrain_type;
amrex::MultiFab* m_z_phys_nd;
bool m_use_real_bcs;
amrex::Vector<amrex::Real*> m_w_bc_data;
};

class ERFPhysBCFunct_w_no_terrain
Expand All @@ -253,14 +277,21 @@ public:
const amrex::Gpu::DeviceVector<amrex::BCRec>& domain_bcs_type_d,
amrex::Array<amrex::Array<amrex::Real,AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> bc_extdir_vals,
amrex::Array<amrex::Array<amrex::Real,AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> bc_neumann_vals,
const bool use_real_bcs)
const bool use_real_bcs,
amrex::Vector<amrex::Real*> w_bc_data)
: m_lev(lev),
m_geom(geom), m_domain_bcs_type(domain_bcs_type),
m_domain_bcs_type_d(domain_bcs_type_d),
m_bc_extdir_vals(bc_extdir_vals),
m_bc_neumann_vals(bc_neumann_vals),
m_use_real_bcs(use_real_bcs)
{}
{
int nval = w_bc_data.size();
m_w_bc_data.resize(nval);
for (int ival(0); ival<nval; ++ival) {
m_w_bc_data[ival] = w_bc_data[ival];
}
}

~ERFPhysBCFunct_w_no_terrain () {}

Expand Down Expand Up @@ -292,6 +323,7 @@ private:
amrex::Array<amrex::Array<amrex::Real, AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> m_bc_extdir_vals;
amrex::Array<amrex::Array<amrex::Real, AMREX_SPACEDIM*2>,AMREX_SPACEDIM+NVAR_max> m_bc_neumann_vals;
bool m_use_real_bcs;
amrex::Vector<amrex::Real*> m_w_bc_data;
};

#endif
9 changes: 9 additions & 0 deletions Source/ERF.H
Original file line number Diff line number Diff line change
Expand Up @@ -580,6 +580,15 @@ private:
// Struct for working with the sponge data we take as an input
InputSpongeData input_sponge_data;

// Vector (6 planes) of DeviceVectors (ncell in plane) for Dirichlet BC data
amrex::Vector<amrex::Gpu::DeviceVector<amrex::Real>> xvel_bc_data;
amrex::Vector<amrex::Gpu::DeviceVector<amrex::Real>> yvel_bc_data;
amrex::Vector<amrex::Gpu::DeviceVector<amrex::Real>> zvel_bc_data;

// Function to read and populate above vectors (if input file exists)
void init_Dirichlet_bc_data (amrex::Orientation ori,
const std::string input_file);

// write checkpoint file to disk
void WriteCheckpointFile () const;

Expand Down
7 changes: 6 additions & 1 deletion Source/ERF.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -190,7 +190,6 @@ ERF::ERF ()
pp_inc.resize(nlevs_max);
#endif


rU_new.resize(nlevs_max);
rV_new.resize(nlevs_max);
rW_new.resize(nlevs_max);
Expand Down Expand Up @@ -306,6 +305,12 @@ ERF::ERF ()
}
#endif

// Dirichlet BC data
int nfaces = 6;
xvel_bc_data.resize(nfaces);
yvel_bc_data.resize(nfaces);
zvel_bc_data.resize(nfaces);

// Initialize tagging criteria for mesh refinement
refinement_criteria_setup();

Expand Down
22 changes: 18 additions & 4 deletions Source/ERF_make_new_arrays.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -480,19 +480,33 @@ ERF::initialize_integrator (int lev, MultiFab& cons_mf, MultiFab& vel_mf)
void
ERF::initialize_bcs (int lev)
{
// Dirichlet BC data only lives on level 0
int nface = 6; // 6 faces
Vector<Real*> u_bc_tmp; u_bc_tmp.resize(nface,nullptr);
Vector<Real*> v_bc_tmp; v_bc_tmp.resize(nface,nullptr);
Vector<Real*> w_bc_tmp; w_bc_tmp.resize(nface,nullptr);
if (lev==0) {
for (int iface(0); iface<nface; ++iface) {
u_bc_tmp[iface] = xvel_bc_data[iface].data();
v_bc_tmp[iface] = yvel_bc_data[iface].data();
w_bc_tmp[iface] = zvel_bc_data[iface].data();
}
}
physbcs_cons[lev] = std::make_unique<ERFPhysBCFunct_cons> (lev, geom[lev], domain_bcs_type, domain_bcs_type_d,
m_bc_extdir_vals, m_bc_neumann_vals,
z_phys_nd[lev], use_real_bcs);
physbcs_u[lev] = std::make_unique<ERFPhysBCFunct_u> (lev, geom[lev], domain_bcs_type, domain_bcs_type_d,
m_bc_extdir_vals, m_bc_neumann_vals,
z_phys_nd[lev], use_real_bcs);
z_phys_nd[lev], use_real_bcs, u_bc_tmp);
physbcs_v[lev] = std::make_unique<ERFPhysBCFunct_v> (lev, geom[lev], domain_bcs_type, domain_bcs_type_d,
m_bc_extdir_vals, m_bc_neumann_vals,
z_phys_nd[lev], use_real_bcs);
z_phys_nd[lev], use_real_bcs, v_bc_tmp);
physbcs_w[lev] = std::make_unique<ERFPhysBCFunct_w> (lev, geom[lev], domain_bcs_type, domain_bcs_type_d,
m_bc_extdir_vals, m_bc_neumann_vals,
solverChoice.terrain_type, z_phys_nd[lev], use_real_bcs);
solverChoice.terrain_type, z_phys_nd[lev],
use_real_bcs, w_bc_tmp);
physbcs_w_no_terrain[lev] = std::make_unique<ERFPhysBCFunct_w_no_terrain>
(lev, geom[lev], domain_bcs_type, domain_bcs_type_d,
m_bc_extdir_vals, m_bc_neumann_vals, use_real_bcs);
m_bc_extdir_vals, m_bc_neumann_vals, use_real_bcs,
w_bc_tmp);
}
Loading
Loading