Skip to content

Commit 0bbddda

Browse files
aliencaocaowilfonba
aliencaocao
and
wilfonba
authored
Fix the azimuthal flux divergence and geometric sources in 3D cylindrical coordinates (#2051)
Co-authored-by: Ben Wilfong <48168887+wilfonba@users.noreply.github.com>
1 parent 07f3d5e commit 0bbddda

10 files changed

Lines changed: 601 additions & 42 deletions

File tree

‎src/simulation/m_rhs.fpp‎

Lines changed: 7 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -1070,7 +1070,7 @@ contains
10701070
integer :: k_loop, l_loop, q_loop !< Standardized spatial loop iterators 0:m, 0:n, 0:p
10711071
integer :: i_fluid_loop
10721072
real(wp) :: inv_ds, flux_face1, flux_face2
1073-
real(wp) :: advected_qty_val, pressure_val, velocity_val
1073+
real(wp) :: advected_qty_val, pressure_val
10741074
real(wp) :: G1_eff, G2_eff, pres_K, alpha_K, alpha_rho_K, blkmod_K
10751075

10761076
G1_eff = 0._wp
@@ -1332,16 +1332,15 @@ contains
13321332
end if
13331333
13341334
if (grid_geometry == 3) then ! Cylindrical Coordinates
1335-
$:GPU_PARALLEL_LOOP(collapse=4,private='[j, k, l, q, inv_ds, velocity_val, flux_face1, flux_face2]')
1335+
$:GPU_PARALLEL_LOOP(collapse=4,private='[j, k, l, q, inv_ds, flux_face1, flux_face2]')
13361336
do j = 1, sys_size
13371337
do k = 0, p
13381338
do q = 0, n
13391339
do l = 0, m
13401340
inv_ds = 1._wp/(dz(k)*y_cc(q))
1341-
velocity_val = q_prim_vf%vf(eqn_idx%cont%end + idir)%sf(l, q, k)
13421341
flux_face1 = flux_n(3)%vf(j)%sf(l, q, k - 1)
13431342
flux_face2 = flux_n(3)%vf(j)%sf(l, q, k)
1344-
rhs_vf(j)%sf(l, q, k) = rhs_vf(j)%sf(l, q, k) + inv_ds*velocity_val*(flux_face1 - flux_face2)
1343+
rhs_vf(j)%sf(l, q, k) = rhs_vf(j)%sf(l, q, k) + inv_ds*(flux_face1 - flux_face2)
13451344
end do
13461345
end do
13471346
end do
@@ -1624,7 +1623,7 @@ contains
16241623
do k_idx = 0, p ! z_extent
16251624
do q_idx = 0, n ! y_extent
16261625
do l_idx = 0, m ! x_extent
1627-
local_inv_ds = 1._wp/dz(k_idx)
1626+
local_inv_ds = 1._wp/(dz(k_idx)*y_cc(q_idx))
16281627
local_term_coeff = q_prim_vf_arg%vf(eqn_idx%cont%end + current_idir)%sf(l_idx, q_idx, k_idx)
16291628
local_flux1 = flux_src_n_vf_arg%vf(j_adv)%sf(l_idx, q_idx, k_idx - 1)
16301629
local_flux2 = flux_src_n_vf_arg%vf(j_adv)%sf(l_idx, q_idx, k_idx)
@@ -1639,7 +1638,7 @@ contains
16391638
$:GPU_PARALLEL_LOOP(collapse=3, private='[k_idx, l_idx, q_idx, local_inv_ds, local_k_term_val, &
16401639
& local_flux1, local_flux2]')
16411640
do k_idx = 0, p; do q_idx = 0, n; do l_idx = 0, m
1642-
local_inv_ds = 1._wp/dz(k_idx)
1641+
local_inv_ds = 1._wp/(dz(k_idx)*y_cc(q_idx))
16431642
local_k_term_val = Kterm_arg(l_idx, q_idx, k_idx)
16441643
local_flux1 = nc_iface_vel_n(3)%vf(3)%sf(l_idx, q_idx, k_idx)
16451644
local_flux2 = nc_iface_vel_n(3)%vf(3)%sf(l_idx, q_idx, k_idx - 1)
@@ -1656,7 +1655,7 @@ contains
16561655
& local_flux1, local_flux2]')
16571656
do j_adv = eqn_idx%adv%beg, eqn_idx%adv%end
16581657
do k_idx = 0, p; do q_idx = 0, n; do l_idx = 0, m
1659-
local_inv_ds = 1._wp/dz(k_idx)
1658+
local_inv_ds = 1._wp/(dz(k_idx)*y_cc(q_idx))
16601659
local_term_coeff = q_cons_vf_arg%vf(j_adv)%sf(l_idx, q_idx, k_idx)
16611660
local_flux1 = flux_src_n_vf_arg%vf(eqn_idx%adv%beg)%sf(l_idx, q_idx, k_idx)
16621661
local_flux2 = flux_src_n_vf_arg%vf(eqn_idx%adv%beg)%sf(l_idx, q_idx, k_idx - 1)
@@ -1669,7 +1668,7 @@ contains
16691668
$:GPU_PARALLEL_LOOP(collapse=3, private='[k_idx, l_idx, q_idx, local_inv_ds, local_k_term_val, &
16701669
& local_flux1, local_flux2]')
16711670
do k_idx = 0, p; do q_idx = 0, n; do l_idx = 0, m
1672-
local_inv_ds = 1._wp/dz(k_idx)
1671+
local_inv_ds = 1._wp/(dz(k_idx)*y_cc(q_idx))
16731672
local_k_term_val = Kterm_arg(l_idx, q_idx, k_idx)
16741673
local_flux1 = flux_src_n_vf_arg%vf(eqn_idx%adv%beg)%sf(l_idx, q_idx, k_idx)
16751674
local_flux2 = flux_src_n_vf_arg%vf(eqn_idx%adv%beg)%sf(l_idx, q_idx, k_idx - 1)

