Skip to content

Commit

Permalink
GRCBC impl in 3 dims
Browse files Browse the repository at this point in the history
  • Loading branch information
Anand committed Nov 6, 2024
1 parent 63c79cb commit 2b79939
Show file tree
Hide file tree
Showing 6 changed files with 167 additions and 9 deletions.
4 changes: 4 additions & 0 deletions src/common/m_derived_types.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,10 @@ module m_derived_types
real(kind(0d0)) :: ve1
real(kind(0d0)) :: ve2
real(kind(0d0)) :: ve3
real(kind(0d0)) :: u_in, v_in, w_in, u_out, v_out, w_out, pres_in, pres_out
real(kind(0d0)), dimension(num_fluids_max) :: alpha_rho_in, alpha_in
logical :: grcbc_in, grcbc_out, grcbc_vel_out

end type int_bounds_info

!> Derived type adding beginning (beg) and end bounds info as attributes
Expand Down
106 changes: 106 additions & 0 deletions src/simulation/m_cbc.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -134,6 +134,31 @@ module m_cbc
integer :: cbc_dir, cbc_loc
!$acc declare create(dj, bcxb, bcxe, bcyb, bcye, bczb, bcze, cbc_dir, cbc_loc)

real(kind(0d0)) :: ux_in, ux_out, vx_in, vx_out, wx_in, wx_out, presx_in, presx_out, Delx_in, Delx_out

real(kind(0d0)) :: uy_in, uy_out, vy_in, vy_out, wy_in, wy_out, presy_in, presy_out, Dely_in, Dely_out

real(kind(0d0)) :: uz_in, uz_out, vz_in, vz_out, wz_in, wz_out, presz_in, presz_out, Delz_in, Delz_out

!$acc declare create(ux_in, ux_out, vx_in, vx_out, wx_in, wx_out, presx_in, presx_out, Delx_in, Delx_out)
!$acc declare create(uy_in, uy_out, vy_in, vy_out, wy_in, wy_out, presy_in, presy_out, Dely_in, Dely_out)
!$acc declare create(uz_in, uz_out, vz_in, vz_out, wz_in, wz_out, presz_in, presz_out, Delz_in, Delz_out)


#ifdef CRAY_ACC_WAR
@:CRAY_DECLARE_GLOBAL(real(kind(0d0)), dimension(:, :, :), alpha_rhox_in, alphax_in)
@:CRAY_DECLARE_GLOBAL(real(kind(0d0)), dimension(:, :, :), alpha_rhoy_in, alphay_in)
@:CRAY_DECLARE_GLOBAL(real(kind(0d0)), dimension(:, :, :), alpha_rhoz_in, alphaz_in)
!$acc declare link(alpha_rhox_in, alphax_in, alpha_rhoy_in, alphay_in, alpha_rhoz_in, alphaz_in)
#else
real(kind(0d0)), allocatable, dimension(:) :: alpha_rhox_in, alphax_in
real(kind(0d0)), allocatable, dimension(:) :: alpha_rhoy_in, alphay_in
real(kind(0d0)), allocatable, dimension(:) :: alpha_rhoz_in, alphaz_in
!$acc declare create(alpha_rhox_in, alphax_in, alpha_rhoy_in, alphay_in, alpha_rhoz_in, alphaz_in)
#endif



