rdb_ice_evp Module

Mechanical transliteration of the validated Python prototype tmp_local_artifacts/proto_evp_core.py (Channel1D) + tmp_local_artifacts/proto_evp_validate.py (Box2D), grounded against SIS_dyn_cgrid.F90 (SIS_C_dynamics :603-1614, limit_stresses :1619-1741, SIS_C_dyn_init :209-270) per SPEC_ice-pr5-evp.md. Field names mirror SIS2 1:1 (str_d, str_t, str_s, sh_Dd, sh_Dt, sh_Ds, zeta, del_sh, mi_ratio_A_q, Tdamp, EC, …).

Index convention (SPEC §1): rdb u-face (i,j) is the WEST face of cell (i,j) (SIS2 face I = east of cell i); rdb v-face (i,j) is the SOUTH face; rdb corner (i,j) is the SW corner of cell (i,j) (SIS2 corner (I,J) = NE of cell (i,j)). The 4 T-cells around rdb corner (ic,jc) are (ic-1,jc-1) (ic,jc-1) (ic-1,jc) (ic,jc) — the wet_q convention.

Ice margins (Lens C, critical — no ice-edge code path anywhere). mi=0/ci=0 cells are ORDINARY wet T-cells: pres_mice*mice=0 => zeta=0 => str_d decays geometrically toward 0 via the I_1pdt_T relaxation — emergent, never a Dirichlet special case. Masks built here (mask_t_w/mask_u_w/mask_v_w/mask_q_w) encode LAND + non-periodic-boundary policy ONLY; they never test ice presence. The momentum-solve denominator’s 0/0 guard is m_neglect; it is not the module’s only division guard — dxharm>0, denom/=0, and the Adcroft i_htot reciprocal each guard their own quotient. With PR 62 a_face_stress=.true. the guard alone is NOT sufficient at an ice-free face: drag_eff = a_u*drag_u is explicitly BRANCHED around (a_fac > 0.0, not floored) in evp_u_momentum_impl/ evp_v_momentum_impl, because m_neglect alone against a generally nonzero dt*fxic_now blows up to O(1e30) on the first substep.

Persistent workspace (never local-allocate scratch in a per-substep kernel on -stdpar=gpu): the EVP scratch lives on the evp_workspace_t slot (ice%evp_ws, in rdb_ice_state), eagerly allocated in ocean_sea_ice_t%init and GPU-mapped via ocean_sea_ice_t%enter_data (which rides ocean_state_enter_data) — matching the rest of src/core/ocean/ (zero module-level save allocatables). ice_evp_dynamics is a thin shim that unpacks the slot’s components into the explicit-shape flat-impl args of ice_evp_dynamics_impl; the whole substep loop is device-resident: no H<->D inside do n = 1, evp_sub_steps.

