Redo Fortran source doc fixes without comment mangling

Previous commit used an overly broad regex that stripped @param,
@brief, and first words from comments across all .fpp files.
Restore all source files from master and re-apply only:
- Remove filename args from @file directives (58 files)
- Convert raw $...$ math to \f$...\f$ in m_riemann_solvers.fpp
- Add @cond/@endcond for #ifdef blocks with vendor attributes
- Fix stale @param name in m_data_input.f90

Co-Authored-By: Claude Opus 4.6 <noreply@anthropic.com>

Author: sbryngelson

src/common/include/1dHardcodedIC.fpp

@@ -1,37 +1,37 @@
 #:def Hardcoded1DVariables()
-    ! any declaration of intermediate variables here
+    ! Place any declaration of intermediate variables here
     real(wp) :: x_mid_diffu, width_sq, profile_shape, temp, molar_mass_inv, y1, y2, y3, y4
 #:enddef
 
 #:def Hardcoded1D()
     select case (patch_icpp(patch_id)%hcid)
     case (150) ! 1D Smooth Alfven Case for MHD
-        !
+        ! velocity
         q_prim_vf(momxb + 1)%sf(i, 0, 0) = 0.1_wp*sin(2._wp*pi*x_cc(i))
         q_prim_vf(momxb + 2)%sf(i, 0, 0) = 0.1_wp*cos(2._wp*pi*x_cc(i))
 
-        ! field
+        ! magnetic field
         q_prim_vf(B_idx%end - 1)%sf(i, 0, 0) = 0.1_wp*sin(2._wp*pi*x_cc(i))
         q_prim_vf(B_idx%end)%sf(i, 0, 0) = 0.1_wp*cos(2._wp*pi*x_cc(i))
 
     case (170)
-        ! hardcoded case can be used to start a simulation with initial conditions given from a known 1D profile (e.g. Cantera, SDtoolbox)
+        ! This hardcoded case can be used to start a simulation with initial conditions given from a known 1D profile (e.g. Cantera, SDtoolbox)
         @: HardcodedReadValues()
 
     case (180)
-        ! is patch is hard-coded for test suite optimization used in the
-        ! cases: "patch_icpp(2)%alpha_rho(1)": "1 + 0.2*sin(5*x)"
+        ! This is patch is hard-coded for test suite optimization used in the
+        ! 1D_shuoser cases: "patch_icpp(2)%alpha_rho(1)": "1 + 0.2*sin(5*x)"
         if (patch_id == 2) then
             q_prim_vf(contxb + 0)%sf(i, 0, 0) = 1 + 0.2*sin(5*x_cc(i))
         end if
 
     case (181)
-        ! is patch is hard-coded for test suite optimization used in the
-        ! cases: "patch_icpp(2)%alpha_rho(1)": "1 + 0.1*sin(20*x*pi)"
+        ! This is patch is hard-coded for test suite optimization used in the
+        ! 1D_titarevtorro cases: "patch_icpp(2)%alpha_rho(1)": "1 + 0.1*sin(20*x*pi)"
         q_prim_vf(contxb + 0)%sf(i, 0, 0) = 1 + 0.1*sin(20*x_cc(i)*pi)
 
     case (182)
-        ! patch is a hard-coded for test suite optimization (multiple component diffusion)
+        ! This patch is a hard-coded for test suite optimization (multiple component diffusion)
         x_mid_diffu = 0.05_wp/2.0_wp
         width_sq = (2.5_wp*10.0_wp**(-3.0_wp))**2
         profile_shape = 1.0_wp - 0.5_wp*exp(-(x_cc(i) - x_mid_diffu)**2/width_sq)

src/common/include/2dHardcodedIC.fpp

@@ -1,5 +1,5 @@
 #:def Hardcoded2DVariables()
-    ! any declaration of intermediate variables here
+    ! Place any declaration of intermediate variables here
     real(wp) :: eps, eps_mhd, C_mhd
     real(wp) :: r, rmax, gam, umax, p0
     real(wp) :: rhoH, rhoL, pRef, pInt, h, lam, wl, amp, intH, intL, alph
@@ -9,9 +9,9 @@
 
     real(wp) :: r_sq
 
-    ! 207
+    ! # 207
     real(wp) :: sigma, gauss1, gauss2
-    ! 208
+    ! # 208
     real(wp) :: ei, d, fsm, alpha_air, alpha_sf6
 
     eps = 1.e-9_wp
