mod_energy Module

Energy field storage, diagnostics, and sensible-enthalpy transport.

This module owns the enthalpy/temperature/radiation-source fields used by the current constant-density energy path. The transported thermodynamic state is h; temperature is recovered from h, composition Y, and the background thermodynamic pressure params%background_press.

In Cantera mode the stored sensible enthalpy is h_abs(T,Y,p0) - h_abs(T_ref,Y,p0), where p0 is the background pressure and T_ref is params%energy_reference_T. After species transport changes Y, Option A is enforced: preserve transported h and recover T(h,Y,p0). Do not rebuild h from the old temperature and new composition.

Density ownership depends on mode. In constant-density mode transport%rho = params%rho and energy%rho_thermo is diagnostic. In guarded variable-density mode with density_eos="cantera", transport%rho <- energy%rho_thermo after Cantera thermo sync, so density comes from the selected Cantera phase and its YAML EOS. Heat conduction uses grad(T), not grad(h). qrad is a volumetric source hook with qrad > 0 adding energy to the gas; external radiation coupling is future work.

Reactions and reaction heat release are not implemented in the supported path. The optional species-enthalpy diffusion correction -div(sum_k h_k J_k) is controlled by enable_species_enthalpy_diffusion.



Variables

Type Visibility Attributes Name Initial
real(kind=rk), private, save :: chemistry_stats_active_cells_local = zero
real(kind=rk), private, save :: chemistry_stats_all_skipped_updates_local = zero
real(kind=rk), private, save :: chemistry_stats_owned_cells_local = zero
real(kind=rk), private, save :: chemistry_stats_reactor_calls_local = zero
real(kind=rk), private, save :: chemistry_stats_reactor_time_local = zero
real(kind=rk), private, save :: chemistry_stats_skipped_cells_local = zero
real(kind=rk), private, save :: chemistry_stats_skipped_legacy_mass_fraction_local = zero
real(kind=rk), private, save :: chemistry_stats_skipped_named_species_local = zero
real(kind=rk), private, save :: chemistry_stats_skipped_temperature_local = zero
real(kind=rk), private, save :: chemistry_stats_total_time_local = zero
real(kind=rk), private, save :: chemistry_stats_updates_local = zero

Local accumulated chemistry-screening statistics. These are reduced across flow ranks at shutdown by write_chemistry_screening_stats.


Interfaces

interface

C-Binding interface for Cantera thermodynamics.

The C++ bridge expects temperatures , thermodynamic pressures , and normalized solver species mass fractions . It returns sensible enthalpy using and recovers temperature by adding the same-composition reference enthalpy before Cantera HP inversion. Projection pressure is never passed to these calls.

  • private subroutine cantera_get_cache_stats_c(stats_out, nstats) bind(c, name="cantera_get_cache_stats_c")

    Arguments

    Type IntentOptional Attributes Name
    real(kind=c_double), intent(out) :: stats_out(*)
    integer(kind=c_int), value :: nstats

interface

C-Binding interface for Cantera thermodynamics.

The C++ bridge expects temperatures , thermodynamic pressures , and normalized solver species mass fractions . It returns sensible enthalpy using and recovers temperature by adding the same-composition reference enthalpy before Cantera HP inversion. Projection pressure is never passed to these calls.

  • private subroutine cantera_heat_release_rate_c(npoints, T, P, nspecies, Y_in, qdot_out, species_names_flat, name_len) bind(c, name="cantera_heat_release_rate_c")

    Arguments

    Type IntentOptional Attributes Name
    integer(kind=c_int), value :: npoints
    real(kind=c_double), intent(in) :: T(npoints)
    real(kind=c_double), intent(in) :: P(npoints)
    integer(kind=c_int), value :: nspecies
    real(kind=c_double), intent(in) :: Y_in(*)
    real(kind=c_double), intent(out) :: qdot_out(npoints)
    character(kind=c_char, len=1), intent(in) :: species_names_flat(*)
    integer(kind=c_int), value :: name_len

interface

C-Binding interface for Cantera thermodynamics.