Documented divergences (D-list, SPEC §4.8): D1 no sea-surface-tilt term (SIS2 PFu/PFv) — v1 assumes flat eta for the ice; a future PR wires the ocean SSH. D2 CLOSED by PR 36 — PROJECT_ICE_CONCENTRATION is ported as &ocean_ice_nml project_ci (SIS2 default .true.; Roundabout default .false. ⇒ byte-identical, the house bit-identity rule). When on, evp_project_ci_impl projects ci (and hence pres_mice) forward each subcycle from the CALL’s initial concentration and cumulative elapsed time — ci_proj = ci*exp(-t_cum*sh_dd), t_cum = n*dt. D3 not ported: landfast (Lemieux/ITD, SIS2 default off), drag_max/MIN_OCN_INTERTIAL_H (default off), vel_underflow/str_underflow (default 0), weak_low_shear (default off), DT_RHEOLOGY (NSTEPS_DYN only), hi-freq/sigI/ sigII diagnostics, drag_bg_vel2 (SIS2 hardwires 0). CFL truncation’s CFL half is CLOSED by PR 36 — &ocean_ice_nml cfl_trunc (SIS2 CFL_TRUNCATE, default 0.5 there, 0.0 here ⇒ byte-identical) clips the FINAL transport velocity to 0.95*cfl_trunc*areaT(donor)/(dt_transport*dy_cu) (evp_truncate_final_impl), counting ice-bearing faces touched into a driver-logged warning (NOT an abort — SIS2 pairs its counter with MAXTRUNC=0, a run-stopper Roundabout does not port; PR-4b transport’s conservation/positivity check stays the backstop of last resort). cfl_trunc_dyn_its (SIS2 CFL_TRUNC_DYN_ITS, default off, matches) additionally clips to the EXACT bound at the bottom of every subcycle. Four documented divergences from SIS2’s port of this feature: (i) the bound uses the dt TRANSPORT will consume (dt_transport, an optional argument threaded through ice_evp_dynamics/ice_evp_step), NOT this call’s dt_slow — SIS2 assumes the two are the same dt, Roundabout decouples EVP (every outer step) from transport (thermo cadence); (ii) ci_proj (D2) is a local() scalar, not an array — SIS2 materialises it only for sigI/sigII/find_ice_strength diagnostics Roundabout does not have; (iii) the truncation count (n_trunc) drives a rank-0 driver WARNING, never an abort — no MAXTRUNC; (iv) SIS2 defaults CFL_TRUNCATE=0.5 / PROJECT_ICE_CONCENTRATION=.true., Roundabout defaults both off (house bit-identity rule) — the shipped polar_freezeup_dynamics.nml example carries SIS2’s defaults instead. Landfast/drag_max/underflows/weak_low_shear/ DT_RHEOLOGY/diagnostics remain unported. D4 v-momentum reads u_tmp (the PRE-update u), per SIS2 :1258-1273 — resolved toward SIS2 (the prototype’s Box2D used the updated u, a defect with no effect on any analytic gate). D5 ncat==1 lumped concentration: ci = 1 where m_ice > 0 (SIS2 has no lumped mode) — see ice_cell_concentration_impl in rdb_ice_state (shared with the tau coupler). D6 domain edges wall-or-periodic only; no OBC, no tripolar fold. Multi-rank (ice_evp_step with a decomposed bc): the ice uses the OCEAN’s decomposition and nghost. ui/vi are halo-exchanged (ocean_halo_face_x/_y, D1 seam ownership) at the top of EVERY subcycle and after the final CFL clip; the local periodic wraps run only on an axis the halo does not own (ocean_halo_is_decomposed_x/_y); mask_t pins a ghost band to land only on a PHYSICAL, non-periodic edge (bc%has_*), so an MPI seam ghost follows wet_T. The stresses are NOT exchanged: every stress kernel runs over the full local array and is point-local (or 4-cell around a corner) in exchanged inputs, so the ghost stresses are recomputed redundantly; only the outermost ring (one-sided corner strain) is wrong, and with nghost >= 3 no physical face reads it. The category inputs’ ghosts (mis/mice/ci) are the caller’s: the engine exchanges the category state at the end of every thermo block and at cold start (ocean_halo_exchange_ice_state). D7 one atmospheric stress field: the ice feels the FULL wind stress snapshot (tau_a_x/tau_a_y); no ice-specific bulk drag law (that part is unchanged — a future ice-specific bulk drag law is PR-55/RESUME #8’s problem). PR 62’s &ocean_ice_nml a_face_stress (default OFF) weights BOTH this wind stress AND the ice-ocean drag in the momentum balance by the face ice concentration a_u = 0.5*(ci(i-1,j)+ci(i,j)) (evp_u_momentum_impl/evp_v_momentum_impl), giving the textbook m du/dt = grad.sigma + a*(tau_a - tau_w) (Hibler 1979 eq. 1) and an EXACTLY closing ice<->ocean momentum budget at every fractional cover a, not just steady free drift. a_face_stress=.false. (default, byte-identical) is the legacy form: the ice absorbs the FULL wind and sheds the FULL drag while ice_ocean_stress_flux hands the ocean (1-a)*tau_a + a*fxoc — leaking (1-a)*(tau_a-fxoc) per face per step at fractional cover (zero at a in {0,1} and at steady free drift fxoc==tau_a, but nonzero in a generic transient — see the F5 caveat in ice_ocean_stress_flux’s docstring, rdb_ice_ocean_coupler). DELIBERATE DIVERGENCE FROM SIS2: SIS_C_dynamics weights NEITHER term (fxat/drag_u are bare against a per-TOTAL-area mis) and carries the same leak; the nearest SIS2 analogue, set_wind_stresses_C’s ice-cover- weighted interpolation, degenerates under Roundabout’s single wind field to a pure a>0 presence gate (kills the ghost-drift artefact below, does not close the budget). a_face_stress is therefore more correct than the reference on Hibler (1979)/CICE conservation grounds — the same call D4 already made toward SIS2. Weighting the wind ALONE (without the drag) is NOT a valid alternative: it converts today’s leak (zero at steady free drift) into a PERMANENT one, -(1-a)*a*tau_a, nonzero at steady state forever — do not “simplify” to a single term. a_face_stress=.true. also kills the ghost free-drift artefact at ice-free wet faces (today’s unweighted form converges a massless slab to the full Nansen free-drift speed sqrt(tau_a/(rho_o*Cdw)), regenerated every substep; with the knob on, a_fac == 0 branches uio_c to exactly 0 => ui==uo). Formerly PRE-EXISTING (not fixed by PR 62): the EVP’s a_u is built from ci_w (masked + periodic-wrapped); the coupler’s a_u (ice_ocean_stress_flux_impl) from the raw halo (ice_cell_concentration_impl), which nothing refreshed, so at the first physical face of a periodic domain the two could differ. CLOSED by the sea-ice MPI exchanges: the category state’s ghosts are exchanged — on one rank, wrapped — at the end of every thermo block and at cold start (ocean_halo_exchange_ice_state), so both gathers see the owner’s cells there (an answer change on periodic + dynamics configurations). D8 tau mediation is one-step-lagged (MEKE/frazil convention); a fresh run’s first outer step drives the ocean with pure wind — literally true as of PR 63 (previously the resume fold’s reconstruct gave a fresh run with ice at configure (1-a)*tau_a on step 1, an artefact of fxoc initialising to 0, not this claim; see PR 63’s F4 note below and its own plan §11.3/§14 Q5). The ELASTIC stress state (str_d/str_t/str_s) + u_ice/v_ice + fxoc/fyoc round-trip bit-exact through the restart; tau_x/tau_y (the field the ocean actually consumes) round-trips bit-exact too, as of PR 63 — ocean_sea_ice_t%tau_ocn_x/y mirror the exact blended value ice_ocean_stress_flux last wrote and ice_ocean_stress_resume_apply (rdb_ice_ocean_coupler) COPIES it back at configure, formula-agnostic. F4 (now historical): the PRE-PR-63 resume fold RECONSTRUCTED tau_x/tau_y from the checkpoint’s POST-thermo ci, which differs from the PRE-thermo/pre-transport ci the uninterrupted run’s blend actually used whenever a checkpoint step’s thermo/transport changed ci after the blend — closed by carrying the blend’s output instead of recomputing it. tau_a_x/tau_a_y remain a fresh configure-time snapshot, not restart-carried (D8 is otherwise unchanged: the mediation is still one-step lagged, only the RESUME path changed).


Uses

  • module~~rdb_ice_evp~~UsesGraph module~rdb_ice_evp rdb_ice_evp module~rdb_constants rdb_constants module~rdb_ice_evp->module~rdb_constants module~rdb_grid rdb_grid module~rdb_ice_evp->module~rdb_grid module~rdb_ice_column rdb_ice_column module~rdb_ice_evp->module~rdb_ice_column module~rdb_ice_state rdb_ice_state module~rdb_ice_evp->module~rdb_ice_state module~rdb_multilayer_state rdb_multilayer_state module~rdb_ice_evp->module~rdb_multilayer_state module~rdb_ocean_boundary_types rdb_ocean_boundary_types module~rdb_ice_evp->module~rdb_ocean_boundary_types module~rdb_ocean_halo rdb_ocean_halo module~rdb_ice_evp->module~rdb_ocean_halo module~rdb_ocean_metrics rdb_ocean_metrics module~rdb_ice_evp->module~rdb_ocean_metrics module~rdb_ocean_periodic rdb_ocean_periodic module~rdb_ice_evp->module~rdb_ocean_periodic pic_types pic_types module~rdb_constants->pic_types module~rdb_grid->module~rdb_constants module~rdb_ice_column->module~rdb_constants module~rdb_ice_enthalpy rdb_ice_enthalpy module~rdb_ice_column->module~rdb_ice_enthalpy module~rdb_ice_mass rdb_ice_mass module~rdb_ice_column->module~rdb_ice_mass module~rdb_ice_optics rdb_ice_optics module~rdb_ice_column->module~rdb_ice_optics module~rdb_ice_state->module~rdb_constants module~rdb_ice_state->module~rdb_grid module~rdb_ice_state->module~rdb_ice_column iso_fortran_env iso_fortran_env module~rdb_ice_state->iso_fortran_env module~rdb_ice_state->module~rdb_ice_enthalpy module~rdb_mem_report rdb_mem_report module~rdb_ice_state->module~rdb_mem_report module~rdb_multilayer_state->module~rdb_constants module~rdb_multilayer_state->module~rdb_grid module~rdb_multilayer_state->iso_fortran_env module~rdb_efp rdb_efp module~rdb_multilayer_state->module~rdb_efp module~rdb_error_ring rdb_error_ring module~rdb_multilayer_state->module~rdb_error_ring module~rdb_multilayer_state->module~rdb_mem_report module~rdb_tracer rdb_tracer module~rdb_multilayer_state->module~rdb_tracer pic_logger pic_logger module~rdb_multilayer_state->pic_logger module~rdb_ocean_boundary_types->module~rdb_constants module~rdb_ocean_boundary_types->module~rdb_grid module~rdb_ocean_boundary_types->iso_fortran_env module~rdb_ocean_boundary_types->module~rdb_mem_report module~rdb_ocean_status rdb_ocean_status module~rdb_ocean_boundary_types->module~rdb_ocean_status module~rdb_ocean_tide_astro rdb_ocean_tide_astro module~rdb_ocean_boundary_types->module~rdb_ocean_tide_astro pic_ascii pic_ascii module~rdb_ocean_boundary_types->pic_ascii module~rdb_ocean_boundary_types->pic_logger module~rdb_ocean_halo->module~rdb_constants module~rdb_ocean_halo->module~rdb_ocean_periodic module~rdb_comm_env rdb_comm_env module~rdb_ocean_halo->module~rdb_comm_env module~rdb_decomp rdb_decomp module~rdb_ocean_halo->module~rdb_decomp module~rdb_ocean_halo->module~rdb_error_ring module~rdb_ocean_halo_counters rdb_ocean_halo_counters module~rdb_ocean_halo->module~rdb_ocean_halo_counters module~rdb_ocean_halo->module~rdb_ocean_status module~rdb_ocean_halo->pic_logger pic_mpi_lib pic_mpi_lib module~rdb_ocean_halo->pic_mpi_lib pic_strings pic_strings module~rdb_ocean_halo->pic_strings module~rdb_ocean_metrics->module~rdb_constants module~rdb_ocean_metrics->module~rdb_grid module~rdb_ocean_metrics->iso_fortran_env module~rdb_ocean_metrics->module~rdb_error_ring module~rdb_io_netcdf rdb_io_netcdf module~rdb_ocean_metrics->module~rdb_io_netcdf module~rdb_ocean_metrics->module~rdb_mem_report module~rdb_ocean_bipolar rdb_ocean_bipolar module~rdb_ocean_metrics->module~rdb_ocean_bipolar module~rdb_ocean_fold rdb_ocean_fold module~rdb_ocean_metrics->module~rdb_ocean_fold module~rdb_ocean_metrics->module~rdb_ocean_status netcdf netcdf module~rdb_ocean_metrics->netcdf module~rdb_ocean_metrics->pic_logger module~rdb_ocean_metrics->pic_strings module~rdb_ocean_periodic->module~rdb_constants module~rdb_ocean_periodic->module~rdb_grid module~rdb_ocean_periodic->module~rdb_multilayer_state module~rdb_ocean_periodic->module~rdb_ocean_boundary_types module~rdb_comm_env->module~rdb_constants module~rdb_comm_env->iso_fortran_env module~rdb_comm_env->pic_mpi_lib module~rdb_config rdb_config module~rdb_decomp->module~rdb_config module~rdb_efp->iso_fortran_env ieee_arithmetic ieee_arithmetic module~rdb_efp->ieee_arithmetic module~rdb_error_ring->pic_logger module~rdb_ice_enthalpy->module~rdb_constants module~rdb_ice_mass->module~rdb_constants module~rdb_ice_mass->module~rdb_ice_enthalpy module~rdb_ice_optics->module~rdb_constants module~rdb_ice_optics->module~rdb_ice_enthalpy module~rdb_io_netcdf->module~rdb_constants module~rdb_io_netcdf->iso_fortran_env module~rdb_io_netcdf->module~rdb_error_ring module~rdb_io_netcdf->netcdf module~rdb_io_netcdf->pic_logger module~rdb_io_netcdf->pic_strings iso_c_binding iso_c_binding module~rdb_io_netcdf->iso_c_binding module~rdb_mem_report->module~rdb_constants module~rdb_mem_report->iso_fortran_env module~rdb_mem_report->pic_logger module~rdb_mem_report->pic_strings module~rdb_ocean_bipolar->module~rdb_constants module~rdb_ocean_fold->module~rdb_constants module~rdb_ocean_halo_counters->iso_fortran_env module~rdb_ocean_halo_counters->pic_strings module~rdb_ocean_tide_astro->module~rdb_constants module~rdb_ocean_tide_astro->iso_fortran_env module~rdb_tracer->module~rdb_constants module~rdb_tracer->module~rdb_grid module~rdb_tracer->iso_fortran_env module~rdb_tracer->module~rdb_mem_report module~rdb_config->module~rdb_constants module~rdb_config->module~rdb_error_ring module~rdb_config->module~rdb_ice_enthalpy module~rdb_config->module~rdb_ocean_status module~rdb_config->pic_ascii module~rdb_config->pic_logger module~rdb_config->pic_strings module~rdb_ice_init rdb_ice_init module~rdb_config->module~rdb_ice_init module~rdb_nml_schema rdb_nml_schema module~rdb_config->module~rdb_nml_schema module~rdb_ice_init->module~rdb_constants module~rdb_ice_init->module~rdb_grid module~rdb_ice_init->module~rdb_ice_column module~rdb_ice_init->module~rdb_ice_state module~rdb_ice_init->module~rdb_multilayer_state module~rdb_ice_init->module~rdb_ocean_metrics module~rdb_ice_init->module~rdb_ice_enthalpy module~rdb_nml_schema->module~rdb_constants module~rdb_nml_schema->module~rdb_error_ring module~rdb_nml_schema->pic_logger

Used by

  • module~~rdb_ice_evp~~UsedByGraph module~rdb_ice_evp rdb_ice_evp module~rdb_ocean_engine rdb_ocean_engine module~rdb_ocean_engine->module~rdb_ice_evp module~rdb_driver rdb_driver module~rdb_driver->module~rdb_ocean_engine module~rdb_handle rdb_handle module~rdb_handle->module~rdb_ocean_engine module~rdb_ocean_api rdb_ocean_api module~rdb_ocean_api->module~rdb_ocean_engine module~rdb_ocean_api->module~rdb_handle

Variables

Type Visibility Attributes Name Initial
real(kind=wp), private, parameter :: EVP_DRAG_LINEARIZE_THRESHOLD = 1.0e8_wp

Large-b_vel0 cutover in the semi-implicit ice-ocean drag solve (SIS2 SIS_C_dynamics): when b_vel0**2 > EVP_DRAG_LINEARIZE_THRESHOLD * I_cdRhoDt * |m_*io_explicit| the quadratic drag term is negligible against the linear one, so the predicted relative velocity linearizes to m_*io_explicit * I_cdRhoDt / b_vel0. Same value in the u- and v-momentum kernels.

real(kind=wp), private, parameter :: M_NEGLECT_FACTOR = 1.0e-30_wp

SIS2 H_subroundoff (:757) — m_neglect = ICE_RHO_ICE*1e-30.

real(kind=wp), private, parameter :: TRUNC_BACKOFF = 0.95_wp

PR 36: back-off factor on the FINAL CFL velocity clip (SIS2 SIS_dyn_cgrid.F90:1456) — the clipped value is set to 0.95*bound, not the bound itself, so a re-check cannot be marginal. The in-loop (cfl_trunc_dyn_its) clip uses the exact bound (backoff = 1.0) instead — it is not the last word on the velocity before transport reads it.


Derived Types

type, public ::  ice_evp_params_t

EVP physical + numerical parameters (SIS2 SIS_C_dyn_CS subset). Passed intent(in) into the core; every scalar is hoisted to a local before the substep loop (no derived-type deref inside a do concurrent). Defaults + meanings mirror the &ocean_ice_nml block in rdb_config.F90 (the config is the authoritative knob set; keep the two in sync).

Components

Type Visibility Attributes Name Initial
logical, public :: a_face_stress = .false.

PR 62: weight the atmospheric stress AND the ice-ocean drag in the momentum balance by the face ice concentration a_u (&ocean_ice_nml a_face_stress). Default off ⇒ byte-identical.

real(kind=wp), public :: c0 = 20.0_wp

SIS2 ICE_STRENGTH_CSTAR — ice-strength exponent constant [nondim].

real(kind=wp), public :: cdw = 3.24e-3_wp

SIS2 ICE_CDRAG_WATER — ice-ocean drag coefficient [nondim].

real(kind=wp), public :: cfl_trunc = 0.0_wp

PR 36: SIS2 CFL_TRUNCATE (SIS2 default 0.5). Transport-CFL ceiling on the final ice velocity; 0 disables the clip. Type-level default is the bit-identity mechanism for every test call site that does not set it.

logical, public :: cfl_trunc_dyn_its = .false.

PR 36: SIS2 CFL_TRUNC_DYN_ITS (SIS2 default .false., matches). Also clip at the bottom of every EVP subcycle.

real(kind=wp), public :: del_sh_min_scale = 2.0_wp

SIS2 ICE_DEL_SH_MIN_SCALE — viscosity-floor scale [nondim].

real(kind=wp), public :: ec = 2.0_wp

SIS2 ICE_YIELD_ELLIPTICITY — yield-curve axis ratio [nondim]. 0 => cavitating-fluid rheology (str_t/str_s stay exactly 0).

integer, public :: evp_sub_steps = 432

SIS2 NSTEPS_DYN — EVP subcycles per slow (outer) step.

real(kind=wp), public :: p0 = 2.75e4_wp

SIS2 ICE_STRENGTH_PSTAR — ice-strength pressure constant [Pa].

logical, public :: project_ci = .false.

PR 36: SIS2 PROJECT_ICE_CONCENTRATION (SIS2 default .true.). Project ci forward along the current divergence each subcycle and recompute pres_mice from it.

real(kind=wp), public :: rho_ocean = 1030.0_wp

SIS2 RHO_OCEAN — ice-drag reference density [kg/m^3]. Deliberately independent of the ocean’s rho0 (usually 1035).

real(kind=wp), public :: tdamp = -0.2_wp

SIS2 ICE_TDAMP_ELASTIC — elastic damping timescale selector. > 0 => seconds; == 0 => max(0.2*dt_slow, 3*dt); < 0 => the special case max(|tdamp|*dt_slow, 3*dt) (i.e. |tdamp| is a fraction of the slow step). Sign-free — no positivity guard.


Functions

public pure function ice_evp_mi_ratio_point(mis_sw, mis_se, mis_nw, mis_ne, mask_u_below, mask_u_above, mask_v_left, mask_v_right, mask_q, area_sw, area_se, area_nw, area_ne, mask_t_sw, mask_t_se, mask_t_nw, mask_t_ne, m_neglect2, m_neglect4) result(mi_ratio)

mi_ratio_A_q at a single corner (SIS2 :926-964), FULL form — all four branches (interior / corner-coast / straight-coast / land). weak_coast_stress=.false. hardwired (SIS2 default): sum_area is the MASKED area sum of the 4 surrounding T-cells. Factored out of the fill kernel so a unit test can pin it directly (SPEC §7 gate 8).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: mis_sw

Ice+snow mass per cell area at the 4 T-cells around the corner (SW/SE/NW/NE, rdb wet_q convention).

real(kind=wp), intent(in) :: mis_se

Ice+snow mass per cell area at the 4 T-cells around the corner (SW/SE/NW/NE, rdb wet_q convention).

real(kind=wp), intent(in) :: mis_nw

Ice+snow mass per cell area at the 4 T-cells around the corner (SW/SE/NW/NE, rdb wet_q convention).

real(kind=wp), intent(in) :: mis_ne

Ice+snow mass per cell area at the 4 T-cells around the corner (SW/SE/NW/NE, rdb wet_q convention).

real(kind=wp), intent(in) :: mask_u_below

u-face masks below/above the corner (SIS2 mask2dCu(I,j)/mask2dCu(I,j+1)).

real(kind=wp), intent(in) :: mask_u_above

u-face masks below/above the corner (SIS2 mask2dCu(I,j)/mask2dCu(I,j+1)).

real(kind=wp), intent(in) :: mask_v_left

v-face masks left/right of the corner (SIS2 mask2dCv(i,J)/mask2dCv(i+1,J)).

real(kind=wp), intent(in) :: mask_v_right

v-face masks left/right of the corner (SIS2 mask2dCv(i,J)/mask2dCv(i+1,J)).

real(kind=wp), intent(in) :: mask_q

mask2dBu at this corner (1 = genuinely interior ocean point).

real(kind=wp), intent(in) :: area_sw

T-cell areas at the 4 surrounding cells.

real(kind=wp), intent(in) :: area_se

T-cell areas at the 4 surrounding cells.

real(kind=wp), intent(in) :: area_nw

T-cell areas at the 4 surrounding cells.

real(kind=wp), intent(in) :: area_ne

T-cell areas at the 4 surrounding cells.

real(kind=wp), intent(in) :: mask_t_sw

T-cell wet masks at the 4 surrounding cells (land => 0).

real(kind=wp), intent(in) :: mask_t_se

T-cell wet masks at the 4 surrounding cells (land => 0).

real(kind=wp), intent(in) :: mask_t_nw

T-cell wet masks at the 4 surrounding cells (land => 0).

real(kind=wp), intent(in) :: mask_t_ne

T-cell wet masks at the 4 surrounding cells (land => 0).

real(kind=wp), intent(in) :: m_neglect2
real(kind=wp), intent(in) :: m_neglect4

Return Value real(kind=wp)

public pure function ice_evp_params_from_config(p0, c0, ec, cdw, rho_ocean, del_sh_min_scale, tdamp, evp_sub_steps, a_face_stress, cfl_trunc, cfl_trunc_dyn_its, project_ci) result(par)

Small constructor — build once from &ocean_ice_nml config. 11 args (> the style guide’s 6): pre-existing deviation, sanctioned by the derived-type-grouping escape hatch (FORTRAN_STYLE.md §Public procedure arguments) — the whole point of ice_evp_params_t is to be this constructor’s one-shot host. Positional (not optional): a knob threaded to the config/schema but not to this constructor reads from the namelist, validates, and does nothing — the dead-knob class the audit indicts. The compiler catches the omission; optional would not.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: p0
real(kind=wp), intent(in) :: c0
real(kind=wp), intent(in) :: ec
real(kind=wp), intent(in) :: cdw
real(kind=wp), intent(in) :: rho_ocean
real(kind=wp), intent(in) :: del_sh_min_scale
real(kind=wp), intent(in) :: tdamp
integer, intent(in) :: evp_sub_steps
logical, intent(in) :: a_face_stress
real(kind=wp), intent(in) :: cfl_trunc
logical, intent(in) :: cfl_trunc_dyn_its
logical, intent(in) :: project_ci

Return Value type(ice_evp_params_t)


Subroutines

public 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).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: wet_t(nx,ny)
real(kind=wp), intent(out) :: mask_t(nx,ny)
real(kind=wp), intent(out) :: mask_u(nx+1,ny)
real(kind=wp), intent(out) :: mask_v(nx,ny+1)
real(kind=wp), intent(out) :: mask_q(nx+1,ny+1)
integer, intent(in) :: nx_phys
integer, intent(in) :: ny_phys
integer, intent(in) :: nghost
logical, intent(in) :: wrap_x
logical, intent(in) :: wrap_y
logical, intent(in) :: land_w
logical, intent(in) :: land_e
logical, intent(in) :: land_s
logical, intent(in) :: land_n
integer, intent(in) :: nx
integer, intent(in) :: ny