#ifndef CRAY_ACC_WAR
!$acc declare create(q_prim_rsx_vf, q_prim_rsy_vf, q_prim_rsz_vf, F_rsx_vf, F_src_rsx_vf,flux_rsx_vf, flux_src_rsx_vf, &
!$acc F_rsy_vf, F_src_rsy_vf,flux_rsy_vf, flux_src_rsy_vf, F_rsz_vf, F_src_rsz_vf,flux_rsz_vf, flux_src_rsz_vf,Re, &
Expand Down Expand Up @@ -426,6 +451,58 @@ contains
!$acc update device(bczb, bcze)
end if

@:ALLOCATE_GLOBAL(alpha_rhox_in(1:num_fluids), alpha_rhoy_in(1:num_fluids), alpha_rhoz_in(1:num_fluids))
@:ALLOCATE_GLOBAL(alphax_in(1:num_fluids), alphay_in(1:num_fluids), alphaz_in(1:num_fluids))

#:for CBC_DIR, XYZ in [(1, 'x'), (2, 'y'), (3, 'z')]

if(${CBC_DIR}$ == 1) then
u${XYZ}$_in = bc_${XYZ}$%u_in
v${XYZ}$_in = bc_${XYZ}$%v_in
w${XYZ}$_in = bc_${XYZ}$%w_in
u${XYZ}$_out = bc_${XYZ}$%u_out
v${XYZ}$_out = bc_${XYZ}$%v_out
w${XYZ}$_out = bc_${XYZ}$%w_out
Del${XYZ}$_in = maxval(dx)
Del${XYZ}$_out = maxval(dx)
else if(${CBC_DIR}$ == 2) then

u${XYZ}$_in = bc_${XYZ}$%v_in
v${XYZ}$_in = bc_${XYZ}$%u_in
w${XYZ}$_in = bc_${XYZ}$%w_in
u${XYZ}$_out = bc_${XYZ}$%v_out
v${XYZ}$_out = bc_${XYZ}$%u_out
w${XYZ}$_out = bc_${XYZ}$%w_out
if(n > 0) then
Del${XYZ}$_in = maxval(dy)
Del${XYZ}$_out = maxval(dy)
end if
else if(${CBC_DIR}$ == 3) then
u${XYZ}$_in = bc_${XYZ}$%w_in
v${XYZ}$_in = bc_${XYZ}$%u_in
w${XYZ}$_in = bc_${XYZ}$%v_in
u${XYZ}$_out = bc_${XYZ}$%w_out
v${XYZ}$_out = bc_${XYZ}$%u_out
w${XYZ}$_out = bc_${XYZ}$%v_out
if(p > 0) then
Del${XYZ}$_in = maxval(dz)
Del${XYZ}$_out = maxval(dz)
end if
end if

pres${XYZ}$_in = bc_${XYZ}$%pres_in
pres${XYZ}$_out = bc_${XYZ}$%pres_out
do i = 1, num_fluids
alpha_rho${XYZ}$_in(i) = bc_${XYZ}$%alpha_rho_in(i)
alpha${XYZ}$_in(i) = bc_${XYZ}$%alpha_in(i)
end do
!$acc update device(u${XYZ}$_in, v${XYZ}$_in, w${XYZ}$_in, u${XYZ}$_out, v${XYZ}$_out, w${XYZ}$_out)
!$acc update device(pres${XYZ}$_in, pres${XYZ}$_out, alpha_rho${XYZ}$_in, alpha${XYZ}$_in)
!$acc update device(Del${XYZ}_in, Del${XYZ}_out)

#:endfor


end subroutine s_initialize_cbc_module

!> Compute CBC coefficients
Expand Down Expand Up @@ -655,6 +732,7 @@ contains
real(kind(0d0)) :: pi_inf !< Cell averaged liquid stiffness
real(kind(0d0)) :: qv !< Cell averaged fluid reference energy
real(kind(0d0)) :: c
real(kind(0d0)) :: Ma

real(kind(0d0)) :: vel_K_sum, vel_dv_dt_sum

Expand Down Expand Up @@ -861,14 +939,42 @@ contains
lambda(2) = vel(dir_idx(1))
lambda(3) = vel(dir_idx(1)) + c

Ma = vel(dir_idx(1)) / c

if ((cbc_loc == -1 .and. bc${XYZ}$b == -5) .or. (cbc_loc == 1 .and. bc${XYZ}$e == -5)) then
call s_compute_slip_wall_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds)
else if ((cbc_loc == -1 .and. bc${XYZ}$b == -6) .or. (cbc_loc == 1 .and. bc${XYZ}$e == -6)) then
call s_compute_nonreflecting_subsonic_buffer_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds)
else if ((cbc_loc == -1 .and. bc${XYZ}$b == -7) .or. (cbc_loc == 1 .and. bc${XYZ}$e == -7)) then
call s_compute_nonreflecting_subsonic_inflow_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds)

if(bc_${XYZ}$%grcbc_in) then
!$acc loop seq
do i = 2, momxb
L(2) = c**3d0 * Ma * (alpha_rho(i-1) - alpha_rho${XYZ}$_in(i-1)) / Del${XYZ}$_in - c * Ma * (pres - pres${XYZ}$_in) / Del${XYZ}$_in
end do
if(n > 0) then
L(momxb + 1) = c * Ma * (vel(dir_idx(2)) - v${XYZ}$_in) / Del${XYZ}$_in
if(p > 0) then
L(momxb + 2) = c * Ma * (vel(dir_idx(3)) - w${XYZ}$_in) / Del${XYZ}$_in
end if
end if
!$acc loop seq
do i = E_idx, advxe - 1
L(i) = c * Ma * (adv(i + 1 - E_idx) - alpha${XYZ}$_in(i + 1 - E_idx)) / Del${XYZ}$_in
end do
L(advxe) = rho*c**2d0*(1d0 + Ma)*(vel(dir_idx(1)) + u${XYZ}$_in * sign(1, cbc_loc)) / Del${XYZ}$_in + c*(1d0 + Ma)*(pres - pres${XYZ}$_in) / Del${XYZ}$_in
end if
else if ((cbc_loc == -1 .and. bc${XYZ}$b == -8) .or. (cbc_loc == 1 .and. bc${XYZ}$e == -8)) then
call s_compute_nonreflecting_subsonic_outflow_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds)