The C++ bridge expects temperatures , thermodynamic pressures , and normalized solver species mass fractions . It returns sensible enthalpy using and recovers temperature by adding the same-composition reference enthalpy before Cantera HP inversion. Projection pressure is never passed to these calls.

  • private subroutine cantera_integrate_const_p_chemistry_c(npoints, dt, T_inout, P, nspecies, Y_inout, h_out, cp_out, lambda_out, rho_out, T_ref, species_names_flat, name_len, rtol, atol, max_steps, energy_enabled) bind(c, name="cantera_integrate_const_p_chemistry_c")

    Arguments

    Type IntentOptional Attributes Name
    integer(kind=c_int), value :: npoints
    real(kind=c_double), value :: dt
    real(kind=c_double), intent(inout) :: T_inout(npoints)
    real(kind=c_double), intent(in) :: P(npoints)
    integer(kind=c_int), value :: nspecies
    real(kind=c_double), intent(inout) :: Y_inout(*)
    real(kind=c_double), intent(out) :: h_out(npoints)
    real(kind=c_double), intent(out) :: cp_out(npoints)
    real(kind=c_double), intent(out) :: lambda_out(npoints)
    real(kind=c_double), intent(out) :: rho_out(npoints)
    real(kind=c_double), value :: T_ref
    character(kind=c_char, len=1), intent(in) :: species_names_flat(*)
    integer(kind=c_int), value :: name_len
    real(kind=c_double), value :: rtol
    real(kind=c_double), value :: atol
    integer(kind=c_int), value :: max_steps
    integer(kind=c_int), value :: energy_enabled

interface

C-Binding interface for Cantera thermodynamics.

The C++ bridge expects temperatures , thermodynamic pressures , and normalized solver species mass fractions . It returns sensible enthalpy using and recovers temperature by adding the same-composition reference enthalpy before Cantera HP inversion. Projection pressure is never passed to these calls.

  • private subroutine cantera_recover_temperature_and_update_thermo_c(ncells, h_in, P, nspecies, Y_in, T_out, cp_out, lambda_out, rho_thermo_out, T_ref, species_names_flat, name_len) bind(c, name="cantera_recover_temperature_and_update_thermo_c")

    Arguments

    Type IntentOptional Attributes Name
    integer(kind=c_int), value :: ncells
    real(kind=c_double), intent(in) :: h_in(ncells)
    real(kind=c_double), intent(in) :: P(ncells)
    integer(kind=c_int), value :: nspecies
    real(kind=c_double), intent(in) :: Y_in(*)
    real(kind=c_double), intent(out) :: T_out(ncells)
    real(kind=c_double), intent(out) :: cp_out(ncells)
    real(kind=c_double), intent(out) :: lambda_out(ncells)
    real(kind=c_double), intent(out) :: rho_thermo_out(ncells)
    real(kind=c_double), value :: T_ref
    character(kind=c_char, len=1), intent(in) :: species_names_flat(*)
    integer(kind=c_int), value :: name_len

interface

C-Binding interface for Cantera thermodynamics.

The C++ bridge expects temperatures , thermodynamic pressures , and normalized solver species mass fractions . It returns sensible enthalpy using and recovers temperature by adding the same-composition reference enthalpy before Cantera HP inversion. Projection pressure is never passed to these calls.

  • private subroutine cantera_recover_temperature_from_h_c(ncells, h_in, P, nspecies, Y_in, T_out, T_ref, species_names_flat, name_len) bind(c, name="cantera_recover_temperature_from_h_c")

    Arguments

    Type IntentOptional Attributes Name
    integer(kind=c_int), value :: ncells
    real(kind=c_double), intent(in) :: h_in(ncells)
    real(kind=c_double), intent(in) :: P(ncells)
    integer(kind=c_int), value :: nspecies
    real(kind=c_double), intent(in) :: Y_in(*)
    real(kind=c_double), intent(out) :: T_out(ncells)
    real(kind=c_double), value :: T_ref
    character(kind=c_char, len=1), intent(in) :: species_names_flat(*)
    integer(kind=c_int), value :: name_len

interface

C-Binding interface for Cantera thermodynamics.

