cavity_comp_apply_impl Subroutine

   pure subroutine cavity_comp_apply_impl(nx, ny, nz, dw, wet_mask, cover_frac, &
                                          h_layer, hTr_S, hTr_T, &
                                          salt_budget, heat_budget, scale, n_thin)
      !! The `volume_compensation = "uniform_open_ocean"` sink: remove
      !! `dw` metres of the top layer from every wet cell the ice does
      !! NOT cover, the removed parcel carrying that cell's own `T` and
      !! `S` so no concentration there changes.
      !!
      !! "Carries its own concentration" is implemented as a RATIO:
      !! `hTr *= h_new/h_old`, which leaves `Tr = hTr/h` algebraically
      !! identical.  The ratio is published in `scale` so every PASSIVE
      !! tracer can be given exactly the same treatment
      !! (`cavity_comp_scale_tracer_impl`) — a tracer left alone here
      !! would be CONCENTRATED by the sink, which is a different physical
      !! statement from the one the knob makes.
      !!
      !! The salt and heat increments ARE mirrored into the surface
      !! budget contributors: a parcel that carries `S` out of the domain
      !! changes the domain salt total, so the budget has to name it.
      !!
      !! Same clamp (`H_CAVITY_FLOOR`, strictly above the vanish marker)
      !! + count + fatal policy as the mass source.
      integer, intent(in) :: nx
         !! First dimension (ghosts included).
      integer, intent(in) :: ny
         !! Second dimension.
      integer, intent(in) :: nz
         !! Layer count; only `k = nz` is touched.
      real(wp), intent(in) :: dw
         !! Thickness (m) to remove from each open-ocean column.
      real(wp), intent(in) :: wet_mask(nx, ny)
         !! Static wet (1) / land (0) mask.
      real(wp), intent(in) :: cover_frac(nx, ny)
         !! Ice-cover fraction (v1 binary).
      real(wp), intent(inout) :: h_layer(nx, ny, nz)
         !! Layer thickness (m).
      real(wp), intent(inout) :: hTr_S(nx, ny, nz)
         !! `h*S`.
      real(wp), intent(inout) :: hTr_T(nx, ny, nz)
         !! `h*T`.
      real(wp), intent(inout) :: salt_budget(nx, ny, nz)
         !! `ms%salt_budget_surface`.
      real(wp), intent(inout) :: heat_budget(nx, ny, nz)
         !! `ms%heat_budget_surface`.
      real(wp), intent(out) :: scale(nx, ny)
         !! `h_new/h_old` where the sink acted, exactly 1 elsewhere.
      integer, intent(out) :: n_thin
         !! Columns whose withdrawal had to be clamped.
      integer :: i, j, c_thin
      real(wp) :: h0, hn, d, f, cell_s, cell_t

      c_thin = 0
      do concurrent(j=1:ny, i=1:nx) local(h0, hn, d, f, cell_s, cell_t) reduce(+:c_thin)
         if (wet_mask(i, j) > 0.5_wp .and. cover_frac(i, j) < 0.5_wp) then
            h0 = h_layer(i, j, nz)
            d = dw
            hn = h0 - d
            if (hn < H_CAVITY_FLOOR) then
               ! Same floor and same refusal rule as the mass source: pin
               ! the RESULT (never the increment — `h0 - (h0 - floor)`
               ! rounds and can land a ulp below), take only what sits
               ! ABOVE the floor, and take nothing at all from a column
               ! already at or below it (`min(h0, ...)`, so `hn <= h0`
               ! always and the sink can never become a deposit).
               hn = min(h0, H_CAVITY_FLOOR)
               c_thin = c_thin + 1
            end if
            if (h0 > H_DIV_EPS) then
               f = hn/h0
            else
               f = 1.0_wp
            end if
            cell_s = hTr_S(i, j, nz)*(f - 1.0_wp)
            cell_t = hTr_T(i, j, nz)*(f - 1.0_wp)
            h_layer(i, j, nz) = hn
            hTr_S(i, j, nz) = hTr_S(i, j, nz) + cell_s
            salt_budget(i, j, nz) = salt_budget(i, j, nz) + cell_s
            hTr_T(i, j, nz) = hTr_T(i, j, nz) + cell_t
            heat_budget(i, j, nz) = heat_budget(i, j, nz) + cell_t
            scale(i, j) = f
         else
            scale(i, j) = 1.0_wp
         end if
      end do
      n_thin = c_thin
   end subroutine cavity_comp_apply_impl