Skip to content

Refactor m_riemann_solvers Module (Non-Solver Subroutines) #908

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.

Sign up for GitHub

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 25 commits into
base: master
Choose a base branch
from
Open

Refactor m_riemann_solvers Module (Non-Solver Subroutines) #908

Show file tree
Hide file tree
Changes from 8 commits
Commits
Show all changes
25 commits
Select commit Hold shift + click to select a range
967c128
s_compute_wave_speed
Malmahrouqi3 Jun 29, 2025
89280ad
low-stakes refactoring & lint/format
Malmahrouqi3 Jun 29, 2025
3bfa9f4
fixed few errors
Malmahrouqi3 Jun 29, 2025
ddaa1b6
quick fix
Malmahrouqi3 Jun 29, 2025
12eabf4
simple fix
Malmahrouqi3 Jun 29, 2025
4e2a439
removed &
Malmahrouqi3 Jun 29, 2025
b684a4e
remove allocatable from pointers
Malmahrouqi3 Jun 29, 2025
6fb8c65
m_riemann_solvers cleanup & format/lint
Malmahrouqi3 Jun 29, 2025
dd158d0
Merge branch 'master' into refactor
sbryngelson Jul 1, 2025
46de417
Merge branch 'master' into refactor
sbryngelson Jul 3, 2025
dcfc3e0
lint-fix
Malmahrouqi3 Jul 11, 2025
7f0c19d
Merge remote-tracking branch 'upstream/master' into refactor
Malmahrouqi3 Aug 14, 2025
6e7281b
format
Aug 14, 2025
3ad2ddb
lint/format + few fixes
Aug 18, 2025
dc8019f
removed raw directive in s_compute_wave_speed
Aug 18, 2025
f3cde3e
corrected something
Aug 18, 2025
156ec55
Merge branch 'MFlowCode:master' into refactor
Malmahrouqi3 Aug 18, 2025
b6a352e
Merge branch 'master' into refactor
Malmahrouqi3 Sep 1, 2025
3c4b4e1
Merge branch 'master' into refactor
Malmahrouqi3 Sep 6, 2025
797c40f
Merge branch 'master' into refactor
Malmahrouqi3 Sep 11, 2025
676082e
undid all changes, retained s_compute_wave_speed
Sep 15, 2025
f5e8c52
minor changes
Sep 15, 2025
e9bc89c
Merge branch 'master' into refactor
Malmahrouqi3 Sep 15, 2025
89064af
Merge branch 'master' into refactor
Malmahrouqi3 Sep 17, 2025
5ba2e99
Merge branch 'master' into refactor
sbryngelson Oct 5, 2025
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
115 changes: 115 additions & 0 deletions src/common/m_variables_conversion.fpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -44,6 +44,7 @@
#ifndef MFC_PRE_PROCESS
s_compute_speed_of_sound, &
s_compute_fast_magnetosonic_speed, &
s_compute_wave_speed, &
#endif
s_finalize_variables_conversion_module

Expand Down Expand Up @@ -1709,4 +1710,118 @@
end subroutine s_compute_fast_magnetosonic_speed
#endif

#ifndef MFC_PRE_PROCESS
subroutine s_compute_wave_speed(wave_speeds, vel_L, vel_R, pres_L, pres_R, rho_L, rho_R, rho_avg, &

Check warning on line 1714 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1714 

Added line #L1714 was not covered by tests
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

surely some of these arguments should be optional since not every riemann solver uses them all?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

yup yup, not all of them will be used in every scenario. I will look into a way to reduce that as the subroutine call looks the same in every solver yet not neat-looking

                            
call s_compute_wave_speed(wave_speeds, vel_L, vel_R, pres_L, pres_R, rho_L, rho_R, rho_avg, &
                     c_L, c_R, c_avg, c_fast%L, c_fast%R, G_L, G_R, &
                     tau_e_L, tau_e_R, gamma_L, gamma_R, pi_inf_L, pi_inf_R, &
                     s_L, s_R, s_S, s_M, s_P, dir_idx(1), dir_idx_tau(1))

c_L, c_R, c_avg, c_fast_L, c_fast_R, G_L, G_R, &
tau_e_L, tau_e_R, gamma_L, gamma_R, pi_inf_L, pi_inf_R, &

Check warning on line 1716 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1716 

Added line #L1716 was not covered by tests
s_L, s_R, s_S, s_M, s_P, idx, idx_tau)

! Computes the wave speeds for the Riemann solver
#ifdef _CRAYFTN
!DIR$ INLINEALWAYS s_compute_wave_speed
#else
!$acc routine seq
#endif

! Input parameters
integer, intent(in) :: wave_speeds
integer, intent(in) :: idx, idx_tau
real(wp), intent(in) :: rho_L, rho_R
real(wp), dimension(:), intent(in) :: vel_L, vel_R, tau_e_L, tau_e_R
real(wp), intent(in) :: pres_L, pres_R, c_L, c_R
real(wp), intent(in) :: gamma_L, gamma_R, pi_inf_L, pi_inf_R
real(wp), intent(in) :: rho_avg, c_avg
real(wp), intent(in) :: c_fast_L, c_fast_R
real(wp), intent(in) :: G_L, G_R