The C++ bridge expects temperatures , thermodynamic pressures , and normalized solver species mass fractions . It returns sensible enthalpy using and recovers temperature by adding the same-composition reference enthalpy before Cantera HP inversion. Projection pressure is never passed to these calls.

  • private subroutine cantera_recover_temperature_update_thermo_and_species_h_c(npoints, h_in, P, nspecies, Y_in, T_out, cp_out, lambda_out, rho_thermo_out, hk_out, T_ref, species_names_flat, name_len) bind(c, name="cantera_recover_temperature_update_thermo_and_species_h_c")

    Arguments

    Type IntentOptional Attributes Name
    integer(kind=c_int), value :: npoints
    real(kind=c_double), intent(in) :: h_in(npoints)
    real(kind=c_double), intent(in) :: P(npoints)
    integer(kind=c_int), value :: nspecies
    real(kind=c_double), intent(in) :: Y_in(*)
    real(kind=c_double), intent(out) :: T_out(npoints)
    real(kind=c_double), intent(out) :: cp_out(npoints)
    real(kind=c_double), intent(out) :: lambda_out(npoints)
    real(kind=c_double), intent(out) :: rho_thermo_out(npoints)
    real(kind=c_double), intent(out) :: hk_out(*)
    real(kind=c_double), value :: T_ref
    character(kind=c_char, len=1), intent(in) :: species_names_flat(*)
    integer(kind=c_int), value :: name_len

interface

C-Binding interface for Cantera thermodynamics.

The C++ bridge expects temperatures , thermodynamic pressures , and normalized solver species mass fractions . It returns sensible enthalpy using and recovers temperature by adding the same-composition reference enthalpy before Cantera HP inversion. Projection pressure is never passed to these calls.

  • private subroutine cantera_species_sensible_enthalpies_c(npoints, T, P, nspecies, Y_in, hk_out, T_ref, species_names_flat, name_len) bind(c, name="cantera_species_sensible_enthalpies_c")

    Arguments

    Type IntentOptional Attributes Name
    integer(kind=c_int), value :: npoints
    real(kind=c_double), intent(in) :: T(npoints)
    real(kind=c_double), intent(in) :: P(npoints)
    integer(kind=c_int), value :: nspecies
    real(kind=c_double), intent(in) :: Y_in(*)
    real(kind=c_double), intent(out) :: hk_out(*)
    real(kind=c_double), value :: T_ref
    character(kind=c_char, len=1), intent(in) :: species_names_flat(*)
    integer(kind=c_int), value :: name_len

interface

C-Binding interface for Cantera thermodynamics.

The C++ bridge expects temperatures , thermodynamic pressures , and normalized solver species mass fractions . It returns sensible enthalpy using and recovers temperature by adding the same-composition reference enthalpy before Cantera HP inversion. Projection pressure is never passed to these calls.

  • private subroutine cantera_update_thermo_c(ncells, T, P, nspecies, Y_in, h_out, cp_out, lambda_out, rho_thermo_out, T_ref, species_names_flat, name_len) bind(c, name="cantera_update_thermo_c")

    Arguments

    Type IntentOptional Attributes Name
    integer(kind=c_int), value :: ncells
    real(kind=c_double), intent(in) :: T(ncells)
    real(kind=c_double), intent(in) :: P(ncells)
    integer(kind=c_int), value :: nspecies
    real(kind=c_double), intent(in) :: Y_in(*)
    real(kind=c_double), intent(out) :: h_out(ncells)
    real(kind=c_double), intent(out) :: cp_out(ncells)
    real(kind=c_double), intent(out) :: lambda_out(ncells)
    real(kind=c_double), intent(out) :: rho_thermo_out(ncells)
    real(kind=c_double), value :: T_ref
    character(kind=c_char, len=1), intent(in) :: species_names_flat(*)
    integer(kind=c_int), value :: name_len

Derived Types

type, public ::  energy_fields_t

Cell-centered energy variables.

Components

