seed_h_layer_uniform_impl Subroutine

   pure subroutine seed_h_layer_uniform_impl(h_layer, b, nz_ml, apply_wetdry_floor)
      !! Even-split layer thickness per column.  Pulled into a flat-impl
      !! so the `do concurrent` body works on plain allocatables — the
      !! outer shim reaches the multilayer + barotropic slot components
      !! once on the host.
      !!
      !! `apply_wetdry_floor` (wet/dry v2): with dynamic wet/dry ON the
      !! static land cutoff moves to `-land_margin`, so the intertidal band
      !! `b ∈ [-land_margin, LAND_DEPTH_THRESHOLD)` seeds LIVE (wet_mask=1).
      !! An emerged rest-bed (`b < 0`, normal for real coastlines / file
      !! bathymetry with intertidal terrain) would then get NEGATIVE layer
      !! thickness `b/nz_ml`, and since `D = Σ h_layer` (derive_bt_from_layers)
      !! that negative column sum becomes a negative barotropic depth at IC —
      !! the positive-definite limiter PRESERVES D≥0, it does not REPAIR a
      !! negative start.  When the flag is set we floor each layer to
      !! `2·H_VANISHED` (never 0 ⇒ no 1/0), so `Σ h_layer ≥ nz·2·H_VANISHED >
      !! 0` ⇒ D≥0 and `bt_eta = Σ h_layer − bt_H_ref` is consistent, and the
      !! emerged column sits at a vanished ≈0 depth ready to flood (its
      !! `wd_wet_dyn` gate seeds DRY since Σ h_layer < dry_depth).
      !!
      !! Why `2·H_VANISHED` and not `H_VANISHED`: the ALE tracer remap drain
      !! (`ocean_remap_tracer_field`) gates the concentration recovery on
      !! `h_old > H_FLOOR` with H_FLOOR == H_VANISHED and a STRICT `>`, so a
      !! layer sitting at exactly `H_VANISHED` fails the gate ⇒ its
      !! concentration is zeroed ⇒ the seeded S/T content (`hTr = const·h`)
      !! is DESTROYED on the first regrid (every outer step on sigma at the
      !! default `dt_therm_ratio=1`).  Flooring at `2·H_VANISHED` clears the
      !! strict drain gate so the seeded tracer survives — the same
      !! inflation floor `rdb_ocean_vcoord` uses for its z*/z-blend layers.
      !! Deep columns (`b/nz_ml > 2·H_VANISHED`) are untouched — only the
      !! near-zero/emerged band is floored.  Flag absent/false ⇒ literal
      !! `b/nz_ml` (byte-identical to the pre-v2 path).
      ! assumed-shape-ok: init routine called once at startup; size(b,1/2) used
      ! to derive loop bounds (flat-impl over registry-dereferenced allocatables).
      real(wp), intent(inout) :: h_layer(:, :, :)
      real(wp), intent(in)    :: b(:, :)  ! assumed-shape-ok: init routine; size(b,1/2) derives loop bounds
      integer, intent(in)    :: nz_ml
      logical, intent(in), optional :: apply_wetdry_floor
         !! When present and .true., floor each seeded layer to `2·H_VANISHED`
         !! so an emerged intertidal rest-bed never seeds a negative column
         !! (and hence a negative barotropic depth) AND survives the strict
         !! `h_old > H_VANISHED` remap-drain gate.  Absent ⇒ byte-identical.
      integer :: i, j, k, nx, ny
      real(wp) :: inv_nz
      logical :: do_floor
      nx = size(b, 1)
      ny = size(b, 2)
      inv_nz = 1.0_wp/real(nz_ml, wp)
      do_floor = .false.
      if (present(apply_wetdry_floor)) do_floor = apply_wetdry_floor
      ! Plain host loop ON PURPOSE: this runs BEFORE enter_data, so a
      ! do concurrent here makes -stdpar=gpu implicitly round-trip the
      ! (unmapped) arrays through the device per loop — measured 36.9 s
      ! of ic_seed at 5.3M cells.  Setup-only code seeds on the host;
      ! enter_data then maps the seeded values once.
      if (do_floor) then
         do k = 1, nz_ml
            do j = 1, ny
               do i = 1, nx
                  h_layer(i, j, k) = max(b(i, j)*inv_nz, 2.0_wp*H_VANISHED)
               end do
            end do
         end do
      else
         do k = 1, nz_ml
            do j = 1, ny
               do i = 1, nx
                  h_layer(i, j, k) = b(i, j)*inv_nz
               end do
            end do
         end do
      end if
   end subroutine seed_h_layer_uniform_impl