m_eos.fpp.f90

Go to the documentation of this file.

# 1 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!>
!! @file
!! @brief Contains module m_eos
 
# 1 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 1
# 1 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 1
# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
! New line at end of file is required for FYPP
# 2 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
# 1 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 1
# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
! New line at end of file is required for FYPP
# 2 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp" 2
 
# 4 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
# 5 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
# 6 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
# 7 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
# 8 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 20 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 43 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 48 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 53 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 58 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 63 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 68 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 81 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 86 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 91 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 96 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 101 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 106 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 111 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 116 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 121 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 126 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 156 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 197 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 211 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 236 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 247 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 249 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
# 260 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 310 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 320 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 330 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 339 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 356 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 366 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 373 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 379 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 385 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 391 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 397 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
 
# 403 "/home/runner/work/MFC/MFC/src/common/include/omp_macros.fpp"
! New line at end of file is required for FYPP
# 3 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
# 1 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 1
# 1 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp" 1
# 2 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 3 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 4 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 5 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 6 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 8 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 9 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 10 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 15 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 16 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
# 17 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 24 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 53 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 65 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 75 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 105 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 117 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 127 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
 
# 174 "/home/runner/work/MFC/MFC/src/common/include/shared_parallel_macros.fpp"
! New line at end of file is required for FYPP
# 2 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp" 2
 
# 7 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 17 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 22 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 27 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 32 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 37 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 42 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 47 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 52 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 57 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 62 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 73 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 78 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 83 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 88 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 103 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 131 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 160 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 175 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 193 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 215 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 244 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 259 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 269 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 278 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 294 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 304 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
 
# 311 "/home/runner/work/MFC/MFC/src/common/include/acc_macros.fpp"
! New line at end of file is required for FYPP
# 4 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp" 2
 
! GPU parallel region (scalar reductions, maxval/minval)
# 23 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! GPU parallel loop over threads (most common GPU macro)
# 43 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Required closing for GPU_PARALLEL_LOOP
# 55 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Mark routine for device compilation
# 112 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Declare device-resident data
# 130 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Inner loop within a GPU parallel region
# 145 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Scoped GPU data region
# 164 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Host code with device pointers (for MPI with GPU buffers)
# 193 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Allocate device memory (unscoped)
# 207 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Free device memory
# 219 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Atomic operation on device
# 231 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! End atomic capture block
# 242 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Copy data between host and device
# 254 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Synchronization barrier
# 266 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Import GPU library module (openacc or omp_lib)
# 275 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Emit code only for AMD compiler
# 282 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Emit code for non-Cray compilers
# 289 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Emit code only for Cray compiler
# 296 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
! Emit code for non-NVIDIA compilers
# 303 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
 
# 305 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
# 306 "/home/runner/work/MFC/MFC/src/common/include/parallel_macros.fpp"
! New line at end of file is required for FYPP
# 2 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp" 2
 
# 14 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
 
! Caution: This macro requires the use of a binding script to set CUDA_VISIBLE_DEVICES, such that we have one GPU device per MPI
! rank. That's because for both cudaMemAdvise (preferred location) and cudaMemPrefetchAsync we use location = device_id = 0. For an
! example see misc/nvidia_uvm/bind.sh.
# 52 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
 
! Allocate and create GPU device memory
# 72 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
 
! Free GPU device memory and deallocate
# 80 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
 
! Cray-specific GPU pointer setup for vector fields
# 104 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
 
! Cray-specific GPU pointer setup for scalar fields
# 120 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
 
! Cray-specific GPU pointer setup for acoustic source spatials
# 145 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
 
# 151 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
 
# 158 "/home/runner/work/MFC/MFC/src/common/include/macros.fpp"
! New line at end of file is required for FYPP
# 6 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp" 2
# 1 "/home/runner/work/MFC/MFC/src/common/include/case.fpp" 1
! This file exists so that Fypp can be run without generating case.fpp files for
! each target. This is useful when generating documentation, for example. This
! should also let MFC be built with CMake directly, without invoking mfc.sh.
 
! For pre-process.
# 8 "/home/runner/work/MFC/MFC/src/common/include/case.fpp"
 
! For moving immersed boundaries in simulation
# 12 "/home/runner/work/MFC/MFC/src/common/include/case.fpp"
# 7 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp" 2
# 1 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp" 1
! AUTO-GENERATED - do not edit directly. Regenerate: cmake reconfigure
!
# 4 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
 
# 6 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
# 9 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
 
# 12 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
# 15 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
 
# 17 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
# 31 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
 
# 34 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
# 58 "/home/runner/work/MFC/MFC/build/include/simulation/generated_eos.fpp"
# 8 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp" 2
 
!> @brief Equations of state in Gamma/Pi form, rho e = Gamma(rho) p + Pi(rho).
!!
!! Stiffened and ideal gas keep constant coefficients, resolved once at start-up. The
!! state-dependent families (Mie-Gruneisen, JWL, Vinet) evaluate theirs per cell from a
!! reference curve.
!!
!! The mixture closure rules that combine the phases -- Wood's law, the six-equation mean, the
!! bubbly branch -- are not themselves equations of state, but s_compute_mixture_coefficients,
!! s_compute_speed_of_sound and their _dt/_avg variants live here rather than in
!! m_variables_conversion, and must stay here. They are the hot-path callers of the phase
!! chain (s_phase_coefficients -> s_eos_coefficients -> s_reference_curve), and on NVHPC that
!! inlining only happens within a single file: the cross-file inliner refuses any device routine
!! with a subroutine call in its call tree. Splitting them from the chain costs ~25% of grind
!! time on NVHPC and nothing on the other backends, so it fails quietly. The solver kernels never
!! inlined these four in the first place, which is why the module boundary is drawn above them and
!! not below. See docs/documentation/gpuParallelization.md, "Module boundaries and NVHPC inlining".
!!
!! This module is a leaf: it directly uses only m_derived_types, m_constants, and
!! m_global_parameters_common. Adding an EOS family means one case in s_reference_curve.
module m_eos
 
    use m_derived_types
    use m_global_parameters_common
    use m_mpi_common, only: s_prohibit_abort
    use m_constants, only: eos_ideal_gas, eos_mie_gruneisen, eos_jwl, eos_vinet, eos_rk4_steps, ode_isentrope, &
        & ode_reference_temperature, sgm_eps, verysmall, dflt_real
 
    implicit none
 
    private
 
    public :: s_compute_mixture_coefficients, s_compute_mixture_coefficients_dt, s_compute_speed_of_sound, &
        & s_compute_speed_of_sound_avg, s_initialize_eos_module, s_finalize_eos_module, f_pressure, f_bulk_modulus, &
        & f_relativistic_enthalpy, f_isentrope_exponent, f_isentrope_pressure, f_sg_thermal, f_mixture_temperature, &
        & f_is_state_dependent, s_phase_coefficients, s_phase_pressure_on_isentrope, s_phase_temperature, &
        & s_phase_density_on_isentrope, s_phase_internal_energy, s_phase_bulk_modulus
 
