ocean_console_stats_report Subroutine

   subroutine ocean_console_stats_report(this, grid, metrics, ms, t, dt, step, &
                                         horiz_adv_budget_valid, &
                                         cfl_vanish_tol, heat_budget_frazil, &
                                         ice_part_size, ice_m_ice, ice_ncat, &
                                         compute_rank, reproducing_sums, &
                                         budget_stage_weight, budget_out)
      !! Compute current totals + means + max-CFL and emit a MOM6-style
      !! console block via the shared `console_stats_report` formatter.
      !!
      !! COLLECTIVE: must be called by ALL compute ranks.  Per-rank local
      !! sums are reduced globally (`halo_allreduce_sum` / `_max`) before any
      !! derived quantity is formed; the shared formatter then prints on rank
      !! 0 only (`is_root`) while its NaN / CFL panic `error stop` stays
      !! collective on every rank.  Single-rank allreduce is an identity ⇒
      !! bit-identical to the serial path.
      !!
      !! Reductions run on device; tracer registry indirection is dereferenced
      !! on the host shim before each flat-impl helper.
      !!
      !! Optional `horiz_adv_budget_valid`:
      !! salt/heat closed-budget fall-back gate (see
      !! `ocean_budget_is_active`).  Absent ⇒ budget active / gate off.
      !! (PR-23 removed the `sponge_relax_tracers` fall-back — both sponge
      !! paths' tracer relaxation are now instrumented into
      !! `ms%salt_budget_sponge`/`heat_budget_sponge`, folded by
      !! `ocean_salt_src_sum`/`ocean_heat_src_sum`, so the gate was dead.)
      !!
      !! Phase-3 optional `cfl_vanish_tol`: when present (> 0) the LOCAL MaxCFL
      !! reduction excludes cells with `h_layer <= vanish_tol` before the
      !! global `allreduce_max`.  Absent ⇒ un-gated (bit-identical).
      type(console_stats_t), intent(inout) :: this
      type(hgrid_t), intent(in) :: grid
      type(ocean_metrics_t), intent(in) :: metrics
      type(multilayer_state_t), intent(in) :: ms
      real(wp), intent(in) :: t, dt
      integer, intent(in) :: step
      logical, intent(in), optional :: horiz_adv_budget_valid
         !! When `.false.` the horizontal-advection accumulator is incomplete
         !! (windowed path); revert to raw drift for salt/heat.
      real(wp), intent(in), optional :: cfl_vanish_tol
         !! When present: vanish-gated MaxCFL (Phase 3). Absent ⇒ un-gated.
      real(wp), intent(in), optional :: heat_budget_frazil(:, :, :)
         !! Sea-ice frazil heat-budget accumulator (`ocean_sea_ice_t%
         !! heat_budget_frazil`, K·m per cell) — folded into `heat_src`
         !! at FULL weight (`ocean_frazil_heat_src`) so the closed heat
         !! budget still closes with ice on.  Absent (ice off — the
         !! driver passes the slot's unallocated array, which an optional
         !! dummy sees as absent) ⇒ bit-identical.
      real(wp), intent(in), optional :: ice_part_size(:, :, 0:)
         !! Sea-ice category area fractions (`ocean_sea_ice_t%part_size`,
         !! category 0 = open water).  Absent (ice off — the driver passes
         !! the slot's unallocated array, seen as absent) ⇒ no ice line ⇒
         !! bit-identical.  Lower bound 0 declared so category indexing
         !! matches the state array.
      real(wp), intent(in), optional :: ice_m_ice(:, :, :)
         !! Sea-ice mass per category (two-mode convention, `rdb_ice_state`).
      integer, intent(in), optional :: ice_ncat
         !! Category count (selects the lumped vs ITD gather mode).
      integer, intent(in) :: compute_rank
         !! Caller's compute-communicator rank (0 = print, others silent).
      logical, intent(in), optional :: reproducing_sums
         !! PR-32: `.true.` routes the `total_h`/`raw_ke`/`raw_heat`/
         !! `raw_salt`/`raw_age`/ice-area reductions and their cross-rank
         !! combine through the order-invariant EFP path (`rdb_efp` +
         !! `halo_allreduce_efp_list`, ONE collective) instead of the FP
         !! `!$acc parallel loop reduction(+:acc)` kernels + SEVEN separate
         !! `halo_allreduce_sum` calls.  Absent / `.false.` (default) ⇒ the
         !! FP path below runs VERBATIM — byte-identical console output.
         !! The salt/heat closed-budget out/src terms (section (d) below)
         !! take the EFP path too (one extra collective per active budget),
         !! so every printed number is independent of the rank count; the
         !! primary totals feeding the `Error` residual also get the
         !! `efp_real_diff` treatment.  The cumulative `mass_out` is
         !! accumulated per step in EFP bins too (`mass_out_efp_on`, set
         !! from this knob); the cavity-only `mass_src` stays an FP running
         !! sum (single-rank path).
      real(wp), intent(in), optional :: budget_stage_weight
         !! Per-outer-step weight for the salt/heat budget accumulators, from
         !! `ocean_budget_stage_weight(is_pc)`.  Absent ⇒ the historical
         !! SSP-RK2 0.5 ⇒ bit-identical.  MUST be 1.0 under
         !! `&ocean_bt_nml split_scheme="pred_corr"`, whose single prognostic
         !! tracer update fills the accumulators once per step rather than
         !! twice; without it the reported `src`/`out` are exactly half and
         !! the closed-budget residual reads ~5e-4 instead of ~1e-13.
      type(conservation_budget_t), intent(out), optional :: budget_out
         !! Hands the computed closed-budget totals (`mass_out`/`salt_out`/
         !! `heat_out`/`*_src`/`*_active`, section (d) below) back to the
         !! caller instead of only formatting them into the printed console
         !! line.  Absent ⇒ inert (no behaviour change; existing callers are
         !! untouched).  Consumer: `tests/mpi/test_ocean_decomp_bitid_mpi`,
         !! which asserts these are bit-identical across every `px x py`
         !! factorisation of each case (`compare_budget`), and, on the
         !! `periodic_sponge` case, against the serial reference too
         !! (`check_periodic_sponge_serial_out`) — the EFP `reproducing_sums` path (default on) makes the boundary
         !! `out` terms order-invariant across rank counts, exactly like the
         !! totals.

      ! (a) per-rank local sums; (b) global after allreduce; (c) derive.
      real(wp) :: total_h, raw_ke, raw_heat, raw_salt, raw_age, max_cfl
      real(wp) :: total_mass, total_ke, total_heat, total_salt
      real(wp) :: mean_S, mean_T, mean_age
      real(wp) :: tmp
      logical :: hav
      type(conservation_budget_t) :: bud
      real(wp) :: b_salt_surf, b_salt_adv, b_salt_hdiff, b_salt_sponge
      real(wp) :: b_heat_surf, b_heat_geo, b_heat_adv, b_heat_hdiff, b_heat_sponge
      real(wp) :: b_heat_frazil
      real(wp) :: bud_w
      logical :: ice_on, use_efp
      real(wp) :: l_wet_area, l_ci_area, l_hi_area
      real(wp) :: g_wet_area, g_ci_area, g_hi_area
      real(wp) :: mean_ci, mean_hi
      character(len=256) :: ice_line
      real(wp) :: g_mass_out, g_mass_src
      ! PR-32 EFP local/global lists — fixed NVAL=9 layout regardless of
      ! gating (ice-off slots stay zero-valued efp_t) so the collective's
      ! size never depends on a per-rank branch.
      integer, parameter :: NVAL_EFP = 10
      integer, parameter :: IX_H = 1, IX_KE = 2, IX_HEAT = 3, IX_SALT = 4, IX_AGE = 5
      integer, parameter :: IX_WET = 6, IX_CI = 7, IX_HI = 8, IX_MOUT = 9
      integer, parameter :: IX_MSRC = 10
      type(efp_t) :: efp_local(NVAL_EFP), efp_global(NVAL_EFP)
      type(efp_t) :: mass_efp_v, salt_efp_v, heat_efp_v
      integer, parameter :: NVAL_BUD = 6
         !! Closed-budget terms combined in ONE EFP collective per tracer
         !! (salt uses 4 slots, heat 5 + frazil); unused slots stay zero.
      type(efp_t) :: efp_bl(NVAL_BUD), efp_bg(NVAL_BUD)

      use_efp = .false.
      if (present(reproducing_sums)) use_efp = reproducing_sums

      ! ---- (a) per-rank local reductions over the subdomain ---------------
      ! Physical integrals are area-weighted: `compute_total_*` fold areaT so
      ! `total_h` is Σ h·areaT [m³] (ghosts excluded ⇒ no double-count on seam).
      ! Phase-3 vanish gate applied to the LOCAL max before allreduce_max —
      ! MaxCFL is untouched by `reproducing_sums` in EITHER branch (max is
      ! already exact/order-invariant in FP; SS2.4 of the plan).
      if (present(cfl_vanish_tol)) then
         max_cfl = compute_max_cfl(ms%u_face_x_layer, ms%v_face_y_layer, &
                                   metrics%idxT, metrics%idyT, dt, grid%nghost, &
                                   ms%h_layer, cfl_vanish_tol)
      else
         max_cfl = compute_max_cfl(ms%u_face_x_layer, ms%v_face_y_layer, &
                                   metrics%idxT, metrics%idyT, dt, grid%nghost)
      end if

      ! Sea-ice conc/thickness area sums (PR ice-diags). `ice_on` is uniform
      ! across ranks (the ice arrays are allocated on every rank iff
      ! ice%enable) so gating the three extra allreduces below on it is
      ! collective-safe.  Pre-zero the global accumulators so the ice-off
      ! path's later `if (ice_on .and. g_wet_area > ...)` read is defined
      ! even if `.and.` does not short-circuit.
      g_wet_area = 0.0_wp
      g_ci_area = 0.0_wp
      g_hi_area = 0.0_wp
      ice_on = present(ice_part_size) .and. present(ice_m_ice) .and. present(ice_ncat)

      if (use_efp) then
         ! ---- EFP path: k-slab-blocked fixed-point reductions, ONE combined
         ! collective (`halo_allreduce_efp_list`) in place of the seven
         ! `halo_allreduce_sum` calls the FP branch below issues for these
         ! same quantities.
         efp_local = efp_t()
         efp_local(IX_H) = compute_total_h_efp(ms%h_layer, metrics%areaT, grid%nghost)
         efp_local(IX_KE) = compute_total_ke_efp(ms%h_layer, ms%u_face_x_layer, &
                                                 ms%v_face_y_layer, metrics%areaT, grid%nghost)
         if (ms%idx_temperature > 0) then
            efp_local(IX_HEAT) = compute_total_tracer_efp(ms%tracers(ms%idx_temperature)%hTr, &
                                                          metrics%areaT, grid%nghost)
         end if
         if (ms%idx_salinity > 0) then
            efp_local(IX_SALT) = compute_total_tracer_efp(ms%tracers(ms%idx_salinity)%hTr, &
                                                          metrics%areaT, grid%nghost)
         end if
         if (ms%idx_age > 0) then
            efp_local(IX_AGE) = compute_total_tracer_efp(ms%tracers(ms%idx_age)%hTr, &
                                                         metrics%areaT, grid%nghost)
         end if
         if (ice_on) then
            call compute_ice_totals_efp(metrics%wet_T, metrics%areaT, ice_part_size, ice_m_ice, &
                                        ice_ncat, grid%nghost, &
                                        efp_local(IX_WET), efp_local(IX_CI), efp_local(IX_HI))
         end if
         if (ms%mass_out_efp_on) then
            efp_local(IX_MOUT)%v = ms%mass_out_efp
            efp_local(IX_MOUT)%poison = ms%mass_out_efp_poison
         else
            efp_local(IX_MOUT) = efp_from_real(real(ms%mass_out, real64))
         end if
         efp_local(IX_MSRC) = efp_from_real(real(ms%mass_src, real64))

         call halo_allreduce_efp_list(efp_local, efp_global, NVAL_EFP)

         total_h = real(efp_to_real(efp_global(IX_H)), wp)
         raw_ke = real(efp_to_real(efp_global(IX_KE)), wp)
         raw_heat = 0.0_wp
         if (ms%idx_temperature > 0) raw_heat = real(efp_to_real(efp_global(IX_HEAT)), wp)
         raw_salt = 0.0_wp
         if (ms%idx_salinity > 0) raw_salt = real(efp_to_real(efp_global(IX_SALT)), wp)
         raw_age = 0.0_wp
         if (ms%idx_age > 0) raw_age = real(efp_to_real(efp_global(IX_AGE)), wp)
         if (ice_on) then
            g_wet_area = real(efp_to_real(efp_global(IX_WET)), wp)
            g_ci_area = real(efp_to_real(efp_global(IX_CI)), wp)
            g_hi_area = real(efp_to_real(efp_global(IX_HI)), wp)
         end if
         g_mass_out = real(efp_to_real(efp_global(IX_MOUT)), wp)
         g_mass_src = real(efp_to_real(efp_global(IX_MSRC)), wp)

         tmp = max_cfl
         call halo_allreduce_max(tmp, max_cfl)
      else
         ! ---- FP path — byte-identical to pre-PR-32 -----------------------
         total_h = compute_total_h(ms%h_layer, metrics%areaT, grid%nghost)
         raw_ke = compute_total_ke(ms%h_layer, ms%u_face_x_layer, &
                                   ms%v_face_y_layer, metrics%areaT, grid%nghost)
         if (ms%idx_temperature > 0) then
            raw_heat = compute_total_tracer(ms%tracers(ms%idx_temperature)%hTr, &
                                            metrics%areaT, grid%nghost)
         else
            raw_heat = 0.0_wp
         end if
         if (ms%idx_salinity > 0) then
            raw_salt = compute_total_tracer(ms%tracers(ms%idx_salinity)%hTr, &
                                            metrics%areaT, grid%nghost)
         else
            raw_salt = 0.0_wp
         end if
         ! Ideal-age tracer: passive, reported whenever registered (volume-mean).
         if (ms%idx_age > 0) then
            raw_age = compute_total_tracer(ms%tracers(ms%idx_age)%hTr, &
                                           metrics%areaT, grid%nghost)
         else
            raw_age = 0.0_wp
         end if
         if (ice_on) then
            call compute_ice_totals(metrics%wet_T, metrics%areaT, ice_part_size, ice_m_ice, &
                                    ice_ncat, grid%nghost, l_wet_area, l_ci_area, l_hi_area)
         end if

         ! ---- (b) global allreduce — COLLECTIVE (all compute ranks) -------
         ! 1-rank ⇒ identities. Distinct in/out temps avoid aliasing.
         tmp = total_h
         call halo_allreduce_sum(tmp, total_h)
         tmp = raw_ke
         call halo_allreduce_sum(tmp, raw_ke)
         tmp = raw_heat
         call halo_allreduce_sum(tmp, raw_heat)
         tmp = raw_salt
         call halo_allreduce_sum(tmp, raw_salt)
         tmp = raw_age
         call halo_allreduce_sum(tmp, raw_age)
         tmp = max_cfl
         call halo_allreduce_max(tmp, max_cfl)
         if (ice_on) then
            tmp = l_wet_area
            call halo_allreduce_sum(tmp, g_wet_area)
            tmp = l_ci_area
            call halo_allreduce_sum(tmp, g_ci_area)
            tmp = l_hi_area
            call halo_allreduce_sum(tmp, g_hi_area)
         end if
         ! `g_mass_out` is computed in section (d) below on this branch (as
         ! before PR-32); the EFP branch above computes it in the combined
         ! list instead, so guard against the section-(d) allreduce running
         ! twice.
      end if

      ! ---- (c) derive conserved totals + means from the GLOBAL sums -------
      total_mass = total_h*RHO_WATER
      total_ke = raw_ke*RHO_WATER
      if (ms%idx_temperature > 0) then
         if (total_h > 0.0_wp) then
            mean_T = raw_heat/total_h
         else
            mean_T = 0.0_wp
         end if
         total_heat = raw_heat*RHO_WATER
      else
         total_heat = 0.0_wp
         mean_T = 0.0_wp
      end if
      if (ms%idx_salinity > 0) then
         if (total_h > 0.0_wp) then
            mean_S = raw_salt/total_h
         else
            mean_S = 0.0_wp
         end if
         total_salt = raw_salt*RHO_WATER
      else
         total_salt = 0.0_wp
         mean_S = 0.0_wp
      end if
      if (ms%idx_age > 0 .and. total_h > 0.0_wp) then
         mean_age = raw_age/total_h
      else
         mean_age = 0.0_wp
      end if
      mean_ci = 0.0_wp
      mean_hi = 0.0_wp
      ! `.and.` is not guaranteed to short-circuit in Fortran, so the
      ! ice-off path could read the g_* accumulators (assigned only inside
      ! the `if (ice_on)` allreduce block); they are pre-zeroed at the top
      ! of section (a) to keep this read defined.
      if (ice_on .and. g_wet_area > tiny(0.0_wp)) then
         mean_ci = g_ci_area/g_wet_area
         mean_hi = g_hi_area/g_wet_area
      end if

      ! ---- (d) closed salt/heat/mass budget — allreduced for multi-rank ---
      hav = .true.
      if (present(horiz_adv_budget_valid)) hav = horiz_adv_budget_valid

      ! Mass out: cumulative open-boundary volume out (scalar accumulator).
      ! `g_mass_out` is already globally combined above on the EFP branch
      ! (packed into the same `halo_allreduce_efp_list` call as the totals);
      ! the FP branch combines it here, as before PR-32.
      if (.not. use_efp) then
         tmp = ms%mass_out
         call halo_allreduce_sum(tmp, g_mass_out)
         tmp = ms%mass_src
         call halo_allreduce_sum(tmp, g_mass_src)
      end if
      bud%mass_out = g_mass_out
      ! Tracked mass SOURCE — the mass twin of `salt_src`/`heat_src`.
      ! Zero on every path but the ice-shelf real-freshwater one
      ! (`&ocean_cavity_melt_nml freshwater="mass"`) and its
      ! `volume_compensation` sink, so the printed budget is unchanged
      ! elsewhere.  Already weighted per RK2 stage at accumulation time,
      ! exactly like `mass_out`, so no `bud_w` appears here.
      bud%mass_src = g_mass_src
      bud%mass_active = ms%mass_out_tracked

      ! Salt / heat closed budget: local per-cell integrals, allreduced, then
      ! the weight+sign+ρ scaling via the pure `ocean_budget_*` helpers.
      ! PR-23: `b_salt_sponge` / `b_heat_sponge` fold in alongside
      ! `b_heat_geo` — zero unless `&ocean_sponge_nml enable=.true.,
      ! relax_tracers=.true.` OR the legacy `sponge_relax_tracers=.true.`.
      bud_w = RK2_STAGE_WEIGHT
      if (present(budget_stage_weight)) bud_w = budget_stage_weight
      if (ocean_budget_is_active(ms%idx_salinity, hav)) then
         if (use_efp) then
            ! Order-invariant budget terms too, so the `out` / `src` columns
            ! print the same digits on every rank count (they used to stay
            ! on the FP path even with `reproducing_sums`).
            efp_bl = efp_t()
            efp_bl(1) = compute_total_tracer_efp(ms%salt_budget_surface, metrics%areaT, grid%nghost)
            efp_bl(2) = compute_total_tracer_efp(ms%salt_budget_sponge, metrics%areaT, grid%nghost)
            efp_bl(3) = compute_total_tracer_efp(ms%salt_budget_horiz_adv, metrics%areaT, grid%nghost)
            efp_bl(4) = compute_total_tracer_efp(ms%salt_budget_hdiff, metrics%areaT, grid%nghost)
            call halo_allreduce_efp_list(efp_bl, efp_bg, NVAL_BUD)
            b_salt_surf = real(efp_to_real(efp_bg(1)), wp)
            b_salt_sponge = real(efp_to_real(efp_bg(2)), wp)
            b_salt_adv = real(efp_to_real(efp_bg(3)), wp)
            b_salt_hdiff = real(efp_to_real(efp_bg(4)), wp)
         else
            b_salt_surf = compute_total_tracer(ms%salt_budget_surface, metrics%areaT, grid%nghost)
            b_salt_sponge = compute_total_tracer(ms%salt_budget_sponge, metrics%areaT, grid%nghost)
            b_salt_adv = compute_total_tracer(ms%salt_budget_horiz_adv, metrics%areaT, grid%nghost)
            b_salt_hdiff = compute_total_tracer(ms%salt_budget_hdiff, metrics%areaT, grid%nghost)
            tmp = b_salt_surf
            call halo_allreduce_sum(tmp, b_salt_surf)
            tmp = b_salt_sponge
            call halo_allreduce_sum(tmp, b_salt_sponge)
            tmp = b_salt_adv
            call halo_allreduce_sum(tmp, b_salt_adv)
            tmp = b_salt_hdiff
            call halo_allreduce_sum(tmp, b_salt_hdiff)
         end if
         bud%salt_src = ocean_budget_src(ocean_salt_src_sum(b_salt_surf, b_salt_sponge), &
                                         stage_weight=bud_w)
         bud%salt_out = ocean_budget_out(b_salt_adv, b_salt_hdiff, stage_weight=bud_w)
         bud%salt_active = .true.
      end if
      if (ocean_budget_is_active(ms%idx_temperature, hav)) then
         ! Sea-ice frazil source (PR 1) — full weight, see
         ! `ocean_frazil_heat_src`.  Absent / ice-off ⇒ adds 0.
         b_heat_frazil = 0.0_wp
         if (use_efp) then
            efp_bl = efp_t()
            efp_bl(1) = compute_total_tracer_efp(ms%heat_budget_surface, metrics%areaT, grid%nghost)
            efp_bl(2) = compute_total_tracer_efp(ms%heat_budget_geothermal, metrics%areaT, grid%nghost)
            efp_bl(3) = compute_total_tracer_efp(ms%heat_budget_sponge, metrics%areaT, grid%nghost)
            efp_bl(4) = compute_total_tracer_efp(ms%heat_budget_horiz_adv, metrics%areaT, grid%nghost)
            efp_bl(5) = compute_total_tracer_efp(ms%heat_budget_hdiff, metrics%areaT, grid%nghost)
            if (present(heat_budget_frazil)) then
               efp_bl(6) = compute_total_tracer_efp(heat_budget_frazil, metrics%areaT, grid%nghost)
            end if
            call halo_allreduce_efp_list(efp_bl, efp_bg, NVAL_BUD)
            b_heat_surf = real(efp_to_real(efp_bg(1)), wp)
            b_heat_geo = real(efp_to_real(efp_bg(2)), wp)
            b_heat_sponge = real(efp_to_real(efp_bg(3)), wp)
            b_heat_adv = real(efp_to_real(efp_bg(4)), wp)
            b_heat_hdiff = real(efp_to_real(efp_bg(5)), wp)
            if (present(heat_budget_frazil)) b_heat_frazil = real(efp_to_real(efp_bg(6)), wp)
         else
            b_heat_surf = compute_total_tracer(ms%heat_budget_surface, metrics%areaT, grid%nghost)
            b_heat_geo = compute_total_tracer(ms%heat_budget_geothermal, metrics%areaT, grid%nghost)
            b_heat_sponge = compute_total_tracer(ms%heat_budget_sponge, metrics%areaT, grid%nghost)
            b_heat_adv = compute_total_tracer(ms%heat_budget_horiz_adv, metrics%areaT, grid%nghost)
            b_heat_hdiff = compute_total_tracer(ms%heat_budget_hdiff, metrics%areaT, grid%nghost)
            tmp = b_heat_surf
            call halo_allreduce_sum(tmp, b_heat_surf)
            tmp = b_heat_geo
            call halo_allreduce_sum(tmp, b_heat_geo)
            tmp = b_heat_sponge
            call halo_allreduce_sum(tmp, b_heat_sponge)
            tmp = b_heat_adv
            call halo_allreduce_sum(tmp, b_heat_adv)
            tmp = b_heat_hdiff
            call halo_allreduce_sum(tmp, b_heat_hdiff)
            if (present(heat_budget_frazil)) then
               b_heat_frazil = compute_total_tracer(heat_budget_frazil, &
                                                    metrics%areaT, grid%nghost)
               tmp = b_heat_frazil
               call halo_allreduce_sum(tmp, b_heat_frazil)
            end if
         end if
         bud%heat_src = ocean_budget_src(ocean_heat_src_sum(b_heat_surf, b_heat_geo, &
                                                            b_heat_sponge), &
                                         stage_weight=bud_w) &
                        + ocean_frazil_heat_src(b_heat_frazil)
         bud%heat_out = ocean_budget_out(b_heat_adv, b_heat_hdiff, stage_weight=bud_w)
         bud%heat_active = .true.
      end if

      if (present(budget_out)) budget_out = bud

      ! ---- (e) hand off to the shared MOM6-style formatter ----------------
      ! COLLECTIVE: all ranks call it (its NaN/CFL panic error-stops on every
      ! rank); `is_root` gates the printed lines to rank 0.
      !
      ! PR-32 SS2.2 fix: the EFP path builds `mass_efp`/`salt_efp`/`heat_efp`
      ! from a SINGLE scalar multiply (`total_mass = total_h*RHO_WATER`,
      ! already formed above) re-decomposed via `efp_from_real` — a lone
      ! multiplication has no accumulation-order issue, so plain FP is fine
      ! for this step; only the SUMMATION that produced `total_h` itself
      ! needed the fixed-point treatment.  `console_stats_report` then
      ! latches these as the EFP reference and forms the `Error` residual
      ! via `efp_real_diff` against it — a difference in FIXED POINT, not a
      ! double subtraction of two already-quantised totals.
      if (use_efp) then
         mass_efp_v = efp_from_real(real(total_mass, real64))
         salt_efp_v = efp_from_real(real(total_salt, real64))
         heat_efp_v = efp_from_real(real(total_heat, real64))
         call console_stats_report(this, t, step, total_mass, total_ke, &
                                   mean_S, mean_T, max_cfl, total_salt, total_heat, &
                                   has_salt=(ms%idx_salinity > 0), &
                                   has_temp=(ms%idx_temperature > 0), &
                                   mean_age=mean_age, has_age=(ms%idx_age > 0), &
                                   budget=bud, is_root=(compute_rank == 0), &
                                   mass_efp=mass_efp_v, salt_efp=salt_efp_v, heat_efp=heat_efp_v)
      else
         call console_stats_report(this, t, step, total_mass, total_ke, &
                                   mean_S, mean_T, max_cfl, total_salt, total_heat, &
                                   has_salt=(ms%idx_salinity > 0), &
                                   has_temp=(ms%idx_temperature > 0), &
                                   mean_age=mean_age, has_age=(ms%idx_age > 0), &
                                   budget=bud, is_root=(compute_rank == 0))
      end if

      ! Sea-ice line — ocean-area-weighted mean conc/thickness (PR ice-diags).
      ! Absent ice arrays (ice off) ⇒ ice_on = .false. ⇒ no line ⇒
      ! bit-identical console output.
      if (ice_on .and. compute_rank == 0) then
         write (ice_line, "('    Ice  : conc ',F7.4,'  thick ',F9.4,' m  (ocean-area mean)')") &
            mean_ci, mean_hi
         call logger%info(trim(ice_line))
      end if

      ! Halo-exchange counter line — cumulative semantic exchange counts,
      ! rank 0 only (the shared formatter owns the physics lines above).
      if (compute_rank == 0) call logger%info("    "//trim(oh_counters_format()))
   end subroutine ocean_console_stats_report