! Local variables
real(wp) :: pres_SL, pres_SR, Ms_L, Ms_R

! Output parameters
real(wp), intent(out) :: s_L, s_R, s_S, s_M, s_P

if (wave_speeds == 1) then
if (elasticity) then
s_L = min(vel_L(idx) - sqrt(c_L*c_L + &
(((4_wp*G_L)/3_wp) + tau_e_L(idx_tau))/rho_L), vel_R(idx) - sqrt(c_R*c_R + &
(((4_wp*G_R)/3_wp) + tau_e_R(idx_tau))/rho_R))
s_R = max(vel_R(idx) + sqrt(c_R*c_R + &
(((4_wp*G_R)/3_wp) + tau_e_R(idx_tau))/rho_R), vel_L(idx) + sqrt(c_L*c_L + &
(((4_wp*G_L)/3_wp) + tau_e_L(idx_tau))/rho_L))
s_S = (pres_R - tau_e_R(idx_tau) - pres_L + &
tau_e_L(idx_tau) + rho_L*vel_L(idx)*(s_L - vel_L(idx)) - &
rho_R*vel_R(idx)*(s_R - vel_R(idx)))/(rho_L*(s_L - vel_L(idx)) - &
rho_R*(s_R - vel_R(idx)))
else if (mhd) then
s_L = min(vel_L(idx) - c_fast_L, vel_R(idx) - c_fast_R)
s_R = max(vel_R(idx) + c_fast_R, vel_L(idx) + c_fast_L)
s_S = (pres_R - pres_L + rho_L*vel_L(idx)* &
(s_L - vel_L(idx)) - rho_R*vel_R(idx)*(s_R - vel_R(idx))) &
/(rho_L*(s_L - vel_L(idx)) - rho_R*(s_R - vel_R(idx)))
else if (hypoelasticity) then
s_L = min(vel_L(idx) - sqrt(c_L*c_L + (((4._wp*G_L)/3._wp) + &
tau_e_L(idx_tau))/rho_L) &
, vel_R(idx) - sqrt(c_R*c_R + (((4._wp*G_R)/3._wp) + &
tau_e_R(idx_tau))/rho_R))
s_R = max(vel_R(idx) + sqrt(c_R*c_R + (((4._wp*G_R)/3._wp) + &
tau_e_R(idx_tau))/rho_R) &
, vel_L(idx) + sqrt(c_L*c_L + (((4._wp*G_L)/3._wp) + &
tau_e_L(idx_tau))/rho_L))
s_S = (pres_R - pres_L + rho_L*vel_L(idx)* &
(s_L - vel_L(idx)) - rho_R*vel_R(idx)*(s_R - vel_R(idx))) &
/(rho_L*(s_L - vel_L(idx)) - rho_R*(s_R - vel_R(idx)))
else if (hyperelasticity) then
s_L = min(vel_L(idx) - sqrt(c_L*c_L + (4._wp*G_L/3._wp)/rho_L) &
, vel_R(idx) - sqrt(c_R*c_R + (4._wp*G_R/3._wp)/rho_R))
s_R = max(vel_R(idx) + sqrt(c_R*c_R + (4._wp*G_R/3._wp)/rho_R) &
, vel_L(idx) + sqrt(c_L*c_L + (4._wp*G_L/3._wp)/rho_L))
s_S = (pres_R - pres_L + rho_L*vel_L(idx)* &
(s_L - vel_L(idx)) - rho_R*vel_R(idx)*(s_R - vel_R(idx))) &
/(rho_L*(s_L - vel_L(idx)) - rho_R*(s_R - vel_R(idx)))

Check warning on line 1780 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1745-L1780 

Added lines #L1745 - L1780 were not covered by tests
else
s_L = min(vel_L(idx) - c_L, vel_R(idx) - c_R)
s_R = max(vel_R(idx) + c_R, vel_L(idx) + c_L)
s_S = (pres_R - pres_L + rho_L*vel_L(idx)* &
(s_L - vel_L(idx)) - rho_R*vel_R(idx)*(s_R - vel_R(idx))) &
/(rho_L*(s_L - vel_L(idx)) - rho_R*(s_R - vel_R(idx)))

Check warning on line 1786 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1782-L1786 

Added lines #L1782 - L1786 were not covered by tests
end if
else if (wave_speeds == 2) then
pres_SL = 5e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL

Check warning on line 1790 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1788-L1790 

Added lines #L1788 - L1790 were not covered by tests
Ms_L = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))

Check warning on line 1793 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1793 

Added line #L1793 was not covered by tests
Ms_R = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))

Check warning on line 1799 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1796-L1799 

Added lines #L1796 - L1799 were not covered by tests
Comment on lines +1803 to +1813
Copy link

Copilot AI Sep 17, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Inconsistent use of floating point literal notation. Line 1804, 1806, 1809, and 1814 use 5e-1_wp while other parts of the codebase use 5.e-1_wp (with decimal point). Should use 5.e-1_wp for consistency with established patterns in the codebase.