contains
 
    !> Resolve every fluid's EOS coefficients once, before any conversion runs.
    impure subroutine s_initialize_eos_module()
 
        integer :: i
        logical :: state_dependent  !< Whether this case's fluids need a density-dependent EOS
 
#ifdef MFC_DEBUG
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:54: ', '@:ALLOCATE(gammas (1:num_fluids))'
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    allocate (gammas(1:num_fluids))
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc enter data create(gammas) 
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target enter data map(always,alloc:gammas) 
# 54 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
#ifdef MFC_DEBUG
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:55: ', '@:ALLOCATE(eoss (1:num_fluids))'
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    allocate (eoss(1:num_fluids))
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc enter data create(eoss) 
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target enter data map(always,alloc:eoss) 
# 55 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
#ifdef MFC_DEBUG
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:56: ', '@:ALLOCATE(isentrope_n (1:num_fluids))'
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    allocate (isentrope_n(1:num_fluids))
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc enter data create(isentrope_n) 
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target enter data map(always,alloc:isentrope_n) 
# 56 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
#ifdef MFC_DEBUG
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:57: ', '@:ALLOCATE(pi_infs(1:num_fluids))'
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    allocate (pi_infs(1:num_fluids))
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc enter data create(pi_infs) 
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target enter data map(always,alloc:pi_infs) 
# 57 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
#ifdef MFC_DEBUG
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:58: ', '@:ALLOCATE(isentrope_B(1:num_fluids))'
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    allocate (isentrope_b(1:num_fluids))
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc enter data create(isentrope_B) 
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target enter data map(always,alloc:isentrope_B) 
# 58 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
#ifdef MFC_DEBUG
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:59: ', '@:ALLOCATE(cvs (1:num_fluids))'
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    allocate (cvs(1:num_fluids))
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc enter data create(cvs) 
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target enter data map(always,alloc:cvs) 
# 59 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
#ifdef MFC_DEBUG
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:60: ', '@:ALLOCATE(qvs (1:num_fluids))'
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    allocate (qvs(1:num_fluids))
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc enter data create(qvs) 
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target enter data map(always,alloc:qvs) 
# 60 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
#ifdef MFC_DEBUG
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:61: ', '@:ALLOCATE(qvps (1:num_fluids))'
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    allocate (qvps(1:num_fluids))
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc enter data create(qvps) 
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target enter data map(always,alloc:qvps) 
# 61 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        state_dependent = .false.
        do i = 1, num_fluids
            gammas(i) = fluid_pp(i)%gamma
            isentrope_n(i) = f_isentrope_exponent(gammas(i))
 
            ! Each EOS supplies its own coefficients. Resolved once here, not per cell: a branch in the mixture loop costs
            ! registers in the Riemann kernels. An EOS whose coefficients depend on state must move to per-cell evaluation.
            select case (fluid_pp(i)%eos)
            case (eos_ideal_gas)
                pi_infs(i) = 0._wp
            case default
                pi_infs(i) = fluid_pp(i)%pi_inf
            end select
            isentrope_b(i) = f_isentrope_pressure(pi_infs(i), gammas(i))
            cvs(i) = fluid_pp(i)%cv
            qvs(i) = fluid_pp(i)%qv
            qvps(i) = fluid_pp(i)%qvp
            eoss(i) = fluid_pp(i)%eos
            ! Every fluid's single-source coefficients, mu_max among them: where a cubic Hugoniot fit turns over.
            ! mu(u_p) peaks where c0 = s2 u_p^2 + 2 s3 u_p^3, and past it no shock state exists, so the Newton below
            ! would wander. Solved once here, on the host.
    eos_coeffs(i)%c0 = fluid_pp(i)%mg_c0
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%s = fluid_pp(i)%mg_s
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%s2 = fluid_pp(i)%mg_s2
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%s3 = fluid_pp(i)%mg_s3
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%mu_max = f_hugoniot_compression_limit(fluid_pp(i)%mg_c0, fluid_pp(i)%mg_s, fluid_pp(i)%mg_s2, &
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        & fluid_pp(i)%mg_s3)
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%a = fluid_pp(i)%jwl_a
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%b = fluid_pp(i)%jwl_b
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%r1 = fluid_pp(i)%jwl_r1
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%r2 = fluid_pp(i)%jwl_r2
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%k0 = fluid_pp(i)%vinet_k0
# 84 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    eos_coeffs(i)%k0p = fluid_pp(i)%vinet_k0p
            ! One reference state and Gruneisen closure for every family; the user-facing names keep their prefix.
    select case (fluid_pp(i)%eos)
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    case (eos_mie_gruneisen)
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%rho0 = fluid_pp(i)%mg_rho0
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%t0 = fluid_pp(i)%mg_t0
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%gruneisen0 = fluid_pp(i)%mg_gruneisen
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%gruneisen_a = fluid_pp(i)%mg_gruneisen_a
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    case (eos_jwl)
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%rho0 = fluid_pp(i)%jwl_rho0
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%t0 = fluid_pp(i)%jwl_t0
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%gruneisen0 = fluid_pp(i)%jwl_omega
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%gruneisen_a = 0._wp
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    case (eos_vinet)
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%rho0 = fluid_pp(i)%vinet_rho0
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%t0 = fluid_pp(i)%vinet_t0
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%gruneisen0 = fluid_pp(i)%vinet_gruneisen
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%gruneisen_a = fluid_pp(i)%vinet_gruneisen_a
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    case default
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%rho0 = dflt_real
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%t0 = dflt_real
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%gruneisen0 = dflt_real
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        eos_coeffs(i)%gruneisen_a = 0._wp
# 86 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end select
            if (f_is_state_dependent(i)) state_dependent = .true.
        end do
