build_thermo_Y Subroutine

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

Build a cellwise thermodynamic composition array for Cantera calls.

If transported species are available, they are used directly and clipped to non-negative values. If their sum is positive, the local composition is normalized before passing to Cantera. If no transported species are available, use the single-species/default mixture prepared during Cantera initialization, usually thermo_default_species.

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

Calls

proc~~build_thermo_y~~CallsGraph proc~build_thermo_y mod_energy::build_thermo_Y proc~fatal_error mod_kinds::fatal_error proc~build_thermo_y->proc~fatal_error

Called by

proc~~build_thermo_y~~CalledByGraph proc~build_thermo_y mod_energy::build_thermo_Y proc~recover_temperature_and_update_thermo_cantera mod_energy::recover_temperature_and_update_thermo_cantera proc~recover_temperature_and_update_thermo_cantera->proc~build_thermo_y proc~recover_temperature_cantera mod_energy::recover_temperature_cantera proc~recover_temperature_cantera->proc~build_thermo_y proc~update_thermo_from_temperature_cantera mod_energy::update_thermo_from_temperature_cantera proc~update_thermo_from_temperature_cantera->proc~build_thermo_y proc~advance_energy_transport mod_energy::advance_energy_transport proc~advance_energy_transport->proc~recover_temperature_and_update_thermo_cantera proc~initialize_energy mod_energy::initialize_energy proc~initialize_energy->proc~update_thermo_from_temperature_cantera program~lowmach_react_hex lowmach_react_hex program~lowmach_react_hex->proc~advance_energy_transport program~lowmach_react_hex->proc~initialize_energy

Source Code

   subroutine build_thermo_Y(params, ncells, Y_local, species_Y)
      type(case_params_t), intent(in) :: params
      integer, intent(in) :: ncells
      real(rk), allocatable, intent(out) :: Y_local(:,:)
      real(rk), intent(in), optional :: species_Y(:,:)

      integer :: c, k, nsp
      real(rk) :: sum_Y

      nsp = max(1, params%nspecies)
      allocate(Y_local(nsp, ncells))
      Y_local = zero

      if (present(species_Y) .and. params%enable_species .and. params%nspecies > 0) then
         if (size(species_Y, 1) < params%nspecies .or. size(species_Y, 2) < ncells) then
            call fatal_error('energy', 'species_Y has incompatible shape for Cantera thermo update')
         end if

         do c = 1, ncells
            sum_Y = zero
            do k = 1, params%nspecies
               Y_local(k, c) = max(zero, species_Y(k, c))
               sum_Y = sum_Y + Y_local(k, c)
            end do

            if (sum_Y > tiny_safe) then
               Y_local(1:params%nspecies, c) = Y_local(1:params%nspecies, c) / sum_Y
            else
               Y_local(1, c) = 1.0_rk
            end if
         end do
      else
         Y_local(1, :) = 1.0_rk
      end if
   end subroutine build_thermo_Y