cavity_melt_columns_2d Subroutine

public pure subroutine cavity_melt_columns_2d(nx, ny, cover, T_w, S_w, p_b, u_far, v_far, S_i, f_cor, cd, u_tide, ustar_min, par, ice, eos, const, u_star, T_b, S_b, m_mass, q_ocean, gamma_t, gamma_s, ierr_col)

Masked 2-D driver: form the friction velocity and solve the interface on every ICE-COVERED column of an (nx, ny) plane, leaving the rest untouched at exactly zero.

This routine exists for the same reason cavity_melt_columns does, and the reason is a linker one. On the GPU build the kernel lives in librdb_core.so and nvlink cannot resolve an !$acc routine seq device symbol out of a shared library into a do concurrent compiled in a DIFFERENT translation unit. So the coupling layer (rdb_ocean_cavity_flux) samples the far field and then calls THIS routine — it does not write its own column loop over cavity_ustar / cavity_melt_point. u* is folded in here for the same reason, and because it keeps ONE status per column: a non-finite far-field velocity comes back as CAVITY_MELT_NONFINITE_INPUT rather than being laundered into ustar_min and then into a small, plausible melt rate.

It is additive: cavity_melt_columns’ 1-D signature is untouched, and so are its tests.

cover is the ICE-COVER MASK (v1 binary, metrics%cover_frac, composed by the caller with the wet mask and with “the far-field sample found mass”), tested > 0.5. An UNCOVERED column is not open ocean with zero melt — it is a column the interface physics does not apply to at all — so it gets u_star = T_b = S_b = m_mass = q_ocean = 0 exactly and CAVITY_MELT_OK, and the caller’s flux components stay zero there. Counting an uncovered column as a solver failure would drown the real failures.

f_cor is a per-column ARRAY (not par%f_cor): CAVITY_LAW_HJ99 divides by |f| and takes ln(.../|f| h_nu), so on a beta plane or a spherical sector the law’s domain is a property of the COLUMN, not of the namelist. It reaches the law as a scalar (cavity_melt_point_gamma_f), and CAVITY_MELT_NO_CORIOLIS comes back for any covered column at f = 0 instead of a plausible number. S_i stays SCALAR: the ice salinity is a namelist constant in v1 (s_ice, 0 by the ISOMIP+ protocol), and a 2-D array of one repeated number is device traffic for nothing.

mem:separate contract, unchanged from the 1-D driver: this routine MOVES NOTHING. Every array must already be device-present, and the caller owns the update self of the outputs it reads on the host.

HOST routine — it LAUNCHES the kernel, so it carries no !$acc routine seq of its own.

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: nx

First dimension (ghosts included — the caller decides).
integer, intent(in) :: ny

Second dimension.
real(kind=wp), intent(in) :: cover(nx,ny)

Ice-cover mask; a column is solved iff cover > 0.5.
real(kind=wp), intent(in) :: T_w(nx,ny)

Far-field temperature (degC).
real(kind=wp), intent(in) :: S_w(nx,ny)

Far-field salinity (g/kg).
real(kind=wp), intent(in) :: p_b(nx,ny)

Interface pressure (Pa) — multilayer_state_t%p_top.
real(kind=wp), intent(in) :: u_far(nx,ny)

Far-field x velocity at the CELL CENTRE (m/s).
real(kind=wp), intent(in) :: v_far(nx,ny)

Far-field y velocity at the CELL CENTRE (m/s).
real(kind=wp), intent(in) :: S_i

Ice salinity (g/kg), >= 0; one scalar for the whole plane.
real(kind=wp), intent(in) :: f_cor(nx,ny)

Coriolis parameter (1/s) per column; read by hj99 only.
real(kind=wp), intent(in) :: cd

Top drag coefficient for the MELT friction velocity.
real(kind=wp), intent(in) :: u_tide

RMS tidal velocity (m/s); melt u* only, never the drag.
real(kind=wp), intent(in) :: ustar_min

Friction-velocity floor (m/s).
type(ocean_cavity_exchange_t), intent(in) :: par

Exchange-law bundle; its f_cor member is OVERRIDDEN per column by the f_cor array above.
type(ocean_cavity_ice_t), intent(in) :: ice

Ice-conduction bundle, shared by every column.
type(eos_t), intent(in) :: eos

Shared EOS handle — the liquidus. Flat POD, by value.
type(ocean_cavity_const_t), intent(in) :: const