# 95 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
            any_state_dependent_eos = state_dependent
# 97 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 97 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 97 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc update device(gammas, isentrope_n, pi_infs, isentrope_B, cvs, qvs, qvps, eoss, eos_coeffs) 
# 97 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 97 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target update to(gammas, isentrope_n, pi_infs, isentrope_B, cvs, qvs, qvps, eoss, eos_coeffs) 
# 97 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 99 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 99 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 99 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc update device(any_state_dependent_eos) 
# 99 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 99 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target update to(any_state_dependent_eos) 
# 99 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 101 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
    end subroutine s_initialize_eos_module
 
    !> Deallocate the fluid property arrays allocated in s_initialize_eos_module.
    impure subroutine s_finalize_eos_module()
 
#ifdef MFC_DEBUG
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    block
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        use iso_fortran_env, only: output_unit
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        print *, 'm_eos.fpp:107: ', '@:DEALLOCATE(gammas, isentrope_n, pi_infs, isentrope_B, cvs, qvs, qvps, eoss)'
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        call flush (output_unit)
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    end block
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc exit data delete(gammas, isentrope_n, pi_infs, isentrope_B, cvs, qvs, qvps, eoss) 
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp target exit data map(release:gammas, isentrope_n, pi_infs, isentrope_B, cvs, qvs, qvps, eoss) 
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 107 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
    deallocate (gammas, isentrope_n, pi_infs, isentrope_b, cvs, qvs, qvps, eoss)
 
    end subroutine s_finalize_eos_module
 
    !> The reference curve of a state-dependent EOS at rho: p_ref, e_ref, their d/drho, and Gamma_G with its d/drho. A new family
    !! adds one case here and nothing else.
    subroutine s_reference_curve(rho, i, p_ref, e_ref, dp_drho, de_drho, G0, dG0)
 
 
# 115 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 115 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 115 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 115 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 115 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 115 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 115 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in)  :: rho
        integer, intent(in)   :: i
        real(wp), intent(out) :: p_ref, e_ref, dp_drho, de_drho, G0, dG0
        real(wp)              :: mu, d, V, ea, eb, up, us, dus, dup_dmu, x, ex, dp_dmu, de_dmu
        integer               :: iter
 
        mu = rho/eos_coeffs(i)%rho0 - 1._wp
        ! Past the fit's turnover there is no shock state to find; clamp rather than let the Newton below wander
        ! off and return a silently wrong pressure. mu_max is huge for the linear fit, so this is a no-op there.
        ! Bounded here so the step finishes and the host-side check in s_write_run_time_information can report it;
        ! past the turnover there is no shock state and the Newton below would wander.
        if (eoss(i) == eos_mie_gruneisen .and. mu > eos_coeffs(i)%mu_max) mu = eos_coeffs(i)%mu_max
        select case (eoss(i))
        case (eos_mie_gruneisen)
            ! Hugoniot reference u_s = c0 + s u_p + s2 u_p^2 + s3 u_p^3, with p_H = rho0 u_s u_p and the Hugoniot
            ! energy e_H = p_H mu/(2 rho0 (1 + mu)); linear on release. Pole at mu = 1/(s - 1) for the linear fit;
            ! the validator refuses initial states outside the EOS.
            if (mu < 0._wp) then
                p_ref = eos_coeffs(i)%rho0*eos_coeffs(i)%c0**2*mu
                dp_dmu = eos_coeffs(i)%rho0*eos_coeffs(i)%c0**2
            else if (eos_coeffs(i)%s2 == 0._wp .and. eos_coeffs(i)%s3 == 0._wp) then
                d = 1._wp - (eos_coeffs(i)%s - 1._wp)*mu
                p_ref = eos_coeffs(i)%rho0*eos_coeffs(i)%c0**2*mu*(1._wp + mu)/(d*d)
                dp_dmu = eos_coeffs(i)%rho0*eos_coeffs(i)%c0**2*((1._wp + 2._wp*mu)*d + 2._wp*(eos_coeffs(i)%s - 1._wp)*mu*(1._wp &
                                    & + mu))/(d*d*d)
            else
                ! u_p solves u_s(u_p) mu = u_p (1 + mu): Newton from the linear fit, then implicit differentiation
                up = eos_coeffs(i)%c0*mu/(1._wp - (eos_coeffs(i)%s - 1._wp)*mu)
 
