ocean_gm_init Subroutine

   subroutine ocean_gm_init(this, grid, nz_ml)
      !! Allocate the 2D face KhTh fields, per-layer face transports, and
      !! the `gm_src` PE-release diagnostic.  Always allocates (configure
      !! runs after init); host allocation (no `do concurrent` before
      !! enter_data).
      class(ocean_gm_t), intent(inout) :: this
      type(hgrid_t), intent(in) :: grid
      integer, intent(in), optional :: nz_ml
      integer :: nx, ny, nz

      nx = grid%nx_total
      ny = grid%ny_total
      nz = 1
      if (present(nz_ml)) nz = nz_ml
      if (nz < 1) nz = 1
      ! Fail loud: the column recurrence uses NZ_STACK_MAX-sized locals.
      if (nz > NZ_STACK_MAX) then
         error stop "ocean_gm_init: nz_ml exceeds NZ_STACK_MAX "// &
            "(raise NZ_STACK_MAX in rdb_constants)"
      end if
      this%nx_total = nx
      this%ny_total = ny
      this%nz_ml = nz

      allocate (this%khth_u(nx + 1, ny), source=0.0_wp)
      allocate (this%khth_v(nx, ny + 1), source=0.0_wp)
      allocate (this%uhD(nx + 1, ny, nz), source=0.0_wp)
      allocate (this%vhD(nx, ny + 1, nz), source=0.0_wp)
      allocate (this%gm_src(nx, ny), source=0.0_wp)
      this%is_init = .true.
   end subroutine ocean_gm_init