Constants bundle, shared by every column.
real(kind=wp), intent(out) :: u_star(nx,ny)

Friction velocity (m/s); 0 where uncovered or refused.
real(kind=wp), intent(out) :: T_b(nx,ny)

Interface temperature (degC); 0 where uncovered.
real(kind=wp), intent(out) :: S_b(nx,ny)

Interface salinity (g/kg); 0 where uncovered.
real(kind=wp), intent(out) :: m_mass(nx,ny)

Melt mass flux (kg/m^2/s), > 0 melting; 0 where uncovered.
real(kind=wp), intent(out) :: q_ocean(nx,ny)

Ocean -> interface heat flux (W/m^2); 0 where uncovered.
real(kind=wp), intent(out) :: gamma_t(nx,ny)

Thermal exchange velocity (m/s) of each column’s converged solve; EXACTLY zero where the column is not solved.
real(kind=wp), intent(out) :: gamma_s(nx,ny)

Haline exchange velocity (m/s), same convention.
integer, intent(out) :: ierr_col(nx,ny)

CAVITY_MELT_* status per column; CAVITY_MELT_OK where uncovered.

Calls

proc~~cavity_melt_columns_2d~~CallsGraph proc~cavity_melt_columns_2d cavity_melt_columns_2d local local proc~cavity_melt_columns_2d->local proc~cavity_melt_point_gamma_f cavity_melt_point_gamma_f proc~cavity_melt_columns_2d->proc~cavity_melt_point_gamma_f proc~cavity_ustar cavity_ustar proc~cavity_melt_columns_2d->proc~cavity_ustar proc~cavity_solve_melt_f cavity_solve_melt_f proc~cavity_melt_point_gamma_f->proc~cavity_solve_melt_f proc~cavity_heat_fluxes cavity_heat_fluxes proc~cavity_solve_melt_f->proc~cavity_heat_fluxes proc~cavity_law_is_implicit cavity_law_is_implicit proc~cavity_solve_melt_f->proc~cavity_law_is_implicit proc~cavity_obukhov_length cavity_obukhov_length proc~cavity_solve_melt_f->proc~cavity_obukhov_length proc~cavity_outer_residual cavity_outer_residual proc~cavity_solve_melt_f->proc~cavity_outer_residual proc~cavity_safe_state cavity_safe_state proc~cavity_solve_melt_f->proc~cavity_safe_state proc~cavity_solution_reset cavity_solution_reset proc~cavity_solve_melt_f->proc~cavity_solution_reset proc~cavity_state_at_x cavity_state_at_x proc~cavity_solve_melt_f->proc~cavity_state_at_x proc~eos_freezing_point eos_freezing_point proc~cavity_solve_melt_f->proc~eos_freezing_point proc~cavity_ice_terms cavity_ice_terms proc~cavity_heat_fluxes->proc~cavity_ice_terms proc~cavity_t_ice cavity_t_ice proc~cavity_heat_fluxes->proc~cavity_t_ice proc~cavity_l_plus_is_neutral cavity_l_plus_is_neutral proc~cavity_outer_residual->proc~cavity_l_plus_is_neutral proc~cavity_safe_state->proc~eos_freezing_point proc~cavity_state_at_x->proc~cavity_safe_state proc~cavity_buoyancy_flux cavity_buoyancy_flux proc~cavity_state_at_x->proc~cavity_buoyancy_flux proc~cavity_exchange_velocities_f cavity_exchange_velocities_f proc~cavity_state_at_x->proc~cavity_exchange_velocities_f proc~cavity_l_plus_from_state cavity_l_plus_from_state proc~cavity_state_at_x->proc~cavity_l_plus_from_state proc~cavity_three_equation cavity_three_equation proc~cavity_state_at_x->proc~cavity_three_equation proc~cavity_gamma_hj99 cavity_gamma_hj99 proc~cavity_exchange_velocities_f->proc~cavity_gamma_hj99 proc~cavity_gamma_yung25 cavity_gamma_yung25 proc~cavity_exchange_velocities_f->proc~cavity_gamma_yung25 proc~cavity_three_equation->proc~cavity_safe_state proc~cavity_three_equation->proc~eos_freezing_point proc~cavity_three_equation->proc~cavity_ice_terms proc~cavity_three_equation->proc~cavity_t_ice

Called by