# 145 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 145 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 145 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 145 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 145 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
                do iter = 1, 8
                    us = eos_coeffs(i)%c0 + up*(eos_coeffs(i)%s + up*(eos_coeffs(i)%s2 + up*eos_coeffs(i)%s3))
                    dus = eos_coeffs(i)%s + up*(2._wp*eos_coeffs(i)%s2 + 3._wp*eos_coeffs(i)%s3*up)
                    up = up - (us*mu - up*(1._wp + mu))/(dus*mu - (1._wp + mu))
                end do
                us = eos_coeffs(i)%c0 + up*(eos_coeffs(i)%s + up*(eos_coeffs(i)%s2 + up*eos_coeffs(i)%s3))
                dus = eos_coeffs(i)%s + up*(2._wp*eos_coeffs(i)%s2 + 3._wp*eos_coeffs(i)%s3*up)
                dup_dmu = (up - us)/(dus*mu - (1._wp + mu))
                p_ref = eos_coeffs(i)%rho0*us*up
                dp_dmu = eos_coeffs(i)%rho0*(dus*up + us)*dup_dmu
            end if
            e_ref = p_ref*mu/(2._wp*eos_coeffs(i)%rho0*(1._wp + mu))
            de_dmu = (dp_dmu*mu*(1._wp + mu) + p_ref)/(2._wp*eos_coeffs(i)%rho0*(1._wp + mu)**2)
            dp_drho = dp_dmu/eos_coeffs(i)%rho0
            de_drho = de_dmu/eos_coeffs(i)%rho0
        case (eos_jwl)
            ! JWL: p_ref = A exp(-R1 V) + B exp(-R2 V), V = rho0/rho. The curve is itself an isentrope, so de_ref = -p_ref d(1/rho).
            v = eos_coeffs(i)%rho0/rho
            ea = eos_coeffs(i)%a*exp(-eos_coeffs(i)%r1*v)
            eb = eos_coeffs(i)%b*exp(-eos_coeffs(i)%r2*v)
            p_ref = ea + eb
            e_ref = (ea/eos_coeffs(i)%r1 + eb/eos_coeffs(i)%r2)/eos_coeffs(i)%rho0
            dp_drho = (eos_coeffs(i)%rho0/rho**2)*(eos_coeffs(i)%r1*ea + eos_coeffs(i)%r2*eb)
            de_drho = p_ref/rho**2
        case (eos_vinet)
            ! Vinet cold curve: p_c = 3 K0 (1 - x)/x^2 exp(eta (1 - x)), x = (rho0/rho)^(1/3), eta = 3 (K0' - 1)/2,
            ! an isentrope like JWL (its energy integrates in closed form).
            d = 1.5_wp*(eos_coeffs(i)%k0p - 1._wp)
            x = (eos_coeffs(i)%rho0/rho)**(1._wp/3._wp)
            ex = exp(d*(1._wp - x))
            p_ref = 3._wp*eos_coeffs(i)%k0*(1._wp - x)/x**2*ex
            e_ref = 9._wp*eos_coeffs(i)%k0/(eos_coeffs(i)%rho0*d**2)*(1._wp - (1._wp - d*(1._wp - x))*ex)
            dp_drho = 3._wp*eos_coeffs(i)%k0*ex*(-1._wp/x**2 - 2._wp*(1._wp - x)/x**3 - d*(1._wp - x)/x**2)*(-x/(3._wp*rho))
            de_drho = p_ref/rho**2
        end select
        g0 = eos_coeffs(i)%gruneisen0 + eos_coeffs(i)%gruneisen_a*mu
        dg0 = eos_coeffs(i)%gruneisen_a/eos_coeffs(i)%rho0
 
    end subroutine s_reference_curve
 
    !> Whether the EOS of fluid i is a family whose coefficients vary with density.
    function f_is_state_dependent(i) result(yes)
 
 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_is_state_dependent
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 189 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        integer, intent(in) :: i
        logical             :: yes
 
        yes = eoss(i) == eos_mie_gruneisen .or. eoss(i) == eos_jwl .or. eoss(i) == eos_vinet
 
    end function f_is_state_dependent
 
    !> True when fluid i's reference curve is itself an isentrope (de_ref = -p_ref d(1/rho), which holds for JWL and Vinet but not
    !! for the Mie-Gruneisen Hugoniot) and its Gruneisen coefficient is constant. Those two together make the isentrope through any
    !! state closed-form, so it never has to be integrated.
    function f_has_isentropic_reference(i) result(yes)
 
 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_has_isentropic_reference
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 203 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        integer, intent(in) :: i
        logical             :: yes
 
        yes = (eoss(i) == eos_jwl .or. eoss(i) == eos_vinet) .and. eos_coeffs(i)%gruneisen_a == 0._wp
 
    end function f_has_isentropic_reference
 
    !> The largest compression a cubic Hugoniot fit can represent. mu(u_p) = u_p/(u_s - u_p) rises, peaks where c0 = s2 u_p^2 + 2 s3
    !! u_p^3, and falls after; only the rising branch is a physical shock. Returns a huge value for the linear fit, which never
    !! turns over. Host-side: called once per fluid at initialization.
    impure function f_hugoniot_compression_limit(c0, s, s2, s3) result(mu_max)
 
        real(wp), intent(in) :: c0, s, s2, s3
        real(wp)             :: mu_max, up, f, df, us
        integer              :: iter
 
        if (s2 == 0._wp .and. s3 == 0._wp) then
            mu_max = huge(1._wp)
            return
        end if
 
        ! Newton on c0 - s2 u^2 - 2 s3 u^3 = 0, from a guess that brackets the physical range
        up = c0
        do iter = 1, 100
            f = c0 - s2*up**2 - 2._wp*s3*up**3
            df = -2._wp*s2*up - 6._wp*s3*up**2
            if (abs(df) < verysmall) exit
            up = max(up - f/df, verysmall)
        end do
        us = c0 + up*(s + up*(s2 + up*s3))
        mu_max = up/max(us - up, verysmall)
 
    end function f_hugoniot_compression_limit
 
    !> Gamma, Pi, dPi/drho and dGamma/drho of fluid i at density rho, the coefficients of rho e = Gamma p + Pi(rho). Stiffened and
    !! ideal gas keep the constants resolved at init, bit for bit.
    subroutine s_eos_coefficients(rho, i, gamma, pi_inf, dpi, dgamma)
 
 
# 243 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 243 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 243 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 243 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 243 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 243 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 243 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in)  :: rho
        integer, intent(in)   :: i
        real(wp), intent(out) :: gamma, pi_inf, dpi, dgamma
        real(wp)              :: p_ref, e_ref, dp_drho, de_drho, G0, dG0
 
        if (.not. f_is_state_dependent(i)) then
            gamma = gammas(i)
            pi_inf = pi_infs(i)
            dpi = 0._wp
            dgamma = 0._wp
            return
        end if
        call s_reference_curve(rho, i, p_ref, e_ref, dp_drho, de_drho, g0, dg0)
        gamma = 1._wp/g0
        pi_inf = rho*e_ref - p_ref/g0
        dpi = e_ref + rho*de_drho - dp_drho/g0 + p_ref*dg0/g0**2
        dgamma = -dg0/g0**2
 
    end subroutine s_eos_coefficients
 
    !> Exponent of the stiffened-gas isentrope p + B = const rho**n. Precomputed per fluid as isentrope_n.
    function f_isentrope_exponent(gamma) result(n)
 
 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_isentrope_exponent
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 268 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: gamma
        real(wp)             :: n
 
        n = 1._wp/gamma + 1._wp
 
    end function f_isentrope_exponent
 
    !> Reference pressure of that isentrope. Precomputed per fluid as isentrope_B.
    function f_isentrope_pressure(pi_inf, gamma) result(B)
 
 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_isentrope_pressure
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 280 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: pi_inf, gamma
        real(wp)             :: b
 
        b = pi_inf/(1._wp + gamma)
 
    end function f_isentrope_pressure
 
    !> Stiffened-gas thermal law p + B = (n - 1)*cv*rho*T. Pass rho to get T, or T to get rho.
    function f_sg_thermal(pres, rho_or_T, n, B, cv) result(T_or_rho)
 
 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_sg_thermal
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 292 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: pres, rho_or_t, n, b, cv
        real(wp)             :: t_or_rho
 
        t_or_rho = (pres + b)/((n - 1._wp)*cv*rho_or_t)
 
    end function f_sg_thermal
 
    !> Thermal-equilibrium mixture temperature for stiffened gas, from primitives. Algebraically identical to the conservative form
    !! in m_phase_change's s_infinite_pt_relaxation_k, T = (rho*e + p - sum(alpha_rho_i*qv_i)) / sum(alpha_rho_i*cv_i*n_i), because
    !! rho*e = gamma_mix*p + pi_inf_mix + sum(alpha_rho_i*qv_i) in MFC's stored variables.
    function f_mixture_temperature(alpha_rho_K, pres, gamma_K, pi_inf_K) result(T)
 
 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_mixture_temperature
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 306 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
# 311 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
            real(wp), dimension(num_fluids), intent(in) :: alpha_rho_k
