evp_build_masks_impl Subroutine

   pure subroutine evp_build_masks_impl(wet_t, mask_t, mask_u, mask_v, mask_q, &
                                        nx_phys, ny_phys, nghost, wrap_x, wrap_y, &
                                        land_w, land_e, land_s, land_n, nx, ny)
      !! `mask_t`: `wet_T` inside the physical domain AND in every ghost
      !! band that is not pinned; a ghost band on a physical, non-periodic
      !! edge (`land_*`) is 0 (SIS2 `mask2dT` semantics — pins a wall
      !! edge's ghosts to land, matching SIS2's own domain-edge convention
      !! even when a driver run's `wet_T` ghost happens to read 1).  An MPI
      !! seam ghost follows `wet_T`, which the ocean setup exchanged; the
      !! local periodic wrap (`wrap_x/_y`, only on an axis the halo does
      !! not own) then overwrites a single-rank periodic band, so on one
      !! rank the result is the same wrap-or-0 mask as before.
      !! `mask_u(i,j) = mask_t(i-1,j)*mask_t(i,j)`, `mask_v` ditto in y,
      !! `mask_q` = product of the 4 surrounding `mask_t` (SIS2
      !! `mask2dBu`).
      integer, intent(in) :: nx_phys, ny_phys, nghost, nx, ny
      real(wp), intent(in) :: wet_t(nx, ny)
      logical, intent(in) :: wrap_x, wrap_y
      logical, intent(in) :: land_w, land_e, land_s, land_n
      real(wp), intent(out) :: mask_t(nx, ny)
      real(wp), intent(out) :: mask_u(nx + 1, ny)
      real(wp), intent(out) :: mask_v(nx, ny + 1)
      real(wp), intent(out) :: mask_q(nx + 1, ny + 1)
      integer :: i, j, i_lo, i_hi, j_lo, j_hi
      real(wp) :: mt_sw, mt_se, mt_nw, mt_ne

      i_lo = nghost + 1
      i_hi = nghost + nx_phys
      j_lo = nghost + 1
      j_hi = nghost + ny_phys

      do concurrent(j=1:ny, i=1:nx)
         if ((i < i_lo .and. land_w) .or. (i > i_hi .and. land_e) .or. &
             (j < j_lo .and. land_s) .or. (j > j_hi .and. land_n)) then
            mask_t(i, j) = 0.0_wp
         else
            mask_t(i, j) = merge(1.0_wp, 0.0_wp, wet_t(i, j) > 0.5_wp)
         end if
      end do
      call ocean_periodic_wrap_centre_2d(mask_t, nx, ny, nx_phys, ny_phys, nghost, &
                                         wrap_x, wrap_y)

      do concurrent(j=1:ny, i=1:nx + 1)
         if (i == 1 .or. i == nx + 1) then
            mask_u(i, j) = 0.0_wp
         else
            mask_u(i, j) = mask_t(i - 1, j)*mask_t(i, j)
         end if
      end do
      do concurrent(j=1:ny + 1, i=1:nx)
         if (j == 1 .or. j == ny + 1) then
            mask_v(i, j) = 0.0_wp
         else
            mask_v(i, j) = mask_t(i, j - 1)*mask_t(i, j)
         end if
      end do
      ! Interior corners, then the array edge (land) in loops of its own.
      ! NOT one loop with a four-way `i == 1 .or. i == nx+1 .or. j == 1
      ! .or. j == ny+1` guard: nvfortran 26.5's CPU vectoriser (-O2 and up)
      ! miscompiles that guard whenever nx+1 is a multiple of the vector
      ! width (AVX2: nx = 3 mod 4) -- the last vector chunk takes the
      ! interior branch at i = nx+1, so the edge column got
      ! mask_t(nx+1, j) = the next row's first cell, and at j = ny a read
      ! past the end of the array (uninitialised memory: the outermost
      ! ghost corner of the EVP stresses then differed run to run, which
      ! the restart round trip caught on a 2x2 decomposition, 19-cell
      ! tiles).  Same split in `evp_q_and_mi_ratio_impl`,
      ! `ice_limit_stresses` and `evp_str_s_relax_impl`.  Bit-identical
      ! wherever the compiler was right.
      do concurrent(j=2:ny, i=2:nx) local(mt_sw, mt_se, mt_nw, mt_ne)
         mt_sw = mask_t(i - 1, j - 1)
         mt_se = mask_t(i, j - 1)
         mt_nw = mask_t(i - 1, j)
         mt_ne = mask_t(i, j)
         mask_q(i, j) = mt_sw*mt_se*mt_nw*mt_ne
      end do
      do concurrent(i=1:nx + 1)
         mask_q(i, 1) = 0.0_wp
         mask_q(i, ny + 1) = 0.0_wp
      end do
      do concurrent(j=2:ny)
         mask_q(1, j) = 0.0_wp
         mask_q(nx + 1, j) = 0.0_wp
      end do
   end subroutine evp_build_masks_impl