Skip to content

Refac helpers #901

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

Open
wants to merge 5 commits into
base: master
Choose a base branch
from
Open
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
8 changes: 6 additions & 2 deletions .github/workflows/bench.yml
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
name: 'Benchmark'

on:
on:
pull_request:
pull_request_review:
types: [submitted]
workflow_dispatch:
Expand All @@ -23,7 +24,10 @@ jobs:

self:
name: Georgia Tech | Phoenix (NVHPC)
if: github.repository == 'MFlowCode/MFC' && needs.file-changes.outputs.checkall == 'true' && ${{ github.event.review.state == 'approved' }}
if: ${{ github.repository == 'MFlowCode/MFC' && needs.file-changes.outputs.checkall == 'true' && (
(github.event_name == 'pull_request_review' && github.event.review.state == 'approved') ||
(github.event_name == 'pull_request' && github.event.pull_request.user.login == 'sbryngelson')
) }}
needs: file-changes
strategy:
matrix:
Expand Down
194 changes: 97 additions & 97 deletions src/simulation/m_sim_helpers.f90
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,65 @@

contains

!> Computes the modified dtheta for Fourier filtering in azimuthal direction
!! @param k y coordinate index
!! @param l z coordinate index
!! @return fltr_dtheta Modified dtheta value for cylindrical coordinates
pure function f_compute_filtered_dtheta(k, l) result(fltr_dtheta)
!$acc routine seq
integer, intent(in) :: k, l
real(wp) :: fltr_dtheta
integer :: Nfq

if (grid_geometry == 3) then
if (k == 0) then
fltr_dtheta = 2._wp*pi*y_cb(0)/3._wp

Check warning on line 29 in src/simulation/m_sim_helpers.f90

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_sim_helpers.f90#L29

Added line #L29 was not covered by tests
elseif (k <= fourier_rings) then
Nfq = min(floor(2._wp*real(k, wp)*pi), (p + 1)/2 + 1)
fltr_dtheta = 2._wp*pi*y_cb(k - 1)/real(Nfq, wp)

Check warning on line 32 in src/simulation/m_sim_helpers.f90

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_sim_helpers.f90#L31-L32

Added lines #L31 - L32 were not covered by tests
else
fltr_dtheta = y_cb(k - 1)*dz(l)

Check warning on line 34 in src/simulation/m_sim_helpers.f90

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_sim_helpers.f90#L34

Added line #L34 was not covered by tests
end if
else
fltr_dtheta = 0._wp
end if
end function f_compute_filtered_dtheta

!> Computes inviscid CFL terms for multi-dimensional cases (2D/3D only)
!! @param vel directional velocities
!! @param c mixture speed of sound
!! @param j x coordinate index
!! @param k y coordinate index
!! @param l z coordinate index
!! @return cfl_terms computed CFL terms for 2D/3D cases
pure function f_compute_multidim_cfl_terms(vel, c, j, k, l) result(cfl_terms)
!$acc routine seq
real(wp), dimension(num_vels), intent(in) :: vel
real(wp), intent(in) :: c
integer, intent(in) :: j, k, l
real(wp) :: cfl_terms
real(wp) :: fltr_dtheta

fltr_dtheta = f_compute_filtered_dtheta(k, l)

if (p > 0) then
!3D
if (grid_geometry == 3) then
cfl_terms = min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c), &
fltr_dtheta/(abs(vel(3)) + c))

Check warning on line 63 in src/simulation/m_sim_helpers.f90

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_sim_helpers.f90#L61-L63

Added lines #L61 - L63 were not covered by tests
else
cfl_terms = min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c), &
dz(l)/(abs(vel(3)) + c))
end if
else
!2D
cfl_terms = min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c))
end if
end function f_compute_multidim_cfl_terms

!> Computes enthalpy
!! @param q_prim_vf cell centered primitive variables
!! @param pres mixture pressure
Expand Down Expand Up @@ -77,7 +136,7 @@

E = gamma*pres + pi_inf + 5.e-1_wp*rho*vel_sum + qv

! ENERGY ADJUSTMENTS FOR HYPERELASTIC ENERGY
! Adjust energy for hyperelasticity
if (hyperelasticity) then
E = E + G*q_prim_vf(xiend + 1)%sf(j, k, l)
end if
Expand All @@ -93,8 +152,8 @@
!! @param j x index
!! @param k y index
!! @param l z index
!! @param icfl_sf cell centered inviscid cfl number
!! @param vcfl_sf (optional) cell centered viscous cfl number
!! @param icfl_sf cell-centered inviscid cfl number
!! @param vcfl_sf (optional) cell-centered viscous CFL number
!! @param Rc_sf (optional) cell centered Rc
pure subroutine s_compute_stability_from_dt(vel, c, rho, Re_l, j, k, l, icfl_sf, vcfl_sf, Rc_sf)
!$acc routine seq
Expand All @@ -105,82 +164,48 @@
real(wp), dimension(2), intent(in) :: Re_l
integer, intent(in) :: j, k, l

real(wp) :: fltr_dtheta !<
!! Modified dtheta accounting for Fourier filtering in azimuthal direction.
integer :: Nfq
real(wp) :: fltr_dtheta

if (grid_geometry == 3) then
if (k == 0) then
fltr_dtheta = 2._wp*pi*y_cb(0)/3._wp
elseif (k <= fourier_rings) then
Nfq = min(floor(2._wp*real(k, wp)*pi), (p + 1)/2 + 1)
fltr_dtheta = 2._wp*pi*y_cb(k - 1)/real(Nfq, wp)
else
fltr_dtheta = y_cb(k - 1)*dz(l)
end if
! Inviscid CFL calculation
if (p > 0 .or. n > 0) then
! 2D/3D
icfl_sf(j, k, l) = dt/f_compute_multidim_cfl_terms(vel, c, j, k, l)
else
! 1D
icfl_sf(j, k, l) = (dt/dx(j))*(abs(vel(1)) + c)
end if