# 313 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        real(wp), intent(in) :: pres, gamma_k, pi_inf_k
        real(wp)             :: t
        real(wp)             :: mcp  !< sum of alpha_rho_i*cp_i; cp_i = n_i*cv_i
        integer              :: i
 
        mcp = 0._wp
 
# 319 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 319 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 319 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 319 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 319 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
        do i = 1, num_fluids
            mcp = mcp + alpha_rho_k(i)*cvs(i)*isentrope_n(i)
        end do
 
        t = ((gamma_k + 1._wp)*pres + pi_inf_k)/max(mcp, sgm_eps)
 
    end function f_mixture_temperature
 
    !> Coefficients of phase i at its own density alpha_rho/alpha: the per-cell dispatch when some fluid's EOS is state dependent,
    !! the constants resolved at init otherwise (bit for bit).
    subroutine s_phase_coefficients(alpha_rho, alpha, i, rho, gamma, pi_inf, dpi, dgamma)
 
 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS s_phase_coefficients
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 332 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in)  :: alpha_rho, alpha
        integer, intent(in)   :: i
        real(wp), intent(out) :: rho, gamma, pi_inf, dpi, dgamma
 
        rho = max(alpha_rho, sgm_eps)/max(alpha, sgm_eps)
        if (any_state_dependent_eos) then
            call s_eos_coefficients(rho, i, gamma, pi_inf, dpi, dgamma)
        else
            gamma = gammas(i)
            pi_inf = pi_infs(i)
            dpi = 0._wp
            dgamma = 0._wp
        end if
 
    end subroutine s_phase_coefficients
 
    !> c^2 = [((Gamma + 1) p + Pi)/rho - dPi/drho - p dGamma/drho]/Gamma, the frozen speed of one phase.
    function f_c2_from_coefficients(rho, pres, gamma, pi_inf, dpi, dgamma) result(c2)
 
 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_c2_from_coefficients
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 353 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: rho, pres, gamma, pi_inf, dpi, dgamma
        real(wp)             :: c2
 
        c2 = (((gamma + 1._wp)*pres + pi_inf)/rho - dpi - pres*dgamma)/gamma
 
    end function f_c2_from_coefficients
 
    !> Frozen sound speed squared of one phase at (rho, p) from its own coefficients. These helpers are subroutines, not functions:
    !! a device function that calls a device subroutine is a pattern no other backend-tested code in MFC uses.
    subroutine s_phase_c2(rho, pres, i, c2)
 
 
# 366 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 366 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 366 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 366 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 366 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 366 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 366 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in)  :: rho, pres
        integer, intent(in)   :: i
        real(wp), intent(out) :: c2
        real(wp)              :: gamma, pi_inf, dpi, dgamma
 
        call s_eos_coefficients(rho, i, gamma, pi_inf, dpi, dgamma)
        c2 = f_c2_from_coefficients(rho, pres, gamma, pi_inf, dpi, dgamma)
 
    end subroutine s_phase_c2
 
    !> Slope of the ODE `kind` for fluid i: dp/drho = c^2 along an isentrope (x = rho, y = p), or the reference temperature dT/dV =
    !! (de_ref/dV + p_ref)/c_v - Gamma_G T/V (x = V, y = T), the Maxwell relation applied to e = e_ref + c_v (T - T_ref).
    subroutine s_ode_slope(kind, i, x, y, dydx)
 
 
# 382 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 382 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 382 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 382 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 382 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 382 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 382 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        integer, intent(in)   :: kind, i
        real(wp), intent(in)  :: x, y
        real(wp), intent(out) :: dydx
        real(wp)              :: p_ref, e_ref, dp_drho, de_drho, G0, dG0
 
        if (kind == ode_isentrope) then
            call s_phase_c2(x, y, i, dydx)
        else
            call s_reference_curve(1._wp/x, i, p_ref, e_ref, dp_drho, de_drho, g0, dg0)
            dydx = (p_ref - de_drho/x**2)/cvs(i) - g0*y/x
        end if
 
    end subroutine s_ode_slope
 
    !> Fixed-step classical RK4 for the ODE `kind` from (x0, y0) to x1.
    subroutine s_rk4(kind, i, x0, y0, x1, y)
 
 
# 401 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 401 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 401 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 401 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 401 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 401 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 401 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        integer, intent(in)   :: kind, i
        real(wp), intent(in)  :: x0, y0, x1
        real(wp), intent(out) :: y
        real(wp)              :: x, h, k1, k2, k3, k4
        integer               :: step
 
        x = x0
        y = y0
        h = (x1 - x0)/eos_rk4_steps
 