‎src/simulation/m_riemann_solver_hll.fpp‎

Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -695,6 +695,17 @@ contains
695695
end if
696696
end if
697697
#:endif
698+
#:if (NORM_DIR == 3)
699+
if (grid_geometry == 3) then
700+
$:GPU_LOOP(parallelism='[seq]')
701+
do i = 1, sys_size
702+
flux_gsrc_rsx_vf(${SF('')}$, i) = 0._wp
703+
end do
704+
flux_gsrc_rsx_vf(${SF('')}$, eqn_idx%mom%beg + 1) = -flux_rsx_vf(${SF('')}$, &
705+
& eqn_idx%mom%end) + (s_M*pres_R - s_P*pres_L)/(s_M - s_P)
706+
flux_gsrc_rsx_vf(${SF('')}$, eqn_idx%mom%end) = flux_rsx_vf(${SF('')}$, eqn_idx%mom%beg + 1)
707+
end if
708+
#:endif
698709
end do
699710
end do
700711
end do

‎src/simulation/m_riemann_solver_lf.fpp‎

Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -315,6 +315,17 @@ contains
315315
end do
316316
end if
317317
#:endif
318+
#:if (NORM_DIR == 3)
319+
if (grid_geometry == 3) then
320+
$:GPU_LOOP(parallelism='[seq]')
321+
do i = 1, sys_size
322+
flux_gsrc_rsx_vf(${SF('')}$, i) = 0._wp
323+
end do
324+
flux_gsrc_rsx_vf(${SF('')}$, eqn_idx%mom%beg + 1) = -flux_rsx_vf(${SF('')}$, &
325+
& eqn_idx%mom%end) + (s_M*pres_R - s_P*pres_L)/(s_M - s_P)
326+
flux_gsrc_rsx_vf(${SF('')}$, eqn_idx%mom%end) = flux_rsx_vf(${SF('')}$, eqn_idx%mom%beg + 1)
327+
end if
328+
#:endif
318329
end do
319330
end do
320331
end do

‎tests/4C1D77D7/golden-metadata.txt‎

Lines changed: 162 additions & 0 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

‎tests/4C1D77D7/golden.txt‎

Lines changed: 16 additions & 0 deletions
Some generated files are not rendered by default. Learn more about customizing how changed files appear on GitHub.

0 commit comments

Comments
 (0)