Suggested change
pres_SL = 5e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL
Ms_L = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))
Ms_R = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))
pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL
Ms_L = max(1._wp, sqrt(1._wp + ((5.e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))
Ms_R = max(1._wp, sqrt(1._wp + ((5.e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5.e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))

Copilot uses AI. Check for mistakes.

Comment on lines +1803 to +1813
Copy link

Copilot AI Sep 17, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Inconsistent use of floating point literal notation. Line 1804, 1806, 1809, and 1814 use 5e-1_wp while other parts of the codebase use 5.e-1_wp (with decimal point). Should use 5.e-1_wp for consistency with established patterns in the codebase.

Suggested change
pres_SL = 5e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL
Ms_L = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))
Ms_R = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))
pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL
Ms_L = max(1._wp, sqrt(1._wp + ((5.e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))
Ms_R = max(1._wp, sqrt(1._wp + ((5.e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5.e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))

Copilot uses AI. Check for mistakes.

Comment on lines +1803 to +1813
Copy link

Copilot AI Sep 17, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Inconsistent use of floating point literal notation. Line 1804, 1806, 1809, and 1814 use 5e-1_wp while other parts of the codebase use 5.e-1_wp (with decimal point). Should use 5.e-1_wp for consistency with established patterns in the codebase.

Suggested change
pres_SL = 5e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL
Ms_L = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))
Ms_R = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))
pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL
Ms_L = max(1._wp, sqrt(1._wp + ((5.e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))
Ms_R = max(1._wp, sqrt(1._wp + ((5.e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5.e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))

Copilot uses AI. Check for mistakes.

Comment on lines +1803 to +1813
Copy link

Copilot AI Sep 17, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Inconsistent use of floating point literal notation. Line 1804, 1806, 1809, and 1814 use 5e-1_wp while other parts of the codebase use 5.e-1_wp (with decimal point). Should use 5.e-1_wp for consistency with established patterns in the codebase.

Suggested change
pres_SL = 5e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL
Ms_L = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))
Ms_R = max(1._wp, sqrt(1._wp + ((5e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))
pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(idx) - vel_R(idx)))
pres_SR = pres_SL
Ms_L = max(1._wp, sqrt(1._wp + ((5.e-1_wp + gamma_L)/(1._wp + gamma_L))* &
(pres_SL/pres_L - 1._wp)*pres_L/ &
((pres_L + pi_inf_L/(1._wp + gamma_L)))))
Ms_R = max(1._wp, sqrt(1._wp + ((5.e-1_wp + gamma_R)/(1._wp + gamma_R))* &
(pres_SR/pres_R - 1._wp)*pres_R/ &
((pres_R + pi_inf_R/(1._wp + gamma_R)))))
s_L = vel_L(idx) - c_L*Ms_L
s_R = vel_R(idx) + c_R*Ms_R
s_S = 5.e-1_wp*((vel_L(idx) + vel_R(idx)) + (pres_L - pres_R)/(rho_avg*c_avg))

Copilot uses AI. Check for mistakes.

end if

! ! follows Einfeldt et al.
Copy link

Copilot AI Sep 17, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Double comment marker ! ! is inconsistent with standard single comment style used elsewhere in the codebase. Should use single ! for consistency.

Suggested change
! ! follows Einfeldt et al.
! follows Einfeldt et al.

Copilot uses AI. Check for mistakes.

! s_M/P = min/max(0.,s_L/R)
s_M = min(0._wp, s_L)
s_P = max(0._wp, s_R)

Check warning on line 1805 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1802-L1805 

Added lines #L1802 - L1805 were not covered by tests

#ifdef DEBUG
! Check for potential issues in wave speed calculation
if (s_R <= s_L) then

Check warning on line 1809 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1809 

Added line #L1809 was not covered by tests
print *, 'WARNING: Wave speed issue detected in s_compute_wave_speed'
print *, 'Left wave speed >= Right wave speed:', s_L, s_R
print *, 'Input velocities :', vel_L(idx), vel_R(idx)

Check warning on line 1812 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1812 

Added line #L1812 was not covered by tests
print *, 'Sound speeds:', c_L, c_R
print *, 'Densities:', rho_L, rho_R
print *, 'Pressures:', pres_L, pres_R

Check warning on line 1815 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1814-L1815 

Added lines #L1814 - L1815 were not covered by tests
print *, 'Wave speeds method:', wave_speeds
if (elasticity .or. hypoelasticity .or. hyperelasticity) then

Check warning on line 1817 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1817 

Added line #L1817 was not covered by tests
print *, 'Shear moduli:', G_L, G_R
end if
call s_mpi_abort('Error: Invalid wave speeds in s_compute_wave_speed')

Check warning on line 1820 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1820 

Added line #L1820 was not covered by tests
end if
#endif

end subroutine s_compute_wave_speed

Check warning on line 1824 in src/common/m_variables_conversion.fpp

View check run for this annotation

Codecov / codecov/patch

src/common/m_variables_conversion.fpp#L1824 

Added line #L1824 was not covered by tests
#endif

end module m_variables_conversion
Loading
Loading