# 412 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 412 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 412 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 412 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 412 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
        do step = 1, eos_rk4_steps
            call s_ode_slope(kind, i, x, y, k1)
            call s_ode_slope(kind, i, x + 0.5_wp*h, y + 0.5_wp*h*k1, k2)
            call s_ode_slope(kind, i, x + 0.5_wp*h, y + 0.5_wp*h*k2, k3)
            call s_ode_slope(kind, i, x + h, y + h*k3, k4)
            y = y + h*(k1 + 2._wp*(k2 + k3) + k4)/6._wp
            x = x + h
        end do
 
    end subroutine s_rk4
 
    !> Pressure of phase i after the isentropic density change rho -> xi rho: closed form for the constant-coefficient families,
    !! integrated for a state-dependent EOS (the star states it serves are close to rho).
    subroutine s_phase_pressure_on_isentrope(pres, rho, xi, i, p_isen)
 
 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS s_phase_pressure_on_isentrope
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 428 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in)  :: pres, rho, xi
        integer, intent(in)   :: i
        real(wp), intent(out) :: p_isen
        real(wp)              :: p_ref_from, p_ref_to, e_ref, dp_drho, de_drho, g0, dg0
 
        if (.not. f_is_state_dependent(i)) then
            p_isen = (pres + isentrope_b(i))*xi**isentrope_n(i) - isentrope_b(i)
        else if (f_has_isentropic_reference(i)) then
            ! Exact: the offset from an isentropic reference obeys dDelta/Delta = Gamma drho/rho, so
            ! p - p_ref scales as (rho'/rho)**(1 + Gamma). Integrating it instead costs a decimal per
            ! doubling of the expansion and turns the pressure negative past roughly twentyfold.
            call s_reference_curve(rho, i, p_ref_from, e_ref, dp_drho, de_drho, g0, dg0)
            call s_reference_curve(xi*rho, i, p_ref_to, e_ref, dp_drho, de_drho, g0, dg0)
            p_isen = p_ref_to + (pres - p_ref_from)*xi**(1._wp + eos_coeffs(i)%gruneisen0)
        else
            call s_rk4(ode_isentrope, i, rho, pres, xi*rho, p_isen)
        end if
 
    end subroutine s_phase_pressure_on_isentrope
 
    !> Temperature of phase i at (rho, p): the stiffened-gas relation, or T_ref(rho) + (e - e_ref)/c_v.
    subroutine s_phase_temperature(rho, pres, i, T)
 
 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS s_phase_temperature
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 453 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in)  :: rho, pres
        integer, intent(in)   :: i
        real(wp), intent(out) :: t
        real(wp)              :: p_ref, e_ref, dp_drho, de_drho, g0, dg0, t0, t_ref
 
        if (f_is_state_dependent(i)) then
            call s_reference_curve(rho, i, p_ref, e_ref, dp_drho, de_drho, g0, dg0)
            t0 = eos_coeffs(i)%t0
            call s_rk4(ode_reference_temperature, i, 1._wp/eos_coeffs(i)%rho0, t0, 1._wp/rho, t_ref)
            t = t_ref + (pres - p_ref)/(rho*g0*cvs(i))
        else
            t = (pres + isentrope_b(i))/((isentrope_n(i) - 1._wp)*cvs(i)*rho)
        end if
 
    end subroutine s_phase_temperature
 
    !> Density of phase i on the isentrope through (rho_from, p_from) at p_to, and c^2 there: Newton on the pressure integrator,
    !! whose slope is c^2. The relaxation's own Newton wraps this, so a few steps suffice.
    subroutine s_phase_density_on_isentrope(i, rho_from, p_from, p_to, rho_to, c2_to)
 
 
# 475 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 475 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 475 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 475 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 475 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 475 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 475 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        integer, intent(in)   :: i
        real(wp), intent(in)  :: rho_from, p_from, p_to
        real(wp), intent(out) :: rho_to, c2_to
        real(wp)              :: p_at, c2_at
        integer               :: iter
 
        rho_to = rho_from
 
# 484 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 484 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 484 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 484 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 484 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
        do iter = 1, 4
            call s_phase_pressure_on_isentrope(p_from, rho_from, rho_to/rho_from, i, p_at)
            call s_phase_c2(rho_to, p_at, i, c2_at)
            rho_to = rho_to - (p_at - p_to)/c2_at
        end do
        call s_phase_c2(rho_to, p_to, i, c2_to)
 
    end subroutine s_phase_density_on_isentrope
 
    !> Internal energy per unit volume of phase i at pressure pres: alpha (Gamma p + Pi) + alpha_rho qv, with the coefficients at
    !! the phase's own density.
    subroutine s_phase_internal_energy(pres, alpha, alpha_rho, i, e_phase)
 
 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS s_phase_internal_energy
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 498 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in)  :: pres, alpha, alpha_rho
        integer, intent(in)   :: i
        real(wp), intent(out) :: e_phase
        real(wp)              :: rho, gamma, pi_inf, dpi, dgamma
 
        call s_phase_coefficients(alpha_rho, alpha, i, rho, gamma, pi_inf, dpi, dgamma)
        e_phase = alpha*(gamma*pres + pi_inf) + alpha_rho*qvs(i)
 
    end subroutine s_phase_internal_energy
 
    !> Bulk modulus rho c^2 of phase i at pressure pres: f_bulk_modulus for a constant-coefficient fluid, bit for bit, minus the
    !! reference-curve terms rho (dPi/drho + p dGamma/drho)/Gamma otherwise.
    subroutine s_phase_bulk_modulus(pres, alpha, alpha_rho, i, blkmod)
 
 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS s_phase_bulk_modulus
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 514 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in)  :: alpha_rho, alpha, pres
        integer, intent(in)   :: i
        real(wp), intent(out) :: blkmod
        real(wp)              :: rho, gamma, pi_inf, dpi, dgamma
 
        call s_phase_coefficients(alpha_rho, alpha, i, rho, gamma, pi_inf, dpi, dgamma)
        blkmod = f_bulk_modulus(pres, gamma, pi_inf) - rho*(dpi + pres*dgamma)/gamma
 
    end subroutine s_phase_bulk_modulus
 
    !> Pressure of a stiffened gas from its internal energy density - the inverse of s_compute_energy. Callers subtract the kinetic,
    !! magnetic and elastic energy first; none of those are equation-of-state terms.
    function f_pressure(e_int, gamma, pi_inf, qv) result(pres)
 
 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_pressure
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 530 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: e_int, gamma, pi_inf, qv
        real(wp)             :: pres
 
        pres = (e_int - pi_inf - qv)/gamma
 
    end function f_pressure
 
    !> Isentropic bulk modulus. Takes coefficients rather than a fluid index, so a mixture - whose effective gamma and pi_inf come
    !! from s_compute_mixture_coefficients - is the same call as a single fluid. Elastic callers add their own shear term.
    function f_bulk_modulus(pres, gamma, pi_inf) result(blkmod)
 
 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_bulk_modulus
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 543 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: pres, gamma, pi_inf
        real(wp)             :: blkmod
 
        blkmod = ((gamma + 1._wp)*pres + pi_inf)/gamma
 
    end function f_bulk_modulus
 
    !> Relativistic specific enthalpy, h = 1 + (Gamma + 1)p/rho. Ideal gas only: the stiffness does not appear, so a fluid with a
    !! nonzero pi_inf is not represented here (the validator refuses that combination).
    function f_relativistic_enthalpy(pres, rho, gamma) result(H)
 
 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS f_relativistic_enthalpy
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 556 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: pres, rho, gamma
        real(wp)             :: h
 
        h = 1._wp + (gamma + 1._wp)*pres/rho
 
    end function f_relativistic_enthalpy
 
    !> Mixture coefficients of one state. Under bubbles_euler with num_fluids == 1 the sole advection slot aliases the void fraction
    !! (eqn_idx%alf == eqn_idx%adv%end), so alpha is not a composition there and the coefficients are the liquid's. Clipping stays
    !! with callers; it differs between solvers and cannot coincide with that case, as mpp_lim requires num_fluids > 1.
    subroutine s_compute_mixture_coefficients(alpha_rho_K, alpha_K, rho_K, gamma_K, pi_inf_K, qv_K)
 
 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS s_compute_mixture_coefficients
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 570 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
# 575 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
            real(wp), dimension(num_fluids), intent(in) :: alpha_rho_k, alpha_k
