write_chemistry_screening_stats Subroutine

   subroutine write_chemistry_screening_stats(flow, params)
      type(flow_mpi_t), intent(in) :: flow
      type(case_params_t), intent(in) :: params

      integer, parameter :: nstats = 11
      real(rk) :: local_stats(nstats), global_stats(nstats)
      real(rk) :: owned_total, active_total, skipped_total, reactor_calls_total
      real(rk) :: skip_fraction, active_fraction, avg_reactor_time_per_call
      real(rk) :: estimated_reactor_time_avoided, update_steps
      logical :: screening_enabled
      integer :: ierr

      if (.not. params%enable_reactions) return

      local_stats = [chemistry_stats_updates_local, &
                     chemistry_stats_owned_cells_local, &
                     chemistry_stats_active_cells_local, &
                     chemistry_stats_skipped_cells_local, &
                     chemistry_stats_reactor_calls_local, &
                     chemistry_stats_skipped_temperature_local, &
                     chemistry_stats_skipped_named_species_local, &
                     chemistry_stats_skipped_legacy_mass_fraction_local, &
                     chemistry_stats_all_skipped_updates_local, &
                     chemistry_stats_total_time_local, &
                     chemistry_stats_reactor_time_local]
      global_stats = zero

      call MPI_Reduce(local_stats, global_stats, nstats, MPI_DOUBLE_PRECISION, MPI_SUM, 0, flow%comm, ierr)
      if (ierr /= MPI_SUCCESS) call fatal_error('chemistry', 'MPI_Reduce failed for chemistry screening statistics')

      if (flow%rank /= 0) return
      if (global_stats(1) <= zero) return

      owned_total = global_stats(2)
      active_total = global_stats(3)
      skipped_total = global_stats(4)
      reactor_calls_total = global_stats(5)
      update_steps = global_stats(1) / real(max(1, flow%nprocs), rk)

      if (owned_total > zero) then
         active_fraction = 100.0_rk * active_total / owned_total
         skip_fraction = 100.0_rk * skipped_total / owned_total
      else
         active_fraction = zero
         skip_fraction = zero
      end if

      if (reactor_calls_total > zero) then
         avg_reactor_time_per_call = global_stats(11) / reactor_calls_total
      else
         avg_reactor_time_per_call = zero
      end if
      estimated_reactor_time_avoided = skipped_total * avg_reactor_time_per_call

      screening_enabled = (params%chemistry_temperature_cutoff > zero) .or. &
                          (params%chemistry_min_reactive_mass_fraction > zero) .or. &
                          (params%chemistry_active_species_threshold > zero .and. &
                           params%chemistry_n_active_species > 0)

      write(output_unit,'(a)') ' ======================================================================'
      write(output_unit,'(a)') '  CHEMISTRY SCREENING STATISTICS'
      write(output_unit,'(a)') '  Counts are summed over flow ranks and chemistry-update calls.'
      write(output_unit,'(a)') ' ======================================================================'
      write(output_unit,'(a,l1)')   '  screening controls enabled        : ', screening_enabled
      write(output_unit,'(a,f12.0)') '  chemistry update steps            : ', update_steps
      write(output_unit,'(a,f12.0)') '  rank-update samples               : ', global_stats(1)
      write(output_unit,'(a,f14.0)') '  owned cell-visits                 : ', owned_total
      write(output_unit,'(a,f14.0)') '  active reactor calls              : ', reactor_calls_total
      write(output_unit,'(a,f14.0)') '  skipped cell-visits               : ', skipped_total
      write(output_unit,'(a,f10.2)') '  active fraction [%]               : ', active_fraction
      write(output_unit,'(a,f10.2)') '  skipped fraction [%]              : ', skip_fraction
      write(output_unit,'(a,f14.0)') '  skipped by temperature cutoff     : ', global_stats(6)
      write(output_unit,'(a,f14.0)') '  skipped by active-species screen  : ', global_stats(7)
      write(output_unit,'(a,f14.0)') '  skipped by legacy max-Y screen    : ', global_stats(8)
      write(output_unit,'(a,f14.0)') '  updates with all owned cells skip : ', global_stats(9)
      write(output_unit,'(a,es12.5)') '  avg reactor time per call [s]     : ', avg_reactor_time_per_call
      write(output_unit,'(a,es12.5)') '  estimated ReactorNet time avoided : ', estimated_reactor_time_avoided
      write(output_unit,'(a,es12.5)') '  accumulated chemistry wall time   : ', global_stats(10)
      write(output_unit,'(a,es12.5)') '  accumulated ReactorNet wall time  : ', global_stats(11)
      write(output_unit,'(a)') ' ======================================================================'
   end subroutine write_chemistry_screening_stats