@@ -23,13 +23,13 @@
 
     case (200)
         if (y_cc(j) <= (-x_cc(i)**3 + 1)**(1._wp/3._wp)) then
-            ! Fractions
+            ! Volume Fractions
             q_prim_vf(advxb)%sf(i, j, 0) = eps
             q_prim_vf(advxe)%sf(i, j, 0) = 1._wp - eps
-            !
+            ! Denssities
             q_prim_vf(contxb)%sf(i, j, 0) = eps*1000._wp
             q_prim_vf(contxe)%sf(i, j, 0) = (1._wp - eps)*1._wp
-            !
+            ! Pressure
             q_prim_vf(E_idx)%sf(i, j, 0) = 1000._wp
         end if
     case (202) ! Gresho vortex (Gouasmi et al 2022 JCP)
@@ -149,7 +149,7 @@
         lam = 1.0_wp
         eps = 1.0e-6_wp
         ei = 5.0_wp
-        ! function to smooth out sharp discontinuity in the interface
+        ! Smoothening function to smooth out sharp discontinuity in the interface
         if (x_cc(i) <= 0.7_wp*lam) then
             d = x_cc(i) - lam*(0.4_wp - 0.1_wp*sin(2.0_wp*pi*(y_cc(j)/lam + 0.25_wp)))
             fsm = 0.5_wp*(1.0_wp + erf(d/(ei*sqrt(dx*dy))))
@@ -162,11 +162,11 @@
         end if
 
     case (250) ! MHD Orszag-Tang vortex
-        ! = 5/3
-        !   = 25/(36*pi)
-        !     = 5/(12*pi)
-        !     = (-sin(2*pi*y), sin(2*pi*x), 0)
-        !     = (-sin(2*pi*y)/sqrt(4*pi), sin(4*pi*x)/sqrt(4*pi), 0)
+        ! gamma = 5/3
+        !   rho = 25/(36*pi)
+        !     p = 5/(12*pi)
+        !     v = (-sin(2*pi*y), sin(2*pi*x), 0)
+        !     B = (-sin(2*pi*y)/sqrt(4*pi), sin(4*pi*x)/sqrt(4*pi), 0)
 
         q_prim_vf(momxb)%sf(i, j, 0) = -sin(2._wp*pi*y_cc(j))
         q_prim_vf(momxb + 1)%sf(i, j, 0) = sin(2._wp*pi*x_cc(i))
@@ -179,7 +179,7 @@
             q_prim_vf(contxb)%sf(i, j, 0) = 0.01
             q_prim_vf(E_idx)%sf(i, j, 0) = 1.0
         elseif (x_cc(i)**2 + y_cc(j)**2 <= 1._wp**2) then
-            ! interpolation between r=0.08 and r=1.0
+            ! Linear interpolation between r=0.08 and r=1.0
             factor = (1.0_wp - sqrt(x_cc(i)**2 + y_cc(j)**2))/(1.0_wp - 0.08_wp)
             q_prim_vf(contxb)%sf(i, j, 0) = 0.01_wp*factor + 1.e-4_wp*(1.0_wp - factor)
             q_prim_vf(E_idx)%sf(i, j, 0) = 1.0_wp*factor + 3.e-5_wp*(1.0_wp - factor)
@@ -188,69 +188,69 @@
             q_prim_vf(E_idx)%sf(i, j, 0) = 3.e-5_wp
         end if
 
-        ! 252 is for the 2D MHD Rotor problem
+        ! case 252 is for the 2D MHD Rotor problem
     case (252) ! 2D MHD Rotor Problem
-        ! conditions are set in the JSON file.
-        ! case imposes the dense, rotating cylinder.
+        ! Ambient conditions are set in the JSON file.
+        ! This case imposes the dense, rotating cylinder.
         !
-        ! = 1.4
-        ! medium (r > 0.1):
-        !   = 1, p = 1, v = 0, B = (1,0,0)
-        ! (r <= 0.1):
-        !   = 10, p = 1
-        !   has angular velocity w=20, giving v_tan=2 at r=0.1
-
-        ! distance squared from the center
+        ! gamma = 1.4
+        ! Ambient medium (r > 0.1):
+        !   rho = 1, p = 1, v = 0, B = (1,0,0)
+        ! Rotor (r <= 0.1):
+        !   rho = 10, p = 1
+        !   v has angular velocity w=20, giving v_tan=2 at r=0.1
+
+        ! Calculate distance squared from the center
         r_sq = (x_cc(i) - 0.5_wp)**2 + (y_cc(j) - 0.5_wp)**2
 
