diag_reduce_stats Subroutine

public pure subroutine diag_reduce_stats(buf, n1, n2, n3, vmin, vmax, vsum, n_valid)

Whole-array min / max / sum of a diagnostic buffer in ONE pass, over the FINITE cells only.

Replaces three separate host passes (minval/maxval/sum) with a single do concurrent ... reduce, so on a GPU build this runs where the buffer already lives and only a few scalars come back. Explicit-shape dummies (never assumed-shape) so NVHPC does not walk a descriptor per launch; index order is (k, j, i) with the contiguous index last.

Missing data. A diagnostic buffer legitimately carries IEEE NaN as the “no water here” sentinel — fill_tracer_impl writes it into every land column and every dynamically vanished layer (see its docstring, and test_fill_vanished_nan). Those cells are not data and must not enter the statistics, so every cell is ieee_is_finite-guarded and n_valid counts the cells that did contribute — the caller divides the sum by THAT, not by the array size. Relying on “comparisons with NaN are false” to make min/max skip them is not enough and not portable: it leaves sum poisoned (the whole field’s mean becomes NaN next to a perfectly finite min/max) and nvfortran’s relaxed-FP default is free to lower an unguarded min/max to a NaN-blind select.

For an all-finite buffer the guard is always taken, so the reduction order — and therefore the result — is bit-identical to the unguarded form. n_valid == 0 (nothing finite anywhere) leaves the +huge / -huge / 0 seeds untouched; the caller substitutes the missing-value sentinel.

The caller MUST only invoke this when the buffer is device-resident on a GPU build — see ocean_diag_t%on_device.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: buf(n1,n2,n3)
integer, intent(in) :: n1
integer, intent(in) :: n2
integer, intent(in) :: n3
real(kind=wp), intent(out) :: vmin
real(kind=wp), intent(out) :: vmax
real(kind=wp), intent(out) :: vsum
integer, intent(out), optional :: n_valid

Number of finite cells folded in. Optional so the pre-existing three-scalar call sites keep working unchanged.

Calls

proc~~diag_reduce_stats~~CallsGraph proc~diag_reduce_stats diag_reduce_stats reduce reduce proc~diag_reduce_stats->reduce

Called by

proc~~diag_reduce_stats~~CalledByGraph proc~diag_reduce_stats diag_reduce_stats proc~diag_field_stats diag_field_stats proc~diag_field_stats->proc~diag_reduce_stats proc~ocean_diag_step ocean_diag_t%ocean_diag_step proc~ocean_diag_step->proc~diag_field_stats proc~engine_step_finalize engine_step_finalize proc~engine_step_finalize->proc~ocean_diag_step proc~driver_run_ocean driver_run_ocean proc~driver_run_ocean->proc~engine_step_finalize proc~rdb_ocean_step rdb_ocean_step proc~rdb_ocean_step->proc~engine_step_finalize proc~driver_run driver_run proc~driver_run->proc~driver_run_ocean

Variables

Type Visibility Attributes Name Initial
integer, private :: i
integer, private :: j
integer, private :: k
integer, private :: nv

Source Code

   pure subroutine diag_reduce_stats(buf, n1, n2, n3, vmin, vmax, vsum, n_valid)
      !! Whole-array min / max / sum of a diagnostic buffer in ONE pass,
      !! over the FINITE cells only.
      !!
      !! Replaces three separate host passes (`minval`/`maxval`/`sum`) with
      !! a single `do concurrent ... reduce`, so on a GPU build this runs
      !! where the buffer already lives and only a few scalars come back.
      !! Explicit-shape dummies (never assumed-shape) so NVHPC does not walk
      !! a descriptor per launch; index order is `(k, j, i)` with the
      !! contiguous index last.
      !!
      !! **Missing data.**  A diagnostic buffer legitimately carries IEEE
      !! NaN as the "no water here" sentinel — `fill_tracer_impl` writes it
      !! into every land column and every dynamically vanished layer (see
      !! its docstring, and `test_fill_vanished_nan`).  Those cells are not
      !! data and must not enter the statistics, so every cell is
      !! `ieee_is_finite`-guarded and `n_valid` counts the cells that did
      !! contribute — the caller divides the sum by THAT, not by the array
      !! size.  Relying on "comparisons with NaN are false" to make
      !! `min`/`max` skip them is not enough and not portable: it leaves
      !! `sum` poisoned (the whole field's mean becomes NaN next to a
      !! perfectly finite min/max) and nvfortran's relaxed-FP default is
      !! free to lower an unguarded `min`/`max` to a NaN-blind select.
      !!
      !! For an all-finite buffer the guard is always taken, so the
      !! reduction order — and therefore the result — is bit-identical to
      !! the unguarded form.  `n_valid == 0` (nothing finite anywhere)
      !! leaves the `+huge` / `-huge` / `0` seeds untouched; the caller
      !! substitutes the missing-value sentinel.
      !!
      !! The caller MUST only invoke this when the buffer is device-resident
      !! on a GPU build — see `ocean_diag_t%on_device`.
      integer, intent(in) :: n1, n2, n3
      real(wp), intent(in) :: buf(n1, n2, n3)
      real(wp), intent(out) :: vmin, vmax, vsum
      integer, intent(out), optional :: n_valid
         !! Number of finite cells folded in. Optional so the pre-existing
         !! three-scalar call sites keep working unchanged.
      integer :: i, j, k, nv
      vmin = huge(1.0_wp)
      vmax = -huge(1.0_wp)
      vsum = 0.0_wp
      nv = 0
      do concurrent(k=1:n3, j=1:n2, i=1:n1) reduce(min:vmin) reduce(max:vmax) &
         reduce(+:vsum) reduce(+:nv)
         if (ieee_is_finite(buf(i, j, k))) then
            vmin = min(vmin, buf(i, j, k))
            vmax = max(vmax, buf(i, j, k))
            vsum = vsum + buf(i, j, k)
            nv = nv + 1
         end if
      end do
      if (present(n_valid)) n_valid = nv
   end subroutine diag_reduce_stats