public pure subroutine evp_q_and_mi_ratio_impl(areaT, f_corner, mask_t, mask_u, mask_v, mask_q, mis, m_neglect, m_neglect2, m_neglect4, q, mi_ratio_a_q, nx, ny)

q(ic,jc) = f_corner*tot_area / (Σ areaT*mis over the 4 cells + tot_area*m_neglect) (:977-982); mi_ratio_A_q via ice_evp_mi_ratio_point (requirement 5). 4 T-cells around corner (ic,jc): (ic-1,jc-1) (ic,jc-1) (ic-1,jc) (ic,jc) (SW/SE/NW/NE, §1). Array-edge corners (no T-cell on one side) get q=0/mi_ratio=0 (land-corner convention — consistent with mask_t=0 beyond the array edge).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: areaT(nx,ny)
real(kind=wp), intent(in) :: f_corner(nx+1,ny+1)
real(kind=wp), intent(in) :: mask_t(nx,ny)
real(kind=wp), intent(in) :: mask_u(nx+1,ny)
real(kind=wp), intent(in) :: mask_v(nx,ny+1)
real(kind=wp), intent(in) :: mask_q(nx+1,ny+1)
real(kind=wp), intent(in) :: mis(nx,ny)
real(kind=wp), intent(in) :: m_neglect
real(kind=wp), intent(in) :: m_neglect2
real(kind=wp), intent(in) :: m_neglect4
real(kind=wp), intent(out) :: q(nx+1,ny+1)
real(kind=wp), intent(out) :: mi_ratio_a_q(nx+1,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny

public pure subroutine evp_truncate_final_impl(areaT, dy_cu, dx_cv, mi_u, mi_v, ui, vi, cfl_trunc, dt_tr, m_neglect, nghost, nx_phys, ny_phys, nx, ny, n_trunc, count_w, count_s)

PR 36: the FINAL CFL clip (SIS2 :1443-1500) – TRUNC_BACKOFF (0.95) back-off instead of the exact bound, PLUS a count of the ice-bearing faces it touched (mi > m_neglect, SIS2 :1466,1469 – massless faces clip silently, matching SIS2: counting them would flood the driver’s warning with meaningless ice-free clips at every margin). Not a do concurrent: reductions use !$acc parallel loop reduction(...) (ice_compress_impl is the local precedent for a reduction that also mutates the arrays it walks). Same bound algebra as evp_truncate_velocity_impl, duplicated rather than shared: the in-loop variant runs evp_sub_steps (432 by default) times per outer step and must NOT carry a reduction (each would be a device->host sync); this variant runs once and must. CLAUDE.md’s “duplicate explicitly” rule – merging the two costs 432 syncs per outer step.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: areaT(nx,ny)
real(kind=wp), intent(in) :: dy_cu(nx+1,ny)
real(kind=wp), intent(in) :: dx_cv(nx,ny+1)
real(kind=wp), intent(in) :: mi_u(nx+1,ny)
real(kind=wp), intent(in) :: mi_v(nx,ny+1)
real(kind=wp), intent(inout) :: ui(nx+1,ny)
real(kind=wp), intent(inout) :: vi(nx,ny+1)
real(kind=wp), intent(in) :: cfl_trunc
real(kind=wp), intent(in) :: dt_tr
real(kind=wp), intent(in) :: m_neglect
integer, intent(in) :: nghost
integer, intent(in) :: nx_phys
integer, intent(in) :: ny_phys
integer, intent(in) :: nx
integer, intent(in) :: ny
integer, intent(out) :: n_trunc
logical, intent(in), optional :: count_w

Count the WEST (i = nghost+1) / SOUTH (j = nghost+1) edge face. That face is clipped either way; it is COUNTED only when this tile owns it — a physical, non-periodic edge. Across an MPI seam the west/south neighbour owns it (D1), and across a periodic seam it is the same face as the east/north edge face, so counting it there would count one face twice in the rank-summed total. Absent => .true. (count every face).

logical, intent(in), optional :: count_s

Count the WEST (i = nghost+1) / SOUTH (j = nghost+1) edge face. That face is clipped either way; it is COUNTED only when this tile owns it — a physical, non-periodic edge. Across an MPI seam the west/south neighbour owns it (D1), and across a periodic seam it is the same face as the east/north edge face, so counting it there would count one face twice in the rank-summed total. Absent => .true. (count every face).

public subroutine ice_evp_dynamics(grid, metrics, f_corner, mis, mice, ci, uo, vo, tau_ax, tau_ay, ui, vi, str_d, str_t, str_s, fxoc, fyoc, dt_slow, par, periodic_x, periodic_y, ws, dt_transport, n_trunc, halo_x, halo_y, land_w, land_e, land_s, land_n)

One outer (slow) EVP call: evp_sub_steps subcycles advancing ui/vi/str_d/str_t/str_s, plus the subcycle-averaged ice->ocean stress fxoc/fyoc.

Read more…

Arguments

Type IntentOptional Attributes Name
type(hgrid_t), intent(in) :: grid
type(ocean_metrics_t), intent(in) :: metrics
real(kind=wp), intent(in) :: f_corner(:,:)
real(kind=wp), intent(in) :: mis(:,:)
real(kind=wp), intent(in) :: mice(:,:)
real(kind=wp), intent(in) :: ci(:,:)
real(kind=wp), intent(in) :: uo(:,:)
real(kind=wp), intent(in) :: vo(:,:)
real(kind=wp), intent(in) :: tau_ax(:,:)
real(kind=wp), intent(in) :: tau_ay(:,:)
real(kind=wp), intent(inout) :: ui(:,:)
real(kind=wp), intent(inout) :: vi(:,:)
real(kind=wp), intent(inout) :: str_d(:,:)
real(kind=wp), intent(inout) :: str_t(:,:)
real(kind=wp), intent(inout) :: str_s(:,:)
real(kind=wp), intent(inout) :: fxoc(:,:)
real(kind=wp), intent(inout) :: fyoc(:,:)
real(kind=wp), intent(in) :: dt_slow
type(ice_evp_params_t), intent(in) :: par
logical, intent(in) :: periodic_x
logical, intent(in) :: periodic_y
type(evp_workspace_t), intent(inout) :: ws
real(kind=wp), intent(in), optional :: dt_transport

PR 36: dt TRANSPORT will use for the CFL bound. Absent => dt_slow.

integer, intent(out), optional :: n_trunc

PR 36: count of ice-bearing faces the final clip touched.

logical, intent(in), optional :: halo_x

The x / y axis is split across ranks: exchange ui/vi through the ocean halo and skip the local periodic wrap on that axis. Absent => .false. (single rank, the unit-test seam).

logical, intent(in), optional :: halo_y

The x / y axis is split across ranks: exchange ui/vi through the ocean halo and skip the local periodic wrap on that axis. Absent => .false. (single rank, the unit-test seam).

logical, intent(in), optional :: land_w

Pin that ghost band of mask_t to land (a physical, non-periodic edge). Absent => .not. periodic_x (W/E) / .not. periodic_y (S/N): the single-rank wall-or-wrap policy.

logical, intent(in), optional :: land_e

Pin that ghost band of mask_t to land (a physical, non-periodic edge). Absent => .not. periodic_x (W/E) / .not. periodic_y (S/N): the single-rank wall-or-wrap policy.

logical, intent(in), optional :: land_s

Pin that ghost band of mask_t to land (a physical, non-periodic edge). Absent => .not. periodic_x (W/E) / .not. periodic_y (S/N): the single-rank wall-or-wrap policy.

logical, intent(in), optional :: land_n

Pin that ghost band of mask_t to land (a physical, non-periodic edge). Absent => .not. periodic_x (W/E) / .not. periodic_y (S/N): the single-rank wall-or-wrap policy.

public subroutine ice_evp_step(grid, metrics, f_corner, ice, ms, dt_slow, par, bc, dt_transport, n_trunc)

Gathers mis/mice/ci from the category state (mode-branched, mirrors PR 4b’s IST->CAS dispatch), pulls the one-step-lagged ocean surface velocity, and calls ice_evp_dynamics on ice%u_ice/v_ice/str_d/str_t/str_s/fxoc/fyoc. No-op when the ice slot is not live or dynamics is off (defence-in-depth; the driver already gates this call on ice%dynamics).

Read more…

Arguments

Type IntentOptional Attributes Name
type(hgrid_t), intent(in) :: grid
type(ocean_metrics_t), intent(in) :: metrics
real(kind=wp), intent(in) :: f_corner(:,:)

Coriolis parameter at corners, shape (nx+1,ny+1) (coriolis_adv_t%f_corner).

type(ocean_sea_ice_t), intent(inout) :: ice
type(multilayer_state_t), intent(in) :: ms
real(kind=wp), intent(in) :: dt_slow
type(ice_evp_params_t), intent(in) :: par
type(ocean_bc_state_t), intent(in) :: bc

Edge policy: periodic_x/_y (wrap) and has_west/_east/_south/ _north (physical edge vs MPI seam). A ghost band is pinned to land only on a physical, non-periodic edge; on a decomposed axis the halo fills it (see D6 in the module docstring).

real(kind=wp), intent(in), optional :: dt_transport

PR 36: the dt the TRANSPORT step will actually consume (ocean_dyn%therm_dt(dt)), NOT this call’s dt_slow — EVP runs every outer step, transport at thermo cadence. Absent => dt_slow (read only when par%cfl_trunc > 0).

integer, intent(out), optional :: n_trunc

PR 36: count of ice-bearing faces the final CFL clip touched (0 when par%cfl_trunc <= 0). Mirrors ice_transport_step’s ok idiom — the caller (driver) logs, this module does not.

private pure subroutine evp_average_stress_impl(mask_u, mask_v, fxoc, fyoc, evp_sub_steps, nx, ny)

fxoc *= mask_u/evp_sub_steps, fyoc *= mask_v/evp_sub_steps (:1415-1441).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: mask_u(nx+1,ny)
real(kind=wp), intent(in) :: mask_v(nx,ny+1)
real(kind=wp), intent(inout) :: fxoc(nx+1,ny)
real(kind=wp), intent(inout) :: fyoc(nx,ny+1)
integer, intent(in) :: evp_sub_steps
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_copy_u_impl(ui, u_tmp, nx, ny)

u_tmp = ui (full array — the v-momentum MUST read pre-update u, D4). Explicit do concurrent element copy rather than a bare whole-array assignment (repo convention for device-resident arrays — see rdb_ml_dynamics’s h_layer0 save pattern).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: ui(nx+1,ny)
real(kind=wp), intent(out) :: u_tmp(nx+1,ny)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_fill_cell_fields_impl(mask_t, mis_in, mice_in, ci_in, mis_out, mice_out, ci_out, nx, ny)

Interior copy of the gathered mis/mice/ci, masked to mask_t (defence-in-depth beyond the caller’s own wet_T gate); ghost rows/cols zeroed (the periodic wrap that follows fills them).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: mask_t(nx,ny)
real(kind=wp), intent(in) :: mis_in(nx,ny)
real(kind=wp), intent(in) :: mice_in(nx,ny)
real(kind=wp), intent(in) :: ci_in(nx,ny)
real(kind=wp), intent(out) :: mis_out(nx,ny)
real(kind=wp), intent(out) :: mice_out(nx,ny)
real(kind=wp), intent(out) :: ci_out(nx,ny)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_mi_face_impl(mis, mi_u, mi_v, nx, ny)

mi_u(i,j) = 0.5*(mis(i-1,j)+mis(i,j)); mi_v(i,j) = 0.5*(mis(i,j-1)+mis(i,j)) (SIS2 :967-974, rdb index translation §1). Array-edge faces (i=1/i=nx+1, j=1/ j=ny+1) have no neighbour on one side; mis at those ghost rows/cols was already periodic-wrapped or zeroed, so a naive mis(i-1,j)/mis(i,j) read is always in-bounds here EXCEPT at the two hard array edges themselves — those faces are handled explicitly.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: mis(nx,ny)
real(kind=wp), intent(out) :: mi_u(nx+1,ny)
real(kind=wp), intent(out) :: mi_v(nx,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_pres_mice_impl(dxT, dyT, ci, p0_rho, c0, del_sh_min_scale, tdamp_eff, dt, pres_mice, del_sh_min_pr, nx, ny)

pres_mice = p0_rho*exp(-c0*max(1-ci,0)) (:878); dxharm = 2*dxT*dyT/(dxT+dyT); del_sh_min_pr = 2*del_sh_min_scale*dt^2 / (Tdamp*dxharm^2) guarded on dxharm > 0 (:880-890).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: dxT(nx,ny)
real(kind=wp), intent(in) :: dyT(nx,ny)
real(kind=wp), intent(in) :: ci(nx,ny)
real(kind=wp), intent(in) :: p0_rho
real(kind=wp), intent(in) :: c0
real(kind=wp), intent(in) :: del_sh_min_scale
real(kind=wp), intent(in) :: tdamp_eff
real(kind=wp), intent(in) :: dt
real(kind=wp), intent(out) :: pres_mice(nx,ny)
real(kind=wp), intent(out) :: del_sh_min_pr(nx,ny)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_project_ci_impl(ci, sh_dd, dt_cum, p0_rho, c0, pres_mice, nx, ny)

PR 36: PROJECT_ICE_CONCENTRATION (SIS2 SIS_dyn_cgrid.F90:1064- 1077). ci_proj = ci*exp(-dt_cum*sh_dd) then pres_mice = p0_rho*exp(-c0*max(1-ci_proj, 0)). ci_proj is a local() scalar, NOT an array: SIS2 materialises it only for the sigI/sigII/find_ice_strength diagnostics Roundabout does not have (documented divergence). del_sh_min_pr is NOT recomputed here (it has no ci dependence, evp_pres_mice_impl above). ci_proj is deliberately unclamped above 1 – max(1-ci_proj, 0) already saturates the effect at p0_rho, and for dt_cum*|sh_dd| > 709 (an unreachable regime in any sane run) exp overflows to +Inf, max(1-Inf, 0) = 0, exp(0) = 1 – IEEE launders the overflow to exactly the correct saturated value, so no guard is needed (SIS2 has none either).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: ci(nx,ny)
real(kind=wp), intent(in) :: sh_dd(nx,ny)
real(kind=wp), intent(in) :: dt_cum
real(kind=wp), intent(in) :: p0_rho
real(kind=wp), intent(in) :: c0
real(kind=wp), intent(inout) :: pres_mice(nx,ny)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_sh_dd_dt_impl(dy_dxT, dx_dyT, iareaT, idyCu, idxCv, dyCu, dxCv, ui, vi, sh_dd, sh_dt, nx, ny)

sh_Dt / sh_Dd at cells (:1053-1061).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: dy_dxT(nx,ny)
real(kind=wp), intent(in) :: dx_dyT(nx,ny)
real(kind=wp), intent(in) :: iareaT(nx,ny)
real(kind=wp), intent(in) :: idyCu(nx+1,ny)
real(kind=wp), intent(in) :: idxCv(nx,ny+1)
real(kind=wp), intent(in) :: dyCu(nx+1,ny)
real(kind=wp), intent(in) :: dxCv(nx,ny+1)
real(kind=wp), intent(in) :: ui(nx+1,ny)
real(kind=wp), intent(in) :: vi(nx,ny+1)
real(kind=wp), intent(out) :: sh_dd(nx,ny)
real(kind=wp), intent(out) :: sh_dt(nx,ny)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_sh_ds_impl(dx_dyBu, dy_dxBu, idxCu, idyCv, mask_q, ui, vi, sh_ds, nx, ny)

sh_Ds at corners (:1045-1050) — requirement (4): the SINGLE scalar no-slip factor (2-mask_q) on the WHOLE combined strain. Computed over the interior+1 ring (ic,jc in [1,nx+1]x[1,ny+1] — the full corner array; out-of-band neighbours contribute 0 via zero ghost velocities at the hard array edges, never per-term mirroring).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: dx_dyBu(nx+1,ny+1)
real(kind=wp), intent(in) :: dy_dxBu(nx+1,ny+1)
real(kind=wp), intent(in) :: idxCu(nx+1,ny)
real(kind=wp), intent(in) :: idyCv(nx,ny+1)
real(kind=wp), intent(in) :: mask_q(nx+1,ny+1)
real(kind=wp), intent(in) :: ui(nx+1,ny)
real(kind=wp), intent(in) :: vi(nx,ny+1)
real(kind=wp), intent(out) :: sh_ds(nx+1,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_str_s_relax_impl(areaT, zeta, sh_ds, mi_ratio_a_q, i_1pdt_t, dt_2tdamp, i_ec2, str_s, nx, ny)

str_s relax (:1137-1143). Corners in [1,nx+1]x[1,ny+1]; the 4 surrounding T-cells at an array-edge corner are handled by zeta’s own ghost values (zero-mass ghost cells => zeta=0 there, contributing nothing) — no special-case branch needed.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: areaT(nx,ny)
real(kind=wp), intent(in) :: zeta(nx,ny)
real(kind=wp), intent(in) :: sh_ds(nx+1,ny+1)
real(kind=wp), intent(in) :: mi_ratio_a_q(nx+1,ny+1)
real(kind=wp), intent(in) :: i_1pdt_t
real(kind=wp), intent(in) :: dt_2tdamp
real(kind=wp), intent(in) :: i_ec2
real(kind=wp), intent(inout) :: str_s(nx+1,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_stress_relax_impl(zeta, sh_dd, sh_dt, pres_mice, mice, i_1pdt_t, dt_2tdamp, i_ec2, str_d, str_t, nx, ny)

str_d/str_t semi-implicit relax (:1124-1134), non-weak_low_shear branch only.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: zeta(nx,ny)
real(kind=wp), intent(in) :: sh_dd(nx,ny)
real(kind=wp), intent(in) :: sh_dt(nx,ny)
real(kind=wp), intent(in) :: pres_mice(nx,ny)
real(kind=wp), intent(in) :: mice(nx,ny)
real(kind=wp), intent(in) :: i_1pdt_t
real(kind=wp), intent(in) :: dt_2tdamp
real(kind=wp), intent(in) :: i_ec2
real(kind=wp), intent(inout) :: str_d(nx,ny)
real(kind=wp), intent(inout) :: str_t(nx,ny)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_truncate_velocity_impl(areaT, dy_cu, dx_cv, ui, vi, cfl_trunc, dt_tr, backoff, nghost, nx_phys, ny_phys, nx, ny)

PR 36: the shared CFL-clip algebra – the transport-CFL bound on the ice velocity (SIS2 SIS_dyn_cgrid.F90:839-870, the in-loop half at :1338-1361, the final half at :1443-1500; this routine is the counting-free, caller-chosen-backoff form both reuse; evp_truncate_final_impl below wraps it with the 0.95 back-off and the mi > m_neglect count).

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: areaT(nx,ny)
real(kind=wp), intent(in) :: dy_cu(nx+1,ny)
real(kind=wp), intent(in) :: dx_cv(nx,ny+1)
real(kind=wp), intent(inout) :: ui(nx+1,ny)
real(kind=wp), intent(inout) :: vi(nx,ny+1)
real(kind=wp), intent(in) :: cfl_trunc
real(kind=wp), intent(in) :: dt_tr
real(kind=wp), intent(in) :: backoff
integer, intent(in) :: nghost
integer, intent(in) :: nx_phys
integer, intent(in) :: ny_phys
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_u_momentum_impl(idxCu, idyCu, dy2h, dx2q, iareaCu, mask_u, mi_u, mi_v, q, str_d, str_t, str_s, uo, vo, tau_ax, ui, vi, fxoc, m_neglect, i_cdrhodt, cdrho, dt, nx_phys, ny_phys, nghost, nx, ny, a_u, a_face_on)

u-momentum (:1172-1231, requirement 1: fxic_now carries the FULL str_t force term). Loop over u-faces ng+1..ng+nxp+1 x ng+1..ng+nyp — each iteration writes only its own face.

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: idxCu(nx+1,ny)
real(kind=wp), intent(in) :: idyCu(nx+1,ny)
real(kind=wp), intent(in) :: dy2h(nx,ny)
real(kind=wp), intent(in) :: dx2q(nx+1,ny+1)
real(kind=wp), intent(in) :: iareaCu(nx+1,ny)
real(kind=wp), intent(in) :: mask_u(nx+1,ny)
real(kind=wp), intent(in) :: mi_u(nx+1,ny)
real(kind=wp), intent(in) :: mi_v(nx,ny+1)
real(kind=wp), intent(in) :: q(nx+1,ny+1)
real(kind=wp), intent(in) :: str_d(nx,ny)
real(kind=wp), intent(in) :: str_t(nx,ny)
real(kind=wp), intent(in) :: str_s(nx+1,ny+1)
real(kind=wp), intent(in) :: uo(nx+1,ny)
real(kind=wp), intent(in) :: vo(nx,ny+1)
real(kind=wp), intent(in) :: tau_ax(nx+1,ny)
real(kind=wp), intent(inout) :: ui(nx+1,ny)
real(kind=wp), intent(in) :: vi(nx,ny+1)
real(kind=wp), intent(inout) :: fxoc(nx+1,ny)
real(kind=wp), intent(in) :: m_neglect
real(kind=wp), intent(in) :: i_cdrhodt
real(kind=wp), intent(in) :: cdrho
real(kind=wp), intent(in) :: dt
integer, intent(in) :: nx_phys
integer, intent(in) :: ny_phys
integer, intent(in) :: nghost
integer, intent(in) :: nx
integer, intent(in) :: ny
real(kind=wp), intent(in) :: a_u(nx+1,ny)

PR 62: face ice concentration. Valid ONLY when a_face_on.

logical, intent(in) :: a_face_on

private pure subroutine evp_v_momentum_impl(idyCv, idxCv, dx2h, dy2q, iareaCv, mask_v, mi_v, mi_u, q, str_d, str_t, str_s, uo, vo, tau_ay, u_tmp, vi, fyoc, m_neglect, i_cdrhodt, cdrho, dt, nx_phys, ny_phys, nghost, nx, ny, a_v, a_face_on)

v-momentum (:1257-1334, mirror of u). D4: reads u_tmp (the PRE-update u), never the just-updated ui. Minus on the str_t divergence term (:1263-1267).

Read more…

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: idyCv(nx,ny+1)
real(kind=wp), intent(in) :: idxCv(nx,ny+1)
real(kind=wp), intent(in) :: dx2h(nx,ny)
real(kind=wp), intent(in) :: dy2q(nx+1,ny+1)
real(kind=wp), intent(in) :: iareaCv(nx,ny+1)
real(kind=wp), intent(in) :: mask_v(nx,ny+1)
real(kind=wp), intent(in) :: mi_v(nx,ny+1)
real(kind=wp), intent(in) :: mi_u(nx+1,ny)
real(kind=wp), intent(in) :: q(nx+1,ny+1)
real(kind=wp), intent(in) :: str_d(nx,ny)
real(kind=wp), intent(in) :: str_t(nx,ny)
real(kind=wp), intent(in) :: str_s(nx+1,ny+1)
real(kind=wp), intent(in) :: uo(nx+1,ny)
real(kind=wp), intent(in) :: vo(nx,ny+1)
real(kind=wp), intent(in) :: tau_ay(nx,ny+1)
real(kind=wp), intent(in) :: u_tmp(nx+1,ny)
real(kind=wp), intent(inout) :: vi(nx,ny+1)
real(kind=wp), intent(inout) :: fyoc(nx,ny+1)
real(kind=wp), intent(in) :: m_neglect
real(kind=wp), intent(in) :: i_cdrhodt
real(kind=wp), intent(in) :: cdrho
real(kind=wp), intent(in) :: dt
integer, intent(in) :: nx_phys
integer, intent(in) :: ny_phys
integer, intent(in) :: nghost
integer, intent(in) :: nx
integer, intent(in) :: ny
real(kind=wp), intent(in) :: a_v(nx,ny+1)

PR 62: face ice concentration. Valid ONLY when a_face_on.

logical, intent(in) :: a_face_on

private pure subroutine evp_wrap_corner_impl(fld, nx_face, ny_face, nx_phys, ny_phys, nghost, wrap_x, wrap_y)

Periodic ghost-wrap for a corner-staggered field (e.g. str_s), shape (nx_total+1, ny_total+1). No corner-wrap helper exists in rdb_ocean_periodic (only centre/face_x/face_y) — this is the EVP-local twin, same two-pass (x-then-y) structure.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(inout) :: fld(nx_face,ny_face)
integer, intent(in) :: nx_face
integer, intent(in) :: ny_face
integer, intent(in) :: nx_phys
integer, intent(in) :: ny_phys
integer, intent(in) :: nghost
logical, intent(in) :: wrap_x
logical, intent(in) :: wrap_y

private pure subroutine evp_zero_massless_velocity_impl(mask_u, mask_v, mis, ui, vi, nx, ny)

SIS2 :899-907 — zero ice velocities where BOTH neighbouring cells are massless (or the face is masked/land).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: mask_u(nx+1,ny)
real(kind=wp), intent(in) :: mask_v(nx,ny+1)
real(kind=wp), intent(in) :: mis(nx,ny)
real(kind=wp), intent(inout) :: ui(nx+1,ny)
real(kind=wp), intent(inout) :: vi(nx,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_zero_stress_impl(fxoc, fyoc, nx, ny)

Zero the subcycle-averaged ice->ocean stress accumulators via an explicit do concurrent device kernel (F1): fxoc/fyoc are copyin-mapped device-resident arrays, so a host = 0.0_wp would zero only the HOST copy and leave the device copy carrying the prior call’s average (the accumulate below would then converge to S/(N-1) instead of S/N). Runs on the device-present arrays; inert no-op on host builds.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(out) :: fxoc(nx+1,ny)
real(kind=wp), intent(out) :: fyoc(nx,ny+1)
integer, intent(in) :: nx
integer, intent(in) :: ny

private pure subroutine evp_zeta_impl(sh_dd, sh_dt, sh_ds, i_ec2, pres_mice, mice, del_sh_min_pr, del_sh, zeta, nx, ny)

del_sh / zeta (:1082-1095). shear_at_T averages the 4 surrounding corner sh_Ds values.

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: sh_dd(nx,ny)
real(kind=wp), intent(in) :: sh_dt(nx,ny)
real(kind=wp), intent(in) :: sh_ds(nx+1,ny+1)
real(kind=wp), intent(in) :: i_ec2
real(kind=wp), intent(in) :: pres_mice(nx,ny)
real(kind=wp), intent(in) :: mice(nx,ny)
real(kind=wp), intent(in) :: del_sh_min_pr(nx,ny)
real(kind=wp), intent(out) :: del_sh(nx,ny)
real(kind=wp), intent(out) :: zeta(nx,ny)
integer, intent(in) :: nx
integer, intent(in) :: ny

private subroutine ice_evp_dynamics_impl(grid, metrics, f_corner, mis, mice, ci, uo, vo, tau_ax, tau_ay, ui, vi, str_d, str_t, str_s, fxoc, fyoc, dt_slow, par, periodic_x, periodic_y, nx, ny, mis_w, mice_w, ci_w, pres_mice_w, del_sh_min_pr_w, sh_dd_w, sh_dt_w, zeta_w, del_sh_w, mask_t_w, mi_u_w, mask_u_w, u_tmp_w, mi_v_w, mask_v_w, a_u_w, a_v_w, sh_ds_w, mi_ratio_a_q_w, q_w, mask_q_w, halo_x, halo_y, land_w, land_e, land_s, land_n, dt_transport, n_trunc)

Flat-impl core of ice_evp_dynamics: the EVP subcycle body with the persistent scratch passed as EXPLICIT-SHAPE dummies (memory: never assumed-shape into a do concurrent feeder — NVHPC would emit descriptor-walk memcpys per launch). The *_w scratch names mirror the retired module workspace 1:1, so the body below is unchanged from the pre-slot version.

Read more…

Arguments

Type IntentOptional Attributes Name
type(hgrid_t), intent(in) :: grid
type(ocean_metrics_t), intent(in) :: metrics
real(kind=wp), intent(in) :: f_corner(:,:)
real(kind=wp), intent(in) :: mis(:,:)
real(kind=wp), intent(in) :: mice(:,:)
real(kind=wp), intent(in) :: ci(:,:)
real(kind=wp), intent(in) :: uo(:,:)
real(kind=wp), intent(in) :: vo(:,:)
real(kind=wp), intent(in) :: tau_ax(:,:)
real(kind=wp), intent(in) :: tau_ay(:,:)
real(kind=wp), intent(inout) :: ui(:,:)
real(kind=wp), intent(inout) :: vi(:,:)
real(kind=wp), intent(inout) :: str_d(:,:)
real(kind=wp), intent(inout) :: str_t(:,:)
real(kind=wp), intent(inout) :: str_s(:,:)
real(kind=wp), intent(inout) :: fxoc(:,:)
real(kind=wp), intent(inout) :: fyoc(:,:)
real(kind=wp), intent(in) :: dt_slow
type(ice_evp_params_t), intent(in) :: par
logical, intent(in) :: periodic_x
logical, intent(in) :: periodic_y
integer, intent(in) :: nx

T-cell extents (declared before the explicit-shape scratch that uses them — decl-order rule).

integer, intent(in) :: ny

T-cell extents (declared before the explicit-shape scratch that uses them — decl-order rule).

real(kind=wp), intent(inout) :: mis_w(nx,ny)
real(kind=wp), intent(inout) :: mice_w(nx,ny)
real(kind=wp), intent(inout) :: ci_w(nx,ny)
real(kind=wp), intent(inout) :: pres_mice_w(nx,ny)
real(kind=wp), intent(inout) :: del_sh_min_pr_w(nx,ny)
real(kind=wp), intent(inout) :: sh_dd_w(nx,ny)
real(kind=wp), intent(inout) :: sh_dt_w(nx,ny)
real(kind=wp), intent(inout) :: zeta_w(nx,ny)
real(kind=wp), intent(inout) :: del_sh_w(nx,ny)
real(kind=wp), intent(inout) :: mask_t_w(nx,ny)
real(kind=wp), intent(inout) :: mi_u_w(nx+1,ny)
real(kind=wp), intent(inout) :: mask_u_w(nx+1,ny)
real(kind=wp), intent(inout) :: u_tmp_w(nx+1,ny)
real(kind=wp), intent(inout) :: mi_v_w(nx,ny+1)
real(kind=wp), intent(inout) :: mask_v_w(nx,ny+1)
real(kind=wp), intent(inout) :: a_u_w(nx+1,ny)

PR 62: face ice concentration, valid ONLY when par%a_face_stress (uninitialised device memory otherwise — never read off-gate).

real(kind=wp), intent(inout) :: a_v_w(nx,ny+1)

PR 62: face ice concentration, valid ONLY when par%a_face_stress (uninitialised device memory otherwise — never read off-gate).

real(kind=wp), intent(inout) :: sh_ds_w(nx+1,ny+1)
real(kind=wp), intent(inout) :: mi_ratio_a_q_w(nx+1,ny+1)
real(kind=wp), intent(inout) :: q_w(nx+1,ny+1)
real(kind=wp), intent(inout) :: mask_q_w(nx+1,ny+1)
logical, intent(in) :: halo_x

Axis split across ranks (see ice_evp_dynamics).

logical, intent(in) :: halo_y

Axis split across ranks (see ice_evp_dynamics).

logical, intent(in) :: land_w

Ghost band pinned to land (see ice_evp_dynamics).

logical, intent(in) :: land_e

Ghost band pinned to land (see ice_evp_dynamics).

logical, intent(in) :: land_s

Ghost band pinned to land (see ice_evp_dynamics).

logical, intent(in) :: land_n

Ghost band pinned to land (see ice_evp_dynamics).

real(kind=wp), intent(in), optional :: dt_transport

PR 36: the dt TRANSPORT will actually consume (ocean_dyn%therm_dt(dt)), NOT this call’s dt_slow — EVP runs every outer step, transport at thermo cadence (rdb_driver.F90). Absent => dt_slow (SIS2’s own assumption: SIS_C_dynamics and SIS_transport share dt_slow). Read only when par%cfl_trunc > 0.

integer, intent(out), optional :: n_trunc

PR 36: count of ice-bearing faces (mi_u/mi_v > m_neglect) the FINAL clip touched. 0 when par%cfl_trunc <= 0.

private pure subroutine ice_limit_stresses(areaT, mask_t, pres_mice, mice, str_d, str_t, str_s, ec, nx, ny)

SIS2 limit_stresses (:1619-1684), lim=1 (no optional arg). Called ONCE per ice_evp_dynamics call, BEFORE the substep loop — requirement (2). Corner clamp uses the MASKED-area-weighted mean pressure of the <=4 wet neighbours — requirement (3).

Arguments

Type IntentOptional Attributes Name
real(kind=wp), intent(in) :: areaT(nx,ny)
real(kind=wp), intent(in) :: mask_t(nx,ny)
real(kind=wp), intent(in) :: pres_mice(nx,ny)
real(kind=wp), intent(in) :: mice(nx,ny)
real(kind=wp), intent(inout) :: str_d(nx,ny)
real(kind=wp), intent(inout) :: str_t(nx,ny)
real(kind=wp), intent(inout) :: str_s(nx+1,ny+1)
real(kind=wp), intent(in) :: ec
integer, intent(in) :: nx
integer, intent(in) :: ny