-        ! radius of 0.1
+        ! inner radius of 0.1
         if (r_sq <= 0.1**2) then
-            ! Inside the rotor --
-            ! density uniformly to 10
+            ! -- Inside the rotor --
+            ! Set density uniformly to 10
             q_prim_vf(contxb)%sf(i, j, 0) = 10._wp
 
-            ! vup constant rotation of rate v=2
-            ! = -omega * (y - y_c)
-            ! =  omega * (x - x_c)
+            ! Set vup constant rotation of rate v=2
+            ! v_x = -omega * (y - y_c)
+            ! v_y =  omega * (x - x_c)
             q_prim_vf(momxb)%sf(i, j, 0) = -20._wp*(y_cc(j) - 0.5_wp)
             q_prim_vf(momxb + 1)%sf(i, j, 0) = 20._wp*(x_cc(i) - 0.5_wp)
 
-            ! width of 0.015
+            ! taper width of 0.015
         else if (r_sq <= 0.115**2) then
-            ! smooth the function between r = 0.1 and 0.115
+            ! linearly smooth the function between r = 0.1 and 0.115
             q_prim_vf(contxb)%sf(i, j, 0) = 1._wp + 9._wp*(0.115_wp - sqrt(r_sq))/(0.015_wp)
 
             q_prim_vf(momxb)%sf(i, j, 0) = -(2._wp/sqrt(r_sq))*(y_cc(j) - 0.5_wp)*(0.115_wp - sqrt(r_sq))/(0.015_wp)
             q_prim_vf(momxb + 1)%sf(i, j, 0) = (2._wp/sqrt(r_sq))*(x_cc(i) - 0.5_wp)*(0.115_wp - sqrt(r_sq))/(0.015_wp)
         end if
 
     case (253) ! MHD Smooth Magnetic Vortex
-        ! 5.2 of
-        ! hybridized discontinuous Galerkin methods for compressible magnetohydrodynamics
-        ! Ciuca, P. Fernandez, A. Christophe, N.C. Nguyen, J. Peraire
+        ! Section 5.2 of
+        ! Implicit hybridized discontinuous Galerkin methods for compressible magnetohydrodynamics
+        ! C. Ciuca, P. Fernandez, A. Christophe, N.C. Nguyen, J. Peraire
 
-        !
+        ! velocity
         q_prim_vf(momxb)%sf(i, j, 0) = 1._wp - (y_cc(j)*exp(1 - (x_cc(i)**2 + y_cc(j)**2))/(2.*pi))
         q_prim_vf(momxb + 1)%sf(i, j, 0) = 1._wp + (x_cc(i)*exp(1 - (x_cc(i)**2 + y_cc(j)**2))/(2.*pi))
 
-        ! field
+        ! magnetic field
         q_prim_vf(B_idx%beg)%sf(i, j, 0) = -y_cc(j)*exp(1 - (x_cc(i)**2 + y_cc(j)**2))/(2.*pi)
         q_prim_vf(B_idx%beg + 1)%sf(i, j, 0) = x_cc(i)*exp(1 - (x_cc(i)**2 + y_cc(j)**2))/(2.*pi)
 
-        !
+        ! pressure
         q_prim_vf(E_idx)%sf(i, j, 0) = 1._wp + (1 - 2._wp*(x_cc(i)**2 + y_cc(j)**2))*exp(1 - (x_cc(i)**2 + y_cc(j)**2))/((2._wp*pi)**3)
 
     case (260)  ! Gaussian Divergence Pulse
-        !  = 1 + C * erf((x-0.5)/σ)
-        !   ∂Bx/∂x = C * (2/√π) * exp[-((x-0.5)/σ)**2] * (1/σ)
-        !  C = ε * σ * √π / 2  ⇒ ∂Bx/∂x = ε * exp[-((x-0.5)/σ)**2]
-        !  is initialized to zero everywhere.
+        !  Bx(x) = 1 + C * erf((x-0.5)/σ)
+        !  ⇒  ∂Bx/∂x = C * (2/√π) * exp[-((x-0.5)/σ)**2] * (1/σ)
+        !  Choose C = ε * σ * √π / 2  ⇒ ∂Bx/∂x = ε * exp[-((x-0.5)/σ)**2]
+        !  ψ is initialized to zero everywhere.
 
         eps_mhd = patch_icpp(patch_id)%a(2)
         sigma = patch_icpp(patch_id)%a(3)
         C_mhd = eps_mhd*sigma*sqrt(pi)*0.5_wp
 