# 577 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        real(wp), intent(out) :: rho_k, gamma_k, pi_inf_k, qv_k
        real(wp)              :: gamma_i, pi_inf_i, dpi_i, dgamma_i
        real(wp)              :: rho_i, alpha_i, alpha_rho_i
        integer               :: i  !< Loop iterator over fluids
 
        ! The bubbly closure is written for one carrier liquid, which keeps its own coefficients
        ! undiluted: Gamma_l*p_l = (E - rho|u|^2/2)/(1 - alf) - Pi_inf_l, the void entering only through
        ! the (1 - alf) that s_compute_pressure applies. There is nothing to sum - the last advection
        ! slot is the void, not a material - and the checker holds num_fluids <= 2 here.
        if (bubbles_euler) then
            rho_k = alpha_rho_k(1)
            gamma_k = gammas(1)
            pi_inf_k = pi_infs(1)
            ! Energy per unit volume, as below: alpha_rho_K(1) is the liquid partial density
            qv_k = alpha_rho_k(1)*qvs(1)
        else
            rho_k = 0._wp
            gamma_k = 0._wp
            pi_inf_k = 0._wp
            qv_k = 0._wp
 
 
# 598 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 598 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 598 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 598 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 598 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
            do i = 1, num_fluids
                rho_k = rho_k + alpha_rho_k(i)
                alpha_rho_i = alpha_rho_k(i)
                alpha_i = alpha_k(i)
                call s_phase_coefficients(alpha_rho_i, alpha_i, i, rho_i, gamma_i, pi_inf_i, dpi_i, dgamma_i)
                gamma_k = gamma_k + alpha_k(i)*gamma_i
                pi_inf_k = pi_inf_k + alpha_k(i)*pi_inf_i
                qv_k = qv_k + alpha_rho_k(i)*qvs(i)
            end do
        end if
 
    end subroutine s_compute_mixture_coefficients
 
    !> Time derivative of the mixture coefficients, mirroring s_compute_mixture_coefficients.
    subroutine s_compute_mixture_coefficients_dt(dalpha_rho_dt, dadv_dt, alpha_rho, adv, drho_dt, dgamma_dt, dpi_inf_dt, dqv_dt)
 
 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#ifdef _CRAYFTN
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#else
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!DIR$ INLINEALWAYS s_compute_mixture_coefficients_dt
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenACC
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 615 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
# 620 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
            real(wp), dimension(num_fluids), intent(in) :: dalpha_rho_dt, dadv_dt, alpha_rho, adv
# 622 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        real(wp), intent(out) :: drho_dt, dgamma_dt, dpi_inf_dt, dqv_dt
        real(wp)              :: rho_i, gamma_i, pi_inf_i, dpi_i, dgamma_i, alpha_i, alpha_rho_i
        integer               :: i  !< Loop iterator over fluids
 
        dgamma_dt = 0._wp
        dpi_inf_dt = 0._wp
        dqv_dt = 0._wp
 
        if (num_fluids == 1 .and. bubbles_euler) then
            ! Fluid 1's coefficients are constants here, so only rho varies.
            drho_dt = dalpha_rho_dt(1)
        else
            drho_dt = 0._wp
 
 
# 636 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 636 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 636 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 636 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 636 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
            do i = 1, num_fluids
                drho_dt = drho_dt + dalpha_rho_dt(i)
                alpha_rho_i = alpha_rho(i)
                alpha_i = adv(i)
                call s_phase_coefficients(alpha_rho_i, alpha_i, i, rho_i, gamma_i, pi_inf_i, dpi_i, dgamma_i)
                ! d(alpha X(rho_i))/dt with rho_i = alpha_rho/alpha; the alpha in dX/dt cancels
                dgamma_dt = dgamma_dt + dadv_dt(i)*gamma_i + dgamma_i*(dalpha_rho_dt(i) - rho_i*dadv_dt(i))
                dpi_inf_dt = dpi_inf_dt + dadv_dt(i)*pi_inf_i + dpi_i*(dalpha_rho_dt(i) - rho_i*dadv_dt(i))
                dqv_dt = dqv_dt + dalpha_rho_dt(i)*qvs(i)
            end do
        end if
 
    end subroutine s_compute_mixture_coefficients_dt
 
    !> Speed of sound of a thermodynamic state. Enthalpy is not an argument: for a real state H, |u|^2 and qv all cancel out of c^2
    !! = ((Gamma + 1)p + Pi)/(Gamma rho). Averaged states, whose enthalpy is a free input, use the _avg variant.
    subroutine s_compute_speed_of_sound(pres, rho, gamma, pi_inf, adv, c, alpha_rho)
 
 