if(bc_${XYZ}$%grcbc_out) then
L(advxe) = c*(1d0 - Ma)*(pres - pres${XYZ}$_out) / Del${XYZ}$_out
if(bc_${XYZ}$%grcbc_vel_out) then
L(advxe) = L(advxe) + rho*c**2d0*(1d0 - Ma)*(vel(dir_idx(1)) + u${XYZ}$_out * sign(1, cbc_loc)) / Del${XYZ}$_out
end if
end if

else if ((cbc_loc == -1 .and. bc${XYZ}$b == -9) .or. (cbc_loc == 1 .and. bc${XYZ}$e == -9)) then
call s_compute_force_free_subsonic_outflow_L(lambda, L, rho, c, mf, dalpha_rho_ds, dpres_ds, dvel_ds, dadv_ds)
else if ((cbc_loc == -1 .and. bc${XYZ}$b == -10) .or. (cbc_loc == 1 .and. bc${XYZ}$e == -10)) then
Expand Down
13 changes: 13 additions & 0 deletions src/simulation/m_global_parameters.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -473,6 +473,7 @@ module m_global_parameters

!$acc declare create(pb_ts, mv_ts)
#endif

! ======================================================================

contains
Expand Down Expand Up @@ -680,6 +681,18 @@ contains
integral(i)%ymax = dflt_real
end do

bc_x%grcbc_in = .false.
bc_x%grcbc_out = .false.
bc_x%grcbc_vel_out = .false.

bc_y%grcbc_in = .false.
bc_y%grcbc_out = .false.
bc_y%grcbc_vel_out = .false.

bc_z%grcbc_in = .false.
bc_z%grcbc_out = .false.
bc_z%grcbc_vel_out = .false.

end subroutine s_assign_default_values_to_user_inputs

!> The computation of parameters, the allocation of memory,
Expand Down
15 changes: 15 additions & 0 deletions src/simulation/m_mpi_proxy.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -203,15 +203,24 @@ contains
& 'polydisperse', 'qbmm', 'acoustic_source', 'probe_wrt', 'integral_wrt', &
& 'prim_vars_wrt', 'weno_avg', 'file_per_process', 'relax', &
& 'adv_n', 'adap_dt', 'ib', 'bodyForces', 'bf_x', 'bf_y', 'bf_z', &
& 'bc_x%grcbc_in', 'bc_x%grcbc_out', 'bc_x%grcbc_vel_out', &
& 'bc_y%grcbc_in', 'bc_y%grcbc_out', 'bc_y%grcbc_vel_out', &
& 'bc_z%grcbc_in', 'bc_z%grcbc_out', 'bc_z%grcbc_vel_out', &
& 'cfl_adap_dt', 'cfl_const_dt', 'cfl_dt' ]
call MPI_BCAST(${VAR}$, 1, MPI_LOGICAL, 0, MPI_COMM_WORLD, ierr)
#:endfor
#:for VAR in [ 'dt','weno_eps','teno_CT','pref','rhoref','R0ref','Web','Ca', 'sigma', &
& 'Re_inv', 'poly_sigma', 'palpha_eps', 'ptgalpha_eps', 'pi_fac', &
& 'bc_x%vb1','bc_x%vb2','bc_x%vb3','bc_x%ve1','bc_x%ve2','bc_x%ve2', &
& 'bc_x%u_in','bc_x%v_in','bc_x%w_in','bc_x%u_out','bc_x%v_out','bc_x%w_out', &
& 'bc_x%pres_in','bc_x%pres_out', &
& 'bc_y%vb1','bc_y%vb2','bc_y%vb3','bc_y%ve1','bc_y%ve2','bc_y%ve3', &
& 'bc_y%u_in','bc_y%v_in','bc_y%w_in','bc_y%u_out','bc_y%v_out','bc_y%w_out', &
& 'bc_y%pres_in','bc_y%pres_out', &
& 'bc_z%vb1','bc_z%vb2','bc_z%vb3','bc_z%ve1','bc_z%ve2','bc_z%ve3', &
& 'bc_z%u_in','bc_z%v_in','bc_z%w_in','bc_z%u_out','bc_z%v_out','bc_z%w_out', &
& 'bc_z%pres_in','bc_z%pres_out', &
& 'x_domain%beg', 'x_domain%end', 'y_domain%beg', 'y_domain%end', &
& 'z_domain%beg', 'z_domain%end', 't_stop', 't_save', 'cfl_target']
call MPI_BCAST(${VAR}$, 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr)
Expand All @@ -234,6 +243,12 @@ contains
call MPI_BCAST(fluid_pp(i)%Re(1), 2, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr)
end do
do i = 1, num_fluids_max
#:for VAR in ['bc_x%alpha_rho_in','bc_x%alpha_in','bc_y%alpha_rho_in','bc_y%alpha_in','bc_z%alpha_rho_in','bc_z%alpha_in']
call MPI_BCAST(${VAR}$(i), 1, MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr )
#:endfor
end do
do i = 1, num_ibs
#:for VAR in [ 'radius', 'length_x', 'length_y', &
& 'x_centroid', 'y_centroid', 'c', 'm', 'p', 't', 'theta', 'slip' ]
Expand Down
3 changes: 3 additions & 0 deletions src/simulation/m_start_up.fpp
Original file line number Diff line number Diff line change
Expand Up @@ -1474,6 +1474,9 @@ contains
!$acc update device(bc_y%vb1, bc_y%vb2, bc_y%vb3, bc_y%ve1, bc_y%ve2, bc_y%ve3)
!$acc update device(bc_z%vb1, bc_z%vb2, bc_z%vb3, bc_z%ve1, bc_z%ve2, bc_z%ve3)

