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.

The energy equation uses constant flow/projection density params%rho. energy%rho_thermo is diagnostic/future-use only and is not used by the projection or momentum equations. Heat conduction uses grad(T), not grad(h). qrad is a volumetric source with qrad > 0 adding energy to the gas.

Missing physics in this path: reactions, reaction heat release, variable-density low-Mach coupling, species-diffusion enthalpy correction -div(sum_k h_k J_k), and external radiation physics.


Uses

  • module~~mod_energy~~UsesGraph module~mod_energy mod_energy iso_c_binding iso_c_binding module~mod_energy->iso_c_binding module~mod_bc mod_bc module~mod_energy->module~mod_bc module~mod_fields mod_fields module~mod_energy->module~mod_fields module~mod_input mod_input module~mod_energy->module~mod_input module~mod_kinds mod_kinds module~mod_energy->module~mod_kinds module~mod_mesh_types mod_mesh_types module~mod_energy->module~mod_mesh_types module~mod_mpi_flow mod_mpi_flow module~mod_energy->module~mod_mpi_flow module~mod_profiler mod_profiler module~mod_energy->module~mod_profiler mpi_f08 mpi_f08 module~mod_energy->mpi_f08 module~mod_bc->module~mod_input module~mod_bc->module~mod_kinds module~mod_bc->module~mod_mesh_types module~mod_fields->module~mod_bc module~mod_fields->module~mod_input module~mod_fields->module~mod_kinds module~mod_fields->module~mod_mesh_types module~mod_input->module~mod_kinds iso_fortran_env iso_fortran_env module~mod_kinds->iso_fortran_env module~mod_mesh_types->module~mod_kinds module~mod_mpi_flow->module~mod_kinds module~mod_mpi_flow->module~mod_mesh_types module~mod_mpi_flow->mpi_f08 module~mod_profiler->module~mod_kinds module~mod_profiler->mpi_f08 module~mod_profiler->iso_fortran_env

Used by

  • module~~mod_energy~~UsedByGraph module~mod_energy mod_energy module~mod_output mod_output module~mod_output->module~mod_energy program~lowmach_react_hex lowmach_react_hex program~lowmach_react_hex->module~mod_energy program~lowmach_react_hex->module~mod_output

Interfaces

interface

C-Binding interface for Cantera thermodynamics.

  • 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.

  • 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.

  • 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 :: cp(:)
real(kind=rk), public, allocatable :: h(:)
real(kind=rk), public, allocatable :: h_old(:)
logical, public :: initialized = .false.
real(kind=rk), public, allocatable :: lambda(:)
real(kind=rk), public, allocatable :: qrad(:)
real(kind=rk), public, allocatable :: rho_thermo(:)

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_energy_transport(mesh, flow, bc, params, fields, energy, species_Y)

Advance transported sensible enthalpy with constant flow density.

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(:,:)

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 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_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 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_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)

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

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 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.

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(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(:,:)

private 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(:,:)

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(:,:)