Type Visibility Attributes Name Initial
real(kind=rk), public, allocatable :: T(:)
real(kind=rk), public, allocatable :: T_old(:)
real(kind=rk), public, allocatable :: chem_active(:)
real(kind=rk), public, allocatable :: chem_dT_dt(:)
real(kind=rk), public, allocatable :: chem_dY_max_dt(:)
real(kind=rk), public, allocatable :: chem_heat_release_rate(:)
real(kind=rk), public, allocatable :: chem_rel_rho_change_dt(:)
real(kind=rk), public, allocatable :: chem_source_alpha(:)
real(kind=rk), public, allocatable :: cp(:)
real(kind=rk), public :: cumulative_boundary_rho_h_advective_flux_out = zero
real(kind=rk), public :: cumulative_boundary_rho_h_conductive_flux_out = zero
integer, public :: cumulative_energy_budget_available = 0
real(kind=rk), public :: cumulative_energy_update_delta_integral = zero
real(kind=rk), public :: cumulative_energy_update_rhs_integral = zero
real(kind=rk), public :: cumulative_qrad_integral = zero
real(kind=rk), public :: cumulative_rho_species_hdiff_integral = zero
real(kind=rk), public, allocatable :: h(:)
real(kind=rk), public, allocatable :: h_old(:)
logical, public :: initialized = .false.
real(kind=rk), public, allocatable :: lambda(:)
integer, public :: last_conductive_boundary_flux_available = 0
real(kind=rk), public :: last_conductive_boundary_flux_out = zero
real(kind=rk), public :: last_energy_update_balance_defect = zero
real(kind=rk), public :: last_energy_update_delta_rate_integral = zero
real(kind=rk), public :: last_energy_update_rhs_integral = zero
real(kind=rk), public :: last_operator_consistent_rho_h_integral = zero
real(kind=rk), public, allocatable :: operator_consistent_rho_h(:)
integer, public :: operator_consistent_rho_h_available = 0
real(kind=rk), public, allocatable :: qrad(:)
real(kind=rk), public :: relative_last_energy_update_balance_defect = zero
real(kind=rk), public, allocatable :: rho_thermo(:)
real(kind=rk), public, allocatable :: rhs_old(:)
logical, public :: rhs_old_valid = .false.
real(kind=rk), public, allocatable :: species_enthalpy_diffusion(:)

Functions

private function boundary_enthalpy_from_temperature(mesh, params, temperature, Y_point) result(h_value)

Convert a boundary temperature to enthalpy using the active thermo model. Boundary enthalpy from a boundary temperature.

Read more…

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(case_params_t), intent(in) :: params
real(kind=rk), intent(in) :: temperature
real(kind=rk), intent(in), optional :: Y_point(:)

Return Value real(kind=rk)

private function energy_face_normal_distance(mesh, bc, face_id, cell_id, nb) result(dist)

Normal distance used for cell-cell and cell-boundary temperature gradients.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(bc_set_t), intent(in) :: bc
integer, intent(in) :: face_id
integer, intent(in) :: cell_id
integer, intent(in) :: nb

Return Value real(kind=rk)

private function energy_outward_normal(mesh, face_id, cell_id) result(nvec)

Outward unit normal from cell_id for face_id.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
integer, intent(in) :: face_id
integer, intent(in) :: cell_id

Return Value real(kind=rk), (3)

private pure function enthalpy_from_temperature_value(params, temperature) result(h_value)

Constant-cp helper for a single boundary/face temperature value.

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
real(kind=rk), intent(in) :: temperature

Return Value real(kind=rk)


Subroutines

public subroutine advance_cantera_chemistry(mesh, flow, params, species_Y, energy, step, chemistry_dt, max_dT_chem, max_dY_chem, max_rel_drho_chem, min_source_alpha, raw_max_dT_chem, raw_max_dY_chem, raw_max_rel_drho_chem)

Advance cell-local Cantera chemistry on owned cells.

Read more…

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(flow_mpi_t), intent(inout) :: flow
type(case_params_t), intent(in) :: params
real(kind=rk), intent(inout) :: species_Y(:,:)
type(energy_fields_t), intent(inout) :: energy
integer, intent(in) :: step
real(kind=rk), intent(in), optional :: chemistry_dt
real(kind=rk), intent(out), optional :: max_dT_chem
real(kind=rk), intent(out), optional :: max_dY_chem
real(kind=rk), intent(out), optional :: max_rel_drho_chem
real(kind=rk), intent(out), optional :: min_source_alpha
real(kind=rk), intent(out), optional :: raw_max_dT_chem
real(kind=rk), intent(out), optional :: raw_max_dY_chem
real(kind=rk), intent(out), optional :: raw_max_rel_drho_chem

public subroutine advance_energy_transport(mesh, flow, bc, params, fields, energy, species_Y, transport, step)

Advance transported sensible enthalpy one explicit finite-volume step.

Read more…

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(flow_mpi_t), intent(inout) :: flow
type(bc_set_t), intent(in) :: bc
type(case_params_t), intent(in) :: params
type(flow_fields_t), intent(in) :: fields
type(energy_fields_t), intent(inout) :: energy
real(kind=rk), intent(in), optional :: species_Y(:,:)
type(transport_properties_t), intent(in), optional :: transport
integer, intent(in), optional :: step

public subroutine allocate_energy(mesh, energy)

Allocate all energy arrays for the mesh.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(energy_fields_t), intent(inout) :: energy

