ocean_surface_flux_assemble Subroutine

   pure subroutine ocean_surface_flux_assemble(grid, sf, ms, active, cover_frac)
      !! **The single gate** that derives `Q_heat`/`Q_salt` from the
      !! component set (§3.1/§3.3 of the PR-12 plan) — the exact analogue
      !! of `vmix_assemble`: fillers contribute components, this routine
      !! alone derives the net fields every downstream kernel reads.  A
      !! no-op unless `sf%use_components` — with components off, `Q_heat`
      !! / `Q_salt` are exactly what `set_surface_flux_const` (or a
      !! field-override path) left them, byte-for-byte.
      !!
      !! Net surface heat into the ocean:
      !!   Q_heat = Q_heat_const + q_sw + q_lw + q_lat + q_sens + heat_added
      !!          + heat_cavity + heat_content_massin + heat_content_massout
      !!   heat_content_massin  = Σ heat_content_{lprec,fprec,vprec,
      !!                            lrunoff,frunoff,seaice_melt}
      !!   heat_content_massout = SEAWATER_CP * T_sst * evap
      !! Net surface salt flux:
      !!   Q_salt = Q_salt_const + salt_flux + salt_cavity
      !! All four outputs multiplied by `ms%wet_mask` (land carries
      !! exactly zero; interior loop bounds are NOT restricted — see
      !! CLAUDE.md's "nghost and the assembler loop bounds" gotcha).
      !!
      !! `has_heat`/`has_salt` are deliberately NOT touched here — a
      !! components-on run with no live filler must reproduce EXACTLY
      !! what `set_surface_flux_const` left them (the "+0.0 bit-identity
      !! trap": forcing them true would change `apply_surface_src_2d_impl`
      !! from an early-return to a `+0.0` stamp, altering a budget sum's
      !! operand COUNT even though the value is unchanged).  Setting
      !! `has_heat`/`has_salt` is the FILLER's job (see the module
      !! docstring's fill contract).
      !!
      !! Outer-shim + flat-impl (`ocean_surfflux_assemble_impl`): the SST
      !! read needs `ms%tracers(idx_temperature)%hTr`, an array-of-DT
      !! registry deref that must happen on the host before the `do
      !! concurrent` (NVHPC device codegen constraint — see
      !! `ocean_surface_flux_apply_tracers` for the identical pattern).
      !! No-op when no temperature tracer is registered (SST is
      !! undefined without one).
      type(hgrid_t), intent(in) :: grid
      type(ocean_surface_flux_t), intent(inout) :: sf
      type(multilayer_state_t), intent(in) :: ms
      logical, intent(in), optional :: active
         !! Optional thermo-cadence gate.  Present-and-false ⇒ early
         !! return; absent ⇒ kernel runs (matches
         !! `ocean_surface_flux_apply_tracers`'s convention).
      real(wp), intent(in), optional :: cover_frac(:, :)
         !! Optional ice-shelf cover fraction (`metrics%cover_frac`,
         !! v1 binary).  Present ⇒ every ATMOSPHERIC contribution
         !! (`Q_heat_const`, `q_sw`, `q_lw`, `q_lat`, `q_sens`,
         !! `heat_added`, both mass-enthalpy terms, `Q_salt_const`,
         !! `salt_flux`) is scaled by `1 - cover_frac`, while the
         !! cavity's OWN `heat_cavity` / `salt_cavity` pass through
         !! unmasked.  This is the single place those two groups are
         !! still distinguishable — see the module docstring for why the
         !! mask lives here and not at apply time.  Absent ⇒ the
         !! original kernel, byte-identical.
      ! assumed-shape-ok: thermo-cadence shim, forwarded to an
      ! explicit-shape `_impl` before the device loop.

      integer :: nx, ny, nz, idx_T
      logical :: masked

      if (.not. sf%use_components) return
      if (present(active)) then
         if (.not. active) return
      end if
      if (.not. allocated(ms%tracers)) return
      idx_T = ms%idx_temperature
      if (idx_T <= 0) return

      nx = grid%nx_total
      ny = grid%ny_total
      nz = ms%nz_ml

      masked = .false.
      if (present(cover_frac)) then
         masked = (size(cover_frac, 1) == nx .and. size(cover_frac, 2) == ny)
      end if
      if (masked) then
         call ocean_surfflux_assemble_cover_impl( &
            sf%heat_content_massin, sf%heat_content_massout, sf%Q_heat, sf%Q_salt, &
            sf%q_sw, sf%q_lw, sf%q_lat, sf%q_sens, sf%heat_added, sf%heat_cavity, &
            sf%heat_content_lprec, sf%heat_content_fprec, sf%heat_content_vprec, &
            sf%heat_content_lrunoff, sf%heat_content_frunoff, sf%heat_content_seaice_melt, &
            sf%evap, sf%salt_flux, sf%salt_cavity, &
            ms%tracers(idx_T)%hTr, ms%h_layer, ms%wet_mask, cover_frac, &
            sf%Q_heat_const, sf%Q_salt_const, sf%cp, sf%h_min, nz, nx, ny)
         sf%q_assembled = 1.0_wp
         return
      end if

      call ocean_surfflux_assemble_impl( &
         sf%heat_content_massin, sf%heat_content_massout, sf%Q_heat, sf%Q_salt, &
         sf%q_sw, sf%q_lw, sf%q_lat, sf%q_sens, sf%heat_added, sf%heat_cavity, &
         sf%heat_content_lprec, sf%heat_content_fprec, sf%heat_content_vprec, &
         sf%heat_content_lrunoff, sf%heat_content_frunoff, sf%heat_content_seaice_melt, &
         sf%evap, sf%salt_flux, sf%salt_cavity, &
         ms%tracers(idx_T)%hTr, ms%h_layer, ms%wet_mask, &
         sf%Q_heat_const, sf%Q_salt_const, sf%cp, sf%h_min, nz, nx, ny)
      sf%q_assembled = 1.0_wp
   end subroutine ocean_surface_flux_assemble