proc~~cavity_melt_columns_2d~~CalledByGraph proc~cavity_melt_columns_2d cavity_melt_columns_2d proc~ocean_cavity_flux_step ocean_cavity_flux_step proc~ocean_cavity_flux_step->proc~cavity_melt_columns_2d proc~engine_step_finalize engine_step_finalize proc~engine_step_finalize->proc~ocean_cavity_flux_step proc~driver_run_ocean driver_run_ocean proc~driver_run_ocean->proc~engine_step_finalize proc~rdb_ocean_step rdb_ocean_step proc~rdb_ocean_step->proc~engine_step_finalize proc~driver_run driver_run proc~driver_run->proc~driver_run_ocean

Variables

Type Visibility Attributes Name Initial
integer, private :: i
integer, private :: ie
integer, private :: j
real(kind=wp), private :: us

Source Code

   pure subroutine cavity_melt_columns_2d(nx, ny, cover, T_w, S_w, p_b, u_far, v_far, &
                                          S_i, f_cor, cd, u_tide, ustar_min, &
                                          par, ice, eos, const, &
                                          u_star, T_b, S_b, m_mass, q_ocean, &
                                          gamma_t, gamma_s, ierr_col)
      !! Masked 2-D driver: form the friction velocity and solve the
      !! interface on every ICE-COVERED column of an `(nx, ny)` plane,
      !! leaving the rest untouched at exactly zero.
      !!
      !! **This routine exists for the same reason `cavity_melt_columns`
      !! does, and the reason is a linker one.**  On the GPU build the
      !! kernel lives in `librdb_core.so` and nvlink cannot resolve an
      !! `!$acc routine seq` device symbol out of a shared library into a
      !! `do concurrent` compiled in a DIFFERENT translation unit.  So the
      !! coupling layer (`rdb_ocean_cavity_flux`) samples the far field
      !! and then calls THIS routine — it does not write its own column
      !! loop over `cavity_ustar` / `cavity_melt_point`.  `u*` is folded
      !! in here for the same reason, and because it keeps ONE status per
      !! column: a non-finite far-field velocity comes back as
      !! `CAVITY_MELT_NONFINITE_INPUT` rather than being laundered into
      !! `ustar_min` and then into a small, plausible melt rate.
      !!
      !! It is additive: `cavity_melt_columns`' 1-D signature is
      !! untouched, and so are its tests.
      !!
      !! `cover` is the ICE-COVER MASK (v1 binary, `metrics%cover_frac`,
      !! composed by the caller with the wet mask and with "the far-field
      !! sample found mass"), tested `> 0.5`.  An UNCOVERED column is not
      !! open ocean with zero melt — it is a column the interface physics
      !! does not apply to at all — so it gets
      !! `u_star = T_b = S_b = m_mass = q_ocean = 0` exactly and
      !! `CAVITY_MELT_OK`, and the caller's flux components stay zero
      !! there.  Counting an uncovered column as a solver failure would
      !! drown the real failures.
      !!
      !! `f_cor` is a per-column ARRAY (not `par%f_cor`):
      !! `CAVITY_LAW_HJ99` divides by `|f|` and takes `ln(.../|f| h_nu)`,
      !! so on a beta plane or a spherical sector the law's domain is a
      !! property of the COLUMN, not of the namelist.
      !! It reaches the law as a scalar (`cavity_melt_point_gamma_f`), and
      !! `CAVITY_MELT_NO_CORIOLIS` comes back for any covered column at
      !! `f = 0` instead of a plausible number.  `S_i` stays SCALAR: the
      !! ice salinity is a namelist constant in v1 (`s_ice`, 0 by the
      !! ISOMIP+ protocol), and a 2-D array of one repeated number is
      !! device traffic for nothing.
      !!
      !! `mem:separate` contract, unchanged from the 1-D driver: this
      !! routine MOVES NOTHING.  Every array must already be
      !! device-present, and the caller owns the `update self` of the
      !! outputs it reads on the host.
      !!
      !! HOST routine — it LAUNCHES the kernel, so it carries no
      !! `!$acc routine seq` of its own.
      integer, intent(in) :: nx
         !! First dimension (ghosts included — the caller decides).
      integer, intent(in) :: ny
         !! Second dimension.
      real(wp), intent(in) :: cover(nx, ny)
         !! Ice-cover mask; a column is solved iff `cover > 0.5`.
      real(wp), intent(in) :: T_w(nx, ny)
         !! Far-field temperature (degC).
      real(wp), intent(in) :: S_w(nx, ny)
         !! Far-field salinity (g/kg).
      real(wp), intent(in) :: p_b(nx, ny)
         !! Interface pressure (Pa) — `multilayer_state_t%p_top`.
      real(wp), intent(in) :: u_far(nx, ny)
         !! Far-field x velocity at the CELL CENTRE (m/s).
      real(wp), intent(in) :: v_far(nx, ny)
         !! Far-field y velocity at the CELL CENTRE (m/s).
      real(wp), intent(in) :: S_i
         !! Ice salinity (g/kg), >= 0; one scalar for the whole plane.
      real(wp), intent(in) :: f_cor(nx, ny)
         !! Coriolis parameter (1/s) per column; read by `hj99` only.
      real(wp), intent(in) :: cd
         !! Top drag coefficient for the MELT friction velocity.
      real(wp), intent(in) :: u_tide
         !! RMS tidal velocity (m/s); melt `u*` only, never the drag.
      real(wp), intent(in) :: ustar_min
         !! Friction-velocity floor (m/s).
      type(ocean_cavity_exchange_t), intent(in) :: par
         !! Exchange-law bundle; its `f_cor` member is OVERRIDDEN per
         !! column by the `f_cor` array above.
      type(ocean_cavity_ice_t), intent(in) :: ice
         !! Ice-conduction bundle, shared by every column.
      type(eos_t), intent(in) :: eos
         !! Shared EOS handle — the liquidus.  Flat POD, by value.
      type(ocean_cavity_const_t), intent(in) :: const
         !! Constants bundle, shared by every column.
      real(wp), intent(out) :: u_star(nx, ny)
         !! Friction velocity (m/s); 0 where uncovered or refused.
      real(wp), intent(out) :: T_b(nx, ny)
         !! Interface temperature (degC); 0 where uncovered.
      real(wp), intent(out) :: S_b(nx, ny)
         !! Interface salinity (g/kg); 0 where uncovered.
      real(wp), intent(out) :: m_mass(nx, ny)
         !! Melt mass flux (kg/m^2/s), > 0 melting; 0 where uncovered.
      real(wp), intent(out) :: q_ocean(nx, ny)
         !! Ocean -> interface heat flux (W/m^2); 0 where uncovered.
      real(wp), intent(out) :: gamma_t(nx, ny)
         !! Thermal exchange velocity (m/s) of each column's converged
         !! solve; EXACTLY zero where the column is not solved.
      real(wp), intent(out) :: gamma_s(nx, ny)
         !! Haline exchange velocity (m/s), same convention.
      integer, intent(out) :: ierr_col(nx, ny)
         !! `CAVITY_MELT_*` status per column; `CAVITY_MELT_OK` where
         !! uncovered.
      integer :: i, j, ie
      real(wp) :: us

      do concurrent(j=1:ny, i=1:nx) local(ie, us)
         if (cover(i, j) > 0.5_wp) then
            call cavity_ustar(u_far(i, j), v_far(i, j), cd, u_tide, ustar_min, us, ie)
            u_star(i, j) = us
            if (ie == CAVITY_MELT_OK) then
               call cavity_melt_point_gamma_f(T_w(i, j), S_w(i, j), p_b(i, j), us, S_i, &
                                              par, f_cor(i, j), ice, eos, const, &
                                              T_b(i, j), S_b(i, j), m_mass(i, j), &
                                              q_ocean(i, j), gamma_t(i, j), gamma_s(i, j), &
                                              ierr_col(i, j))
            else
               T_b(i, j) = 0.0_wp
               S_b(i, j) = 0.0_wp
               m_mass(i, j) = 0.0_wp
               q_ocean(i, j) = 0.0_wp
               gamma_t(i, j) = 0.0_wp
               gamma_s(i, j) = 0.0_wp
               ierr_col(i, j) = ie
            end if
         else
            u_star(i, j) = 0.0_wp
            T_b(i, j) = 0.0_wp
            S_b(i, j) = 0.0_wp
            m_mass(i, j) = 0.0_wp
            q_ocean(i, j) = 0.0_wp
            gamma_t(i, j) = 0.0_wp
            gamma_s(i, j) = 0.0_wp
            ierr_col(i, j) = CAVITY_MELT_OK
         end if
      end do
   end subroutine cavity_melt_columns_2d