if (p > 0) then
!3D
if (grid_geometry == 3) then
icfl_sf(j, k, l) = dt/min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c), &
fltr_dtheta/(abs(vel(3)) + c))
else
icfl_sf(j, k, l) = dt/min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c), &
dz(l)/(abs(vel(3)) + c))
end if

if (viscous) then

! Viscous calculations
if (viscous) then
if (p > 0) then
!3D
if (grid_geometry == 3) then
fltr_dtheta = f_compute_filtered_dtheta(k, l)
vcfl_sf(j, k, l) = maxval(dt/Re_l/rho) &
/min(dx(j), dy(k), fltr_dtheta)**2._wp

Rc_sf(j, k, l) = min(dx(j)*(abs(vel(1)) + c), &
dy(k)*(abs(vel(2)) + c), &
fltr_dtheta*(abs(vel(3)) + c)) &
/maxval(1._wp/Re_l)
else
vcfl_sf(j, k, l) = maxval(dt/Re_l/rho) &
/min(dx(j), dy(k), dz(l))**2._wp

Rc_sf(j, k, l) = min(dx(j)*(abs(vel(1)) + c), &
dy(k)*(abs(vel(2)) + c), &
dz(l)*(abs(vel(3)) + c)) &
/maxval(1._wp/Re_l)
end if

end if

elseif (n > 0) then
!2D
icfl_sf(j, k, l) = dt/min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c))

if (viscous) then

elseif (n > 0) then
!2D
vcfl_sf(j, k, l) = maxval(dt/Re_l/rho)/min(dx(j), dy(k))**2._wp

Rc_sf(j, k, l) = min(dx(j)*(abs(vel(1)) + c), &
dy(k)*(abs(vel(2)) + c)) &
/maxval(1._wp/Re_l)

end if

else
!1D
icfl_sf(j, k, l) = (dt/dx(j))*(abs(vel(1)) + c)

if (viscous) then

else
!1D
vcfl_sf(j, k, l) = maxval(dt/Re_l/rho)/dx(j)**2._wp

Rc_sf(j, k, l) = dx(j)*(abs(vel(1)) + c)/maxval(1._wp/Re_l)

end if

end if

end subroutine s_compute_stability_from_dt
Expand All @@ -202,64 +227,39 @@
integer, intent(in) :: j, k, l

real(wp) :: icfl_dt, vcfl_dt
real(wp) :: fltr_dtheta !<
!! Modified dtheta accounting for Fourier filtering in azimuthal direction.

integer :: Nfq
real(wp) :: fltr_dtheta

if (grid_geometry == 3) then
if (k == 0) then
fltr_dtheta = 2._wp*pi*y_cb(0)/3._wp
elseif (k <= fourier_rings) then
Nfq = min(floor(2._wp*real(k, wp)*pi), (p + 1)/2 + 1)
fltr_dtheta = 2._wp*pi*y_cb(k - 1)/real(Nfq, wp)
else
fltr_dtheta = y_cb(k - 1)*dz(l)
end if
! Inviscid CFL calculation
if (p > 0 .or. n > 0) then
! 2D/3D cases
icfl_dt = cfl_target*f_compute_multidim_cfl_terms(vel, c, j, k, l)
else
! 1D case
icfl_dt = cfl_target*(dx(j)/(abs(vel(1)) + c))

Check warning on line 238 in src/simulation/m_sim_helpers.f90

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_sim_helpers.f90#L238

Added line #L238 was not covered by tests
end if

if (p > 0) then
!3D
if (grid_geometry == 3) then
icfl_dt = cfl_target*min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c), &
fltr_dtheta/(abs(vel(3)) + c))
else
icfl_dt = cfl_target*min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c), &
dz(l)/(abs(vel(3)) + c))
end if

if (viscous) then
! Viscous calculations
if (viscous) then
if (p > 0) then
!3D
if (grid_geometry == 3) then
fltr_dtheta = f_compute_filtered_dtheta(k, l)

Check warning on line 246 in src/simulation/m_sim_helpers.f90

View check run for this annotation

Codecov / codecov/patch

src/simulation/m_sim_helpers.f90#L246

Added line #L246 was not covered by tests
vcfl_dt = cfl_target*(min(dx(j), dy(k), fltr_dtheta)**2._wp) &
/minval(1/(rho*Re_l))
else
vcfl_dt = cfl_target*(min(dx(j), dy(k), dz(l))**2._wp) &
/minval(1/(rho*Re_l))
end if
end if

elseif (n > 0) then
!2D
icfl_dt = cfl_target*min(dx(j)/(abs(vel(1)) + c), &
dy(k)/(abs(vel(2)) + c))

if (viscous) then
elseif (n > 0) then
!2D
vcfl_dt = cfl_target*(min(dx(j), dy(k))**2._wp)/maxval((1/Re_l)/rho)
end if

else
!1D
icfl_dt = cfl_target*(dx(j)/(abs(vel(1)) + c))

if (viscous) then
else
!1D
vcfl_dt = cfl_target*(dx(j)**2._wp)/minval(1/(rho*Re_l))
end if

end if

if (any(re_size > 0)) then
if (any(Re_size > 0)) then
max_dt(j, k, l) = min(icfl_dt, vcfl_dt)
else
max_dt(j, k, l) = icfl_dt
Expand Down
Loading