-        !
+        ! B-field
         q_prim_vf(B_idx%beg)%sf(i, j, 0) = 1._wp + C_mhd*erf((x_cc(i) - 0.5_wp)/sigma)
 
     case (261)  ! Blob
@@ -260,21 +260,21 @@
         alpha = r/r0
         if (alpha < 1) then
             q_prim_vf(B_idx%beg)%sf(i, j, 0) = 1._wp/sqrt(4._wp*pi)*(alpha**8 - 2._wp*alpha**4 + 1._wp)
-            ! = 1._wp/sqrt(4000._wp*pi) * (4096._wp*r2**4 - 128._wp*r2**2 + 1._wp)
-            ! = 1._wp/(4._wp*pi) * (alpha**8 - 2._wp*alpha**4 + 1._wp)
-            ! = 6._wp - q_prim_vf(B_idx%beg)%sf(i,j,0)**2/2._wp
+            ! q_prim_vf(B_idx%beg)%sf(i,j,0) = 1._wp/sqrt(4000._wp*pi) * (4096._wp*r2**4 - 128._wp*r2**2 + 1._wp)
+            ! q_prim_vf(B_idx%beg)%sf(i,j,0) = 1._wp/(4._wp*pi) * (alpha**8 - 2._wp*alpha**4 + 1._wp)
+            ! q_prim_vf(E_idx)%sf(i,j,0) = 6._wp - q_prim_vf(B_idx%beg)%sf(i,j,0)**2/2._wp
         end if
 
     case (262)  ! Tilted 2D MHD shock‐tube at α = arctan2 (≈63.4°)
-        ! by α = atan(2)
+        ! rotate by α = atan(2)
         alpha = atan(2._wp)
         cosA = cos(alpha)
         sinA = sin(alpha)
-        ! along shock normal
+        ! projection along shock normal
         r = x_cc(i)*cosA + y_cc(j)*sinA
 
         if (r <= 0.5_wp) then
-            ! state: ρ=1, v∥=+10, v⊥=0, p=20, B∥=B⊥=5/√(4π)
+            ! LEFT state: ρ=1, v∥=+10, v⊥=0, p=20, B∥=B⊥=5/√(4π)
             q_prim_vf(contxb)%sf(i, j, 0) = 1._wp
             q_prim_vf(momxb)%sf(i, j, 0) = 10._wp*cosA
             q_prim_vf(momxb + 1)%sf(i, j, 0) = 10._wp*sinA
@@ -284,7 +284,7 @@
             q_prim_vf(B_idx%beg + 1)%sf(i, j, 0) = (5._wp/sqrt(4._wp*pi))*sinA &
                                                    + (5._wp/sqrt(4._wp*pi))*cosA
         else
-            ! state: ρ=1, v∥=−10, v⊥=0, p=1, B∥=B⊥=5/√(4π)
+            ! RIGHT state: ρ=1, v∥=−10, v⊥=0, p=1, B∥=B⊥=5/√(4π)
             q_prim_vf(contxb)%sf(i, j, 0) = 1._wp
             q_prim_vf(momxb)%sf(i, j, 0) = -10._wp*cosA
             q_prim_vf(momxb + 1)%sf(i, j, 0) = -10._wp*sinA
@@ -294,16 +294,16 @@
             q_prim_vf(B_idx%beg + 1)%sf(i, j, 0) = (5._wp/sqrt(4._wp*pi))*sinA &
                                                    + (5._wp/sqrt(4._wp*pi))*cosA
         end if
-        ! and B^z remain zero by default
+        ! v^z and B^z remain zero by default
 
     case (270)
-        ! hardcoded case extrudes a 1D profile to initialize a 2D simulation domain
+        ! This hardcoded case extrudes a 1D profile to initialize a 2D simulation domain
         @: HardcodedReadValues()
 
     case (280)
