recover_temperature_and_update_thermo_cantera_owned Subroutine

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

Owned-cell thermo sync for the timestep hot path.

The legacy helper above updates the replicated full mesh. During a timestep each MPI rank only needs to evaluate Cantera on its owned cells; the caller performs the normal halo exchanges for stencil fields. This preserves the separate radiation MPI communicator because all exchanges still go through flow_mpi_t.

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

Source Code

   subroutine recover_temperature_and_update_thermo_cantera_owned(mesh, flow, params, energy, species_Y)
      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(rk), intent(in), optional :: species_Y(:,:)

      integer :: n_len, nloc, i, c, k
      real(rk) :: sum_Y
      real(rk), allocatable :: h_local(:), P_arr(:), Y_local(:,:)
      real(rk), allocatable :: T_local(:), cp_local(:), lambda_local(:), rho_local(:)
      character(kind=c_char), allocatable :: c_names_flat(:)

      if (.not. params%enable_cantera_thermo) return
      if (params%nspecies <= 0) then
         call fatal_error('energy', 'enable_cantera_thermo requires at least one inert thermo species')
      end if

      if (mesh%ncells < flow%last_cell) then
         call fatal_error('energy', 'flow ownership exceeds mesh size in owned Cantera thermo sync')
      end if

      nloc = max(0, flow%last_cell - flow%first_cell + 1)
      if (nloc <= 0) return

      allocate(h_local(nloc), P_arr(nloc), Y_local(max(1, params%nspecies), nloc))
      allocate(T_local(nloc), cp_local(nloc), lambda_local(nloc), rho_local(nloc))
      P_arr = params%background_press

      do i = 1, nloc
         c = flow%first_cell + i - 1
         h_local(i) = energy%h(c)
         if (present(species_Y) .and. params%nspecies > 0) then
            sum_Y = zero
            do k = 1, params%nspecies
               Y_local(k, i) = max(zero, species_Y(k, c))
               sum_Y = sum_Y + Y_local(k, i)
            end do
            if (sum_Y > tiny_safe) then
               Y_local(1:params%nspecies, i) = Y_local(1:params%nspecies, i) / sum_Y
            else
               Y_local(:, i) = zero
               Y_local(1, i) = 1.0_rk
            end if
         else
            Y_local(:, i) = zero
            Y_local(1, i) = 1.0_rk
         end if
      end do

      call build_c_species_names(params, c_names_flat, n_len)
      call cantera_recover_temperature_and_update_thermo_c(nloc, h_local, P_arr, params%nspecies, &
                                   Y_local, T_local, cp_local, lambda_local, rho_local, &
                                   params%energy_reference_T, c_names_flat, n_len)

      do i = 1, nloc
         c = flow%first_cell + i - 1
         energy%T(c) = T_local(i)
         energy%cp(c) = cp_local(i)
         energy%lambda(c) = lambda_local(i)
         energy%rho_thermo(c) = rho_local(i)
      end do

      deallocate(h_local, P_arr, Y_local, T_local, cp_local, lambda_local, rho_local, c_names_flat)
   end subroutine recover_temperature_and_update_thermo_cantera_owned