!$acc update device(bc_x%grcbc_in, bc_x%grcbc_out, bc_x%grcbc_vel_out)
!$acc update device(bc_y%grcbc_in, bc_y%grcbc_out, bc_y%grcbc_vel_out)
!$acc update device(bc_z%grcbc_in, bc_z%grcbc_out, bc_z%grcbc_vel_out)

!$acc update device(relax, relax_model)
if (relax) then
Expand Down
35 changes: 26 additions & 9 deletions toolchain/mfc/run/case_dicts.py
Original file line number Diff line number Diff line change
Expand Up @@ -237,15 +237,31 @@ def analytic(self):
SIMULATION[f'patch_ib({ib_id})%{cmp}_centroid'] = ParamType.REAL
SIMULATION[f'patch_ib({ib_id})%length_{cmp}'] = ParamType.REAL

for cmp in ["x", "y", "z"]:
SIMULATION[f'bc_{cmp}%beg'] = ParamType.INT
SIMULATION[f'bc_{cmp}%end'] = ParamType.INT
SIMULATION[f'bc_{cmp}%vb1'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%vb2'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%vb3'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%ve1'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%ve2'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%ve3'] = ParamType.REAL
for cmp in ["x", "y", "z"]:
SIMULATION[f'bc_{cmp}%beg'] = ParamType.INT
SIMULATION[f'bc_{cmp}%end'] = ParamType.INT
SIMULATION[f'bc_{cmp}%vb1'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%vb2'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%vb3'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%ve1'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%ve2'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%ve3'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%u_in'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%v_in'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%w_in'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%u_out'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%v_out'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%w_out'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%pres_in'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%pres_out'] = ParamType.REAL
SIMULATION[f'bc_{cmp}%grcbc_in'] = ParamType.LOG
SIMULATION[f'bc_{cmp}%grcbc_out'] = ParamType.LOG
SIMULATION[f'bc_{cmp}%grcbc_vel_out'] = ParamType.LOG

for int_id in range(1, 10+1):
for cmp in ["x", "y", "z"]:
SIMULATION[f"bc_{cmp}%alpha_rho_in({int_id})"] = ParamType.REAL
SIMULATION[f"bc_{cmp}%alpha_in({int_id})"] = ParamType.REAL

for var in ["k", "w", "p", "g"]:
SIMULATION[f'{var}_{cmp}'] = ParamType.REAL
Expand Down Expand Up @@ -295,6 +311,7 @@ def analytic(self):
SIMULATION[f"integral({int_id})%{cmp}max"] = ParamType.REAL



# Removed: 'fourier_modes%beg', 'fourier_modes%end', 'chem_wrt'
# Feel free to return them if they are needed once more.
POST_PROCESS = COMMON.copy()
Expand Down

0 comments on commit 2b79939

Please sign in to comment.