# 655 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 655 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 655 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 655 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 655 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 655 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 655 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: pres, rho, gamma, pi_inf
# 661 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
            real(wp), dimension(num_fluids), intent(in) :: adv
# 663 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        real(wp), intent(out) :: c
# 667 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
            real(wp), dimension(num_fluids), intent(in), optional :: alpha_rho
# 669 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        real(wp) :: alf  !< Subgrid void fraction; dilute by construction
        real(wp) :: blkmod_q, alpha_q, alpha_rho_q, gamma_q, pi_inf_q
        integer  :: q
 
        if (chemistry) then  ! Reacting mixture sound speed
            c = sqrt((1.0_wp + 1.0_wp/gamma)*pres/rho)
        else if (relativity) then  ! Relativistic sound speed, whose enthalpy is 1 + (Gamma + 1)p/rho
            c = sqrt((1._wp + 1._wp/gamma)*pres/rho/f_relativistic_enthalpy(pres, rho, gamma))
        else
            ! Every case below is a bulk modulus over a density. The equation of state enters
            ! only through f_bulk_modulus; the cases differ in how the phases are mixed.
            if (any_state_dependent_eos .and. present(alpha_rho)) then  ! frozen mixing: each phase's modulus at its own density
                c = 0._wp
 
# 682 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 682 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 682 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 682 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 682 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
                do q = 1, num_fluids
                    alpha_q = adv(q)
                    alpha_rho_q = alpha_rho(q)
                    call s_phase_bulk_modulus(pres, alpha_q, alpha_rho_q, q, blkmod_q)
                    if (alt_soundspeed) then
                        c = c + adv(q)/blkmod_q
                    else
                        c = c + adv(q)*blkmod_q
                    end if
                end do
                if (alt_soundspeed) then
                    c = 1._wp/(rho*c)
                else
                    c = c/rho
                end if
            else if (alt_soundspeed) then  ! Wood's law: volume-weighted harmonic mean
                c = 0._wp
 
# 700 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 700 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 700 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 700 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 700 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
                do q = 1, num_fluids
                    gamma_q = gammas(q)
                    pi_inf_q = pi_infs(q)
                    c = c + adv(q)/f_bulk_modulus(pres, gamma_q, pi_inf_q)
                end do
                c = 1._wp/(rho*c)
            else if (model_eqns == model_eqns_6eq) then  ! volume-weighted arithmetic mean
                c = 0._wp
 
# 709 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if defined(MFC_OpenACC)
# 709 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc loop seq 
# 709 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif defined(MFC_OpenMP)
# 709 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 709 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
                do q = 1, num_fluids
                    gamma_q = gammas(q)
                    pi_inf_q = pi_infs(q)
                    c = c + adv(q)*f_bulk_modulus(pres, gamma_q, pi_inf_q)
                end do
                c = c/rho
            else  ! the mixture coefficients already carry the mixing
                c = f_bulk_modulus(pres, gamma, pi_inf)/rho
 
                ! Subgrid bubbles: c = c_l/(1 - alf), the carrier-liquid speed with an O(alf) void
                ! correction. alf is dilute by construction; near one means a wrong index or an
                ! out-of-regime case, which the toolchain warns about at case load (#1793).
                if (model_eqns == model_eqns_5eq .and. bubbles_euler .and. .not. (mpp_lim .and. num_fluids > 1)) then
                    alf = adv(num_fluids)
                    c = c/(1._wp - alf)
                end if
            end if
 
            if (mixture_err .and. c < 0._wp) then
                c = 100._wp*sgm_eps
            else
                c = sqrt(c)
            end if
        end if
 
    end subroutine s_compute_speed_of_sound
 
    !> Speed of sound of an interface-averaged state. An average of two states is not a state - its enthalpy is not the one its
    !! pressure and density imply - so the caller supplies H, |u|^2 and qv. Only the enthalpy-reading branches differ from
    !! s_compute_speed_of_sound; keep the condition below in step with the branch list there.
    subroutine s_compute_speed_of_sound_avg(pres, rho, gamma, pi_inf, qv, vel_sum, H, c_c, adv, c, alpha_rho)
 
 
# 742 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#if MFC_OpenACC
# 742 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$acc routine seq 
# 742 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#elif MFC_OpenMP
# 742 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 742 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
# 742 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
!$omp declare target device_type(any) 
# 742 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
#endif
 
        real(wp), intent(in) :: pres, rho, gamma, pi_inf, qv, vel_sum, h, c_c
# 748 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
            real(wp), dimension(num_fluids), intent(in) :: adv
# 750 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
        real(wp), intent(out) :: c
# 754 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
            real(wp), dimension(num_fluids), intent(in), optional :: alpha_rho
# 756 "/home/runner/work/MFC/MFC/src/common/m_eos.fpp"
 
        if (chemistry) then  ! Reacting mixture sound speed
            if (avg_state == avg_state_roe .and. abs(c_c) > verysmall) then
                c = sqrt(c_c - (gamma - 1.0_wp)*(vel_sum - h))
            else
                call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, adv, c, alpha_rho)
            end if
        else if (relativity) then  ! Relativistic sound speed
            c = sqrt((1._wp + 1._wp/gamma)*pres/rho/h)
        else if (alt_soundspeed .or. model_eqns == model_eqns_6eq .or. (model_eqns == model_eqns_5eq .and. bubbles_euler) &
                 & .or. any_state_dependent_eos) then
            call s_compute_speed_of_sound(pres, rho, gamma, pi_inf, adv, c, alpha_rho)
        else  ! Stiffened-gas mixture, the one branch where the averaged enthalpy survives
            c = (h - 5.e-1*vel_sum - qv/rho)/gamma
 
            if (mixture_err .and. c < 0._wp) then
                c = 100._wp*sgm_eps
            else
                c = sqrt(c)
            end if
        end if
 
    end subroutine s_compute_speed_of_sound_avg
 
end module m_eos