public subroutine compute_species_sensible_enthalpies_cantera(mesh, params, T_state, species_Y, hk_species)

Compute species sensible enthalpies h_k(T)-h_k(T_ref) [J/kg_k] for all cells.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(case_params_t), intent(in) :: params
real(kind=rk), intent(in) :: T_state(:)
real(kind=rk), intent(in) :: species_Y(:,:)
real(kind=rk), intent(out) :: hk_species(:,:)

public subroutine finalize_energy(energy)

Deallocate all energy arrays.

Arguments

Type IntentOptional Attributes Name
type(energy_fields_t), intent(inout) :: energy

public subroutine initialize_energy(mesh, params, energy, species_Y)

Initialize energy fields from case parameters.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy
real(kind=rk), intent(in), optional :: species_Y(:,:)

public subroutine recover_temperature_and_update_thermo_cantera(mesh, params, energy, species_Y)

Recover T from h and refresh cp/lambda/rho_thermo in one Cantera sync.

Read more…

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy
real(kind=rk), intent(in), optional :: species_Y(:,:)

public subroutine recover_temperature_cantera(mesh, params, energy, species_Y)

Recover T from h using Cantera HPY inversion.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy
real(kind=rk), intent(in), optional :: species_Y(:,:)

public subroutine recover_temperature_constant_cp(params, energy)

Recover T from h using the constant-cp thermodynamic model.

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy

public subroutine refresh_energy_from_temperature(mesh, params, energy, species_Y)

Refresh h/cp/lambda/rho_thermo from the current temperature field.

Read more…

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy
real(kind=rk), intent(in), optional :: species_Y(:,:)

public subroutine update_enthalpy_from_temperature_constant_cp(params, energy)

Update h from T using the constant-cp thermodynamic model.

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy

public subroutine write_cantera_cache_stats(flow)

Print bridge-level Cantera cache statistics, summed over flow ranks.

Read more…

Arguments

Type IntentOptional Attributes Name
type(flow_mpi_t), intent(in) :: flow

public subroutine write_chemistry_diagnostics_header(params, flow)

Initialize chemistry workload diagnostics CSV.

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
type(flow_mpi_t), intent(in) :: flow

public subroutine write_chemistry_screening_stats(flow, params)

Print accumulated chemistry workload and screening effectiveness statistics.

Read more…

Arguments

Type IntentOptional Attributes Name
type(flow_mpi_t), intent(in) :: flow
type(case_params_t), intent(in) :: params

public subroutine write_energy_diagnostics_header(params, flow)

Writes the CSV header for energy diagnostics.

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
type(flow_mpi_t), intent(in) :: flow

public subroutine write_energy_diagnostics_row(mesh, flow, params, energy, step, time, write_file)

Appends one row of global energy diagnostics.

Read more…

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(flow_mpi_t), intent(in) :: flow
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(in) :: energy
integer, intent(in) :: step
real(kind=rk), intent(in) :: time
logical, intent(in), optional :: write_file

public subroutine zero_radiation_source(energy)

Reset the volumetric energy source to zero.

Arguments

Type IntentOptional Attributes Name
type(energy_fields_t), intent(inout) :: energy

private subroutine accumulate_chemistry_screening_stats(nowned, nactive, skipped_temperature, skipped_named_species, skipped_legacy_mass_fraction, total_time, reactor_time, reactor_calls)

Accumulate local chemistry workload/screening statistics.

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: nowned
integer, intent(in) :: nactive
integer, intent(in) :: skipped_temperature
integer, intent(in) :: skipped_named_species
integer, intent(in) :: skipped_legacy_mass_fraction
real(kind=rk), intent(in) :: total_time
real(kind=rk), intent(in) :: reactor_time
integer, intent(in), optional :: reactor_calls

private subroutine build_boundary_thermo_Y(mesh, bc, params, face_id, cell_id, species_Y, Y_point)

Build a boundary thermodynamic composition vector for fixed-T boundaries.

Read more…

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(bc_set_t), intent(in) :: bc
type(case_params_t), intent(in) :: params
integer, intent(in) :: face_id
integer, intent(in) :: cell_id
real(kind=rk), intent(in) :: species_Y(:,:)
real(kind=rk), intent(out) :: Y_point(:)

private subroutine build_c_species_names(params, c_names_flat, n_len)

Build a flattened C-compatible species-name buffer.

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
character(kind=c_char, len=1), intent(out), allocatable :: c_names_flat(:)
integer, intent(out) :: n_len