-        ! is patch is hard-coded for test suite optimization used in the
-        ! case:
-        ! analytic patch uses geometry 2
+        ! This is patch is hard-coded for test suite optimization used in the
+        ! 2D_isentropicvortex case:
+        ! This analytic patch uses geometry 2
         if (patch_id == 1) then
             q_prim_vf(E_idx)%sf(i, j, 0) = 1.0*(1.0 - (1.0/1.0)*(5.0/(2.0*pi))*(5.0/(8.0*1.0*(1.4 + 1.0)*pi))*exp(2.0*1.0*(1.0 - (x_cc(i) - patch_icpp(1)%x_centroid)**2.0 - (y_cc(j) - patch_icpp(1)%y_centroid)**2.0)))**(1.4 + 1.0)
             q_prim_vf(contxb + 0)%sf(i, j, 0) = 1.0*(1.0 - (1.0/1.0)*(5.0/(2.0*pi))*(5.0/(8.0*1.0*(1.4 + 1.0)*pi))*exp(2.0*1.0*(1.0 - (x_cc(i) - patch_icpp(1)%x_centroid)**2.0 - (y_cc(j) - patch_icpp(1)%y_centroid)**2.0)))**1.4
@@ -312,18 +312,18 @@
         end if
 
     case (281)
-        ! is patch is hard-coded for test suite optimization used in the
-        ! case:
-        ! analytic patch uses geometry 2
+        ! This is patch is hard-coded for test suite optimization used in the
+        ! 2D_acoustic_pulse case:
+        ! This analytic patch uses geometry 2
         if (patch_id == 2) then
             q_prim_vf(E_idx)%sf(i, j, 0) = 101325*(1 - 0.5*(1.4 - 1)*(0.4)**2*exp(0.5*(1 - sqrt(x_cc(i)**2 + y_cc(j)**2))))**(1.4/(1.4 - 1))
             q_prim_vf(contxb + 0)%sf(i, j, 0) = 1*(1 - 0.5*(1.4 - 1)*(0.4)**2*exp(0.5*(1 - sqrt(x_cc(i)**2 + y_cc(j)**2))))**(1/(1.4 - 1))
         end if
 
     case (282)
-        ! is patch is hard-coded for test suite optimization used in the
-        ! case:
-        ! analytic patch uses geometry 2
+        ! This is patch is hard-coded for test suite optimization used in the
+        ! 2D_zero_circ_vortex case:
+        ! This analytic patch uses geometry 2
         if (patch_id == 2) then
             q_prim_vf(E_idx)%sf(i, j, 0) = 101325*(1 - 0.5*(1.4 - 1)*(0.1/0.3)**2*exp(0.5*(1 - sqrt(x_cc(i)**2 + y_cc(j)**2))))**(1.4/(1.4 - 1))
             q_prim_vf(contxb + 0)%sf(i, j, 0) = 1*(1 - 0.5*(1.4 - 1)*(0.1/0.3)**2*exp(0.5*(1 - sqrt(x_cc(i)**2 + y_cc(j)**2))))**(1/(1.4 - 1))

src/common/include/3dHardcodedIC.fpp

@@ -1,12 +1,12 @@
 #:def Hardcoded3DVariables()
-    ! any declaration of intermediate variables here
+    ! Place any declaration of intermediate variables here
     real(wp) :: rhoH, rhoL, pRef, pInt, h, lam, wl, amp, intH, alph, Mach
     real(wp) :: eps
 
-    ! Jets
-    ! to stor position and radii of jets from input file
+    ! IGR Jets
+    ! Arrays to stor position and radii of jets from input file
     real(wp), dimension(:), allocatable :: y_th_arr, z_th_arr, r_th_arr
-    ! to describe initial condition of jet
+    ! Variables to describe initial condition of jet
     real(wp) :: r, ux_th, ux_am, p_th, p_am, rho_th, rho_am, y_th, z_th, r_th, eps_smooth
     real(wp) :: rcut, xcut ! Intermediate variables for creating smooth initial condition
 
@@ -175,13 +175,13 @@
         q_prim_vf(E_idx)%sf(i, j, k) = p_th*rcut*xcut + p_am
 
     case (370)
-        ! hardcoded case extrudes a 2D profile to initialize a 3D simulation domain
+        ! This hardcoded case extrudes a 2D profile to initialize a 3D simulation domain
         @: HardcodedReadValues()
 
     case (380)
-        ! is patch is hard-coded for test suite optimization used in the
-        ! case:
-        ! analytic patch used geometry 9
+        ! This is patch is hard-coded for test suite optimization used in the
+        ! 3D_TaylorGreenVortex case:
+        ! This analytic patch used geometry 9
         Mach = 0.1
         if (patch_id == 1) then
             q_prim_vf(E_idx)%sf(i, j, k) = 101325 + (Mach**2*376.636429464809**2/16)*(cos(2*x_cc(i)/1) + cos(2*y_cc(j)/1))*(cos(2*z_cc(k)/1) + 2)