private subroutine build_thermo_Y(params, ncells, Y_local, species_Y)

Build a cellwise thermodynamic composition array for Cantera calls.

Read more…

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
integer, intent(in) :: ncells
real(kind=rk), intent(out), allocatable :: Y_local(:,:)
real(kind=rk), intent(in), optional :: species_Y(:,:)

private subroutine enthalpy_from_temperature_cantera_point(params, temperature, h_value, Y_point)

Compute h(T,Y,p0) for one boundary state using Cantera.

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
real(kind=rk), intent(in) :: temperature
real(kind=rk), intent(out) :: h_value
real(kind=rk), intent(in), optional :: Y_point(:)

private subroutine recover_temperature_and_update_thermo_cantera_owned(mesh, flow, params, energy, species_Y)

Owned-cell thermo sync for the timestep hot path.

Read more…

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(flow_mpi_t), intent(in) :: flow
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy
real(kind=rk), intent(in), optional :: species_Y(:,:)

private subroutine recover_temperature_update_thermo_and_species_h_cantera_owned(mesh, flow, params, energy, species_Y, hk_species)

Owned-cell fused thermo sync plus species sensible enthalpies.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(flow_mpi_t), intent(inout) :: flow
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy
real(kind=rk), intent(in) :: species_Y(:,:)
real(kind=rk), intent(out) :: hk_species(:,:)

private subroutine resolve_chemistry_active_species(params, nactive_species, active_species_indices)

Resolve optional named chemistry active-species screening names to current species indices.

Read more…

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
integer, intent(out) :: nactive_species
integer, intent(out) :: active_species_indices(:)

private subroutine species_sensible_enthalpies_from_temperature_cantera_point(params, temperature, hk_value, Y_point)

Compute species sensible enthalpies for one boundary/face state.

Arguments

Type IntentOptional Attributes Name
type(case_params_t), intent(in) :: params
real(kind=rk), intent(in) :: temperature
real(kind=rk), intent(out) :: hk_value(:)
real(kind=rk), intent(in) :: Y_point(:)

private subroutine update_thermo_from_temperature_cantera(mesh, params, energy, species_Y)

Update h, cp, lambda, and diagnostic thermo density from current T.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(case_params_t), intent(in) :: params
type(energy_fields_t), intent(inout) :: energy
real(kind=rk), intent(in), optional :: species_Y(:,:)

private subroutine write_chemistry_diagnostics_row(flow, params, step, chemistry_dt, owned_local, active_local, skipped_local, total_time_local, reactor_time_local, max_dT_local, max_dY_local, max_rel_drho_local, min_source_alpha_local, raw_max_dT_local, raw_max_dY_local, raw_max_rel_drho_local, max_qdot_local, integrated_qdot_local)

Arguments

Type IntentOptional Attributes Name
type(flow_mpi_t), intent(in) :: flow
type(case_params_t), intent(in) :: params
integer, intent(in) :: step
real(kind=rk), intent(in) :: chemistry_dt
integer, intent(in) :: owned_local
integer, intent(in) :: active_local
integer, intent(in) :: skipped_local
real(kind=rk), intent(in) :: total_time_local
real(kind=rk), intent(in) :: reactor_time_local
real(kind=rk), intent(in), optional :: max_dT_local
real(kind=rk), intent(in), optional :: max_dY_local
real(kind=rk), intent(in), optional :: max_rel_drho_local
real(kind=rk), intent(in), optional :: min_source_alpha_local
real(kind=rk), intent(in), optional :: raw_max_dT_local
real(kind=rk), intent(in), optional :: raw_max_dY_local
real(kind=rk), intent(in), optional :: raw_max_rel_drho_local
real(kind=rk), intent(in), optional :: max_qdot_local
real(kind=rk), intent(in), optional :: integrated_qdot_local

private subroutine write_variable_density_energy_debug(mesh, flow, params, fields, energy, label, species_Y, transport, step)

Targeted diagnostics for experimental variable-density energy/Cantera HP failures.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(flow_mpi_t), intent(in) :: flow
type(case_params_t), intent(in) :: params
type(flow_fields_t), intent(in) :: fields
type(energy_fields_t), intent(in) :: energy
character(len=*), intent(in) :: label
real(kind=rk), intent(in), optional :: species_Y(:,:)
type(transport_properties_t), intent(in), optional :: transport
integer, intent(in), optional :: step