Skip to content

Use SST and SSS consistent with MARBL tracers #374

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 4 commits into from
Aug 26, 2025
Merged

Use SST and SSS consistent with MARBL tracers #374

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
930af7f
Use prediabatic SST and SSS for MARBL surface flux
mnlevy1981 Jul 17, 2025
c4bbe6b
Clean up comments
mnlevy1981 Jul 17, 2025
92e33be
Use prediabatic T & S for MARBL interior tendency
mnlevy1981 Jul 17, 2025
f75be1d
Clean up following code review
mnlevy1981 Jul 30, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -387,7 +387,7 @@ subroutine ocean_model_init(Ocean_sfc, OS, Time_init, Time_in, gas_fields_ocn, i
! vertical integrals, since the related 3-d sums are not negligible in cost.
call allocate_surface_state(OS%sfc_state, OS%grid, use_temperature, &
do_integrals=.true., gas_fields_ocn=gas_fields_ocn, &
use_meltpot=use_melt_pot, use_marbl_tracers=OS%use_MARBL)
use_meltpot=use_melt_pot, use_MARBL_tracers=OS%use_MARBL)

call surface_forcing_init(Time_in, OS%grid, OS%US, param_file, OS%diag, &
OS%forcing_CSp, OS%restore_salinity, OS%restore_temp, OS%use_waves)
Expand Down
10 changes: 5 additions & 5 deletions config_src/drivers/nuopc_cap/mom_surface_forcing_nuopc.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -81,7 +81,7 @@ module MOM_surface_forcing_nuopc
!! pressure limited by max_p_surf instead of the
!! full atmospheric pressure. The default is true.
logical :: use_CFC !< enables the MOM_CFC_cap tracer package.
logical :: use_marbl_tracers !< enables the MARBL tracer package.
logical :: use_MARBL_tracers !< enables the MARBL tracer package.
logical :: enthalpy_cpl !< Controls if enthalpy terms are provided by the coupler or computed
!! internally.
real :: gust_const !< constant unresolved background gustiness for ustar [R L Z T-2 ~> Pa]
Expand Down Expand Up @@ -326,7 +326,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,
if (fluxes%dt_buoy_accum < 0) then
call allocate_forcing_type(G, fluxes, water=.true., heat=.true., ustar=.true., &
press=.true., fix_accum_bug=.not.CS%ustar_gustless_bug, &
cfc=CS%use_CFC, marbl=CS%use_marbl_tracers, hevap=CS%enthalpy_cpl, &
cfc=CS%use_CFC, marbl=CS%use_MARBL_tracers, hevap=CS%enthalpy_cpl, &
tau_mag=.true., ice_ncat=IOB%ice_ncat)
call safe_alloc_ptr(fluxes%omega_w2x,isd,ied,jsd,jed)
call safe_alloc_ptr(fluxes%sw_vis_dir,isd,ied,jsd,jed)
Expand Down Expand Up @@ -615,7 +615,7 @@ subroutine convert_IOB_to_fluxes(IOB, fluxes, index_bounds, Time, valid_time, G,

! Copy MARBL-specific IOB fields into fluxes; also set some MARBL-specific forcings to other values
! (constants, values from netCDF, etc)
if (CS%use_marbl_tracers) &
if (CS%use_MARBL_tracers) &
call convert_driver_fields_to_forcings(IOB%atm_fine_dust_flux, IOB%atm_coarse_dust_flux, &
IOB%seaice_dust_flux, IOB%atm_bc_flux, IOB%seaice_bc_flux, &
IOB%nhx_dep, IOB%noy_dep, IOB%atm_co2_prog, IOB%atm_co2_diag, &
Expand Down Expand Up @@ -1272,7 +1272,7 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, restore_salt,
call get_param(param_file, mdl, "USE_CFC_CAP", CS%use_CFC, &
default=.false., do_not_log=.true.)

call get_param(param_file, mdl, "USE_MARBL_TRACERS", CS%use_marbl_tracers, &
call get_param(param_file, mdl, "USE_MARBL_TRACERS", CS%use_MARBL_tracers, &
default=.false., do_not_log=.true.)

call get_param(param_file, mdl, "ENTHALPY_FROM_COUPLER", CS%enthalpy_cpl, &
Expand Down Expand Up @@ -1483,7 +1483,7 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, restore_salt,
endif

! Set up MARBL forcing control structure
call MARBL_forcing_init(G, US, param_file, diag, Time, CS%inputdir, CS%use_marbl_tracers, &
call MARBL_forcing_init(G, US, param_file, diag, Time, CS%inputdir, CS%use_MARBL_tracers, &
CS%marbl_forcing_CSp)

if (present(restore_salt)) then ; if (restore_salt) then
Expand Down
10 changes: 5 additions & 5 deletions config_src/drivers/solo_driver/MOM_surface_forcing.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -116,7 +116,7 @@ module MOM_surface_forcing
!! rotationally invariant and more likely to be the same between compilers.
logical :: ustar_gustless_bug !< If true, include a bug in the time-averaging of the
!! gustless wind friction velocity.
logical :: use_marbl_tracers !< If true, allocate memory for forcing needed by MARBL
logical :: use_MARBL_tracers !< If true, allocate memory for forcing needed by MARBL
! if WIND_CONFIG=='scurves' then use the following to define a piecewise scurve profile
real :: scurves_ydata(20) = 90. !< Latitudes of scurve nodes [degreesN]
real :: scurves_taux(20) = 0. !< Zonal wind stress values at scurve nodes [R L Z T-2 ~> Pa]
Expand Down Expand Up @@ -287,7 +287,7 @@ subroutine set_forcing(sfc_state, forces, fluxes, day_start, day_interval, G, US
! Allocate memory for the mechanical and thermodynamic forcing fields.
call allocate_mech_forcing(G, forces, stress=.true., ustar=.not.CS%nonBous, press=.true., tau_mag=CS%nonBous)

call allocate_forcing_type(G, fluxes, ustar=.not.CS%nonBous, marbl=CS%use_marbl_tracers, tau_mag=CS%nonBous, &
call allocate_forcing_type(G, fluxes, ustar=.not.CS%nonBous, marbl=CS%use_MARBL_tracers, tau_mag=CS%nonBous, &
fix_accum_bug=.not.CS%ustar_gustless_bug)
if (trim(CS%buoy_config) /= "NONE") then
if ( CS%use_temperature ) then
Expand Down Expand Up @@ -384,7 +384,7 @@ subroutine set_forcing(sfc_state, forces, fluxes, day_start, day_interval, G, US
endif
endif

if (CS%use_marbl_tracers) then
if (CS%use_MARBL_tracers) then
call MARBL_forcing_from_data_override(fluxes, day_center, G, US, CS)
endif

Expand Down Expand Up @@ -2148,7 +2148,7 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, tracer_flow_C
call read_netCDF_data(filename, 'gustiness', CS%gust, G%Domain, &
rescale=US%Pa_to_RLZ_T2*US%L_to_Z) ! units in file should be [Pa]
endif
call get_param(param_file, mdl, "USE_MARBL_TRACERS", CS%use_marbl_tracers, &
call get_param(param_file, mdl, "USE_MARBL_TRACERS", CS%use_MARBL_tracers, &
default=.false., do_not_log=.true.)

! All parameter settings are now known.
Expand Down Expand Up @@ -2181,7 +2181,7 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, tracer_flow_C
endif

! Set up MARBL forcing control structure
call MARBL_forcing_init(G, US, param_file, diag, Time, CS%inputdir, CS%use_marbl_tracers, &
call MARBL_forcing_init(G, US, param_file, diag, Time, CS%inputdir, CS%use_MARBL_tracers, &
CS%marbl_forcing_CSp)

call register_forcing_type_diags(Time, diag, US, CS%use_temperature, CS%handles)
Expand Down
2 changes: 1 addition & 1 deletion src/core/MOM_forcing_type.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -3524,7 +3524,7 @@ subroutine allocate_forcing_by_group(G, fluxes, water, heat, ustar, press, &

if (present(fix_accum_bug)) fluxes%gustless_accum_bug = .not.fix_accum_bug

!These fields should only be allocated when USE_MARBL is activated.
!These fields should only be allocated when USE_MARBL_TRACERS is activated.
call myAlloc(fluxes%ice_fraction,isd,ied,jsd,jed, marbl)
call myAlloc(fluxes%u10_sqr,isd,ied,jsd,jed, marbl)
call myAlloc(fluxes%noy_dep,isd,ied,jsd,jed, marbl)
Expand Down
6 changes: 3 additions & 3 deletions src/core/MOM_variables.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -364,7 +364,7 @@ module MOM_variables
!! the ocean model. Unused fields are unallocated.
subroutine allocate_surface_state(sfc_state, G, use_temperature, do_integrals, &
gas_fields_ocn, use_meltpot, use_iceshelves, &
omit_frazil, sfc_state_in, turns, use_marbl_tracers)
omit_frazil, sfc_state_in, turns, use_MARBL_tracers)
type(ocean_grid_type), intent(in) :: G !< ocean grid structure
type(surface), intent(inout) :: sfc_state !< ocean surface state type to be allocated.
logical, optional, intent(in) :: use_temperature !< If true, allocate the space for thermodynamic variables.
Expand All @@ -390,7 +390,7 @@ subroutine allocate_surface_state(sfc_state, G, use_temperature, do_integrals, &
!! is present, it is used and tr_fields_in is ignored.
integer, optional, intent(in) :: turns !< If present, the number of counterclockwise quarter
!! turns to use on the new grid.
logical, optional, intent(in) :: use_marbl_tracers !< If true, allocate the space for CO2 flux from MARBL
logical, optional, intent(in) :: use_MARBL_tracers !< If true, allocate the space for CO2 flux from MARBL

! local variables
logical :: use_temp, alloc_integ, use_melt_potential, alloc_iceshelves, alloc_frazil, alloc_fco2
Expand All @@ -408,7 +408,7 @@ subroutine allocate_surface_state(sfc_state, G, use_temperature, do_integrals, &
use_melt_potential = .false. ; if (present(use_meltpot)) use_melt_potential = use_meltpot
alloc_iceshelves = .false. ; if (present(use_iceshelves)) alloc_iceshelves = use_iceshelves
alloc_frazil = .true. ; if (present(omit_frazil)) alloc_frazil = .not.omit_frazil
alloc_fco2 = .false. ; if (present(use_marbl_tracers)) alloc_fco2 = use_marbl_tracers
alloc_fco2 = .false. ; if (present(use_MARBL_tracers)) alloc_fco2 = use_MARBL_tracers

if (sfc_state%arrays_allocated) return

Expand Down
56 changes: 51 additions & 5 deletions src/parameterizations/vertical/MOM_diabatic_driver.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -67,7 +67,7 @@ module MOM_diabatic_driver
use MOM_sponge, only : apply_sponge, sponge_CS
use MOM_ALE_sponge, only : apply_ALE_sponge, ALE_sponge_CS
use MOM_time_manager, only : time_type, real_to_time, operator(-), operator(<=)
use MOM_tracer_flow_control, only : call_tracer_column_fns, tracer_flow_control_CS
use MOM_tracer_flow_control, only : call_tracer_column_fns, tracer_flow_control_CS, extract_tracer_flow_member
use MOM_tracer_diabatic, only : tracer_vertdiff, tracer_vertdiff_Eulerian
use MOM_unit_scaling, only : unit_scale_type
use MOM_variables, only : thermo_var_ptrs, vertvisc_type, accel_diag_ptrs
Expand Down Expand Up @@ -183,6 +183,10 @@ module MOM_diabatic_driver
logical :: Use_KdWork_diag = .false. !< Logical flag to indicate if any Kd_work diagnostics are on.
logical :: Use_N2_diag = .false. !< Logical flag to indicate if any N2 diagnostics are on.

! MARBL needs T & S from before the tracer_vertdiff call
real, allocatable, dimension(:,:,:) :: prediabatic_T !< Temperature prior to calling diabatic driver [C ~> degC]
real, allocatable, dimension(:,:,:) :: prediabatic_S !< Salinity prior to calling diabatic driver [S ~> ppt]

!>@{ Diagnostic IDs
integer :: id_ea = -1, id_eb = -1 ! used by layer diabatic
integer :: id_ea_t = -1, id_eb_t = -1, id_ea_s = -1, id_eb_s = -1
Expand Down Expand Up @@ -637,6 +641,22 @@ subroutine diabatic_ALE_legacy(u, v, h, tv, BLD, fluxes, visc, ADp, CDp, dt, Tim
showCallTree = callTree_showQuery()
if (showCallTree) call callTree_enter("diabatic_ALE_legacy(), MOM_diabatic_driver.F90")

! Some tracer packages require T & S from the beginning of the diabatic step to
! provide forcing consistent with the passive tracer values. The initialization
! routine will allocate prediabatic_T and prediabatic_S if the tracer flow control
! structure indicates it is necessary. If these arrays are allocated, they will store
! a copy of tv%T & tv%S before this subroutine modifies the tv structure.
if (allocated(CS%prediabatic_T)) then
do k=1,nz ; do j=js,je ; do i=is,ie
CS%prediabatic_T(i,j,k) = tv%T(i,j,k)
enddo ; enddo ; enddo
endif
if (allocated(CS%prediabatic_S)) then
do k=1,nz ; do j=js,je ; do i=is,ie
CS%prediabatic_S(i,j,k) = tv%S(i,j,k)
enddo ; enddo ; enddo
endif

! For all other diabatic subroutines, the averaging window should be the entire diabatic timestep
call enable_averages(dt, Time_end, CS%diag)

Expand Down Expand Up @@ -1192,7 +1212,8 @@ subroutine diabatic_ALE_legacy(u, v, h, tv, BLD, fluxes, visc, ADp, CDp, dt, Tim
KPP_CSp=CS%KPP_CSp, &
nonLocalTrans=KPP_NLTscalar, &
evap_CFL_limit=CS%evap_CFL_limit, &
minimum_forcing_depth=CS%minimum_forcing_depth, h_BL=visc%h_ML)
minimum_forcing_depth=CS%minimum_forcing_depth, h_BL=visc%h_ML, &
prediabatic_T=CS%prediabatic_T, prediabatic_S=CS%prediabatic_S)

call cpu_clock_end(id_clock_tracers)

Expand Down Expand Up @@ -1352,6 +1373,22 @@ subroutine diabatic_ALE(u, v, h, tv, BLD, fluxes, visc, ADp, CDp, dt, Time_end,
if (.not. (CS%useALEalgorithm)) call MOM_error(FATAL, "MOM_diabatic_driver: "// &
"The ALE algorithm must be enabled when using MOM_diabatic_driver.")

! Some tracer packages require T & S from the beginning of the diabatic step to
! provide forcing consistent with the passive tracer values. The initialization
! routine will allocate prediabatic_T and prediabatic_S if the tracer flow control
! structure indicates it is necessary. If these arrays are allocated, they will store
! a copy of tv%T & tv%S before this subroutine modifies the tv structure.
if (allocated(CS%prediabatic_T)) then
do k=1,nz ; do j=js,je ; do i=is,ie
CS%prediabatic_T(i,j,k) = tv%T(i,j,k)
enddo ; enddo ; enddo
endif
if (allocated(CS%prediabatic_S)) then
do k=1,nz ; do j=js,je ; do i=is,ie
CS%prediabatic_S(i,j,k) = tv%S(i,j,k)
enddo ; enddo ; enddo
endif

! For all other diabatic subroutines, the averaging window should be the entire diabatic timestep
call enable_averages(dt, Time_end, CS%diag)

Expand Down Expand Up @@ -1818,7 +1855,8 @@ subroutine diabatic_ALE(u, v, h, tv, BLD, fluxes, visc, ADp, CDp, dt, Time_end,
KPP_CSp=CS%KPP_CSp, &
nonLocalTrans=KPP_NLTscalar, &
evap_CFL_limit=CS%evap_CFL_limit, &
minimum_forcing_depth=CS%minimum_forcing_depth, h_BL=visc%h_ML)
minimum_forcing_depth=CS%minimum_forcing_depth, h_BL=visc%h_ML, &
prediabatic_T=CS%prediabatic_T, prediabatic_S=CS%prediabatic_S)

call cpu_clock_end(id_clock_tracers)

Expand Down Expand Up @@ -3227,7 +3265,7 @@ subroutine diabatic_driver_init(Time, G, GV, US, param_file, useALEalgorithm, di

! Local variables
real :: Kd ! A diffusivity used in the default for other tracer diffusivities [Z2 T-1 ~> m2 s-1]
logical :: use_temperature
logical :: use_temperature, use_MARBL_tracers
character(len=20) :: EN1, EN2, EN3

! This "include" declares and sets the variable "version".
Expand All @@ -3246,7 +3284,12 @@ subroutine diabatic_driver_init(Time, G, GV, US, param_file, useALEalgorithm, di
CS%diag => diag
CS%Time => Time

if (associated(tracer_flow_CSp)) CS%tracer_flow_CSp => tracer_flow_CSp
if (associated(tracer_flow_CSp)) then
CS%tracer_flow_CSp => tracer_flow_CSp
call extract_tracer_flow_member(tracer_flow_CSp, use_MARBL_tracers=use_MARBL_tracers)
if (use_MARBL_tracers) &
allocate(CS%prediabatic_T(SZI_(G),SZJ_(G), SZK_(G)), CS%prediabatic_S(SZI_(G),SZJ_(G), SZK_(G)))
end if
if (associated(sponge_CSp)) CS%sponge_CSp => sponge_CSp
if (associated(ALE_sponge_CSp)) CS%ALE_sponge_CSp => ALE_sponge_CSp
if (associated(oda_incupd_CSp)) CS%oda_incupd_CSp => oda_incupd_CSp
Expand Down Expand Up @@ -3842,6 +3885,9 @@ subroutine diabatic_driver_end(CS)
deallocate(CS%optics)
endif

if (allocated(CS%prediabatic_T)) deallocate(CS%prediabatic_T)
if (allocated(CS%prediabatic_S)) deallocate(CS%prediabatic_S)

if (CS%debug_energy_req) &
call diapyc_energy_req_end(CS%diapyc_en_rec_CSp)

Expand Down
14 changes: 7 additions & 7 deletions src/tracer/MARBL_forcing_mod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -55,7 +55,7 @@ module MARBL_forcing_mod

type(marbl_forcing_diag_ids) :: diag_ids !< used for registering and posting some MARBL forcing fields as diagnostics

logical :: use_marbl_tracers !< most functions can return immediately
logical :: use_MARBL_tracers !< most functions can return immediately
!! MARBL tracers are turned off
integer :: atm_co2_iopt !< Integer version of atm_co2_opt, which determines source of atm_co2
integer :: atm_alt_co2_iopt !< Integer version of atm_alt_co2_opt, which determines source of atm_alt_co2
Expand All @@ -69,14 +69,14 @@ module MARBL_forcing_mod

contains

subroutine MARBL_forcing_init(G, US, param_file, diag, day, inputdir, use_marbl, CS)
subroutine MARBL_forcing_init(G, US, param_file, diag, day, inputdir, use_MARBL_tracers, CS)
type(ocean_grid_type), intent(in) :: G !< The ocean's grid structure
type(unit_scale_type), intent(in) :: US !< A dimensional unit scaling type
type(param_file_type), intent(in) :: param_file !< A structure to parse for run-time parameters
type(diag_ctrl), target, intent(in) :: diag !< Structure used to regulate diagnostic output.
type(time_type), target, intent(in) :: day !< Time of the start of the run.
character(len=*), intent(in) :: inputdir !< Directory containing input files
logical, intent(in) :: use_marbl !< Is MARBL tracer package active?
logical, intent(in) :: use_MARBL_tracers !< Is MARBL tracer package active?
type(marbl_forcing_CS), pointer, intent(inout) :: CS !< A pointer that is set to point to control
!! structure for MARBL forcing

Expand All @@ -92,9 +92,9 @@ subroutine MARBL_forcing_init(G, US, param_file, diag, day, inputdir, use_marbl,
allocate(CS)
CS%diag => diag

CS%use_marbl_tracers = .true.
if (.not. use_marbl) then
CS%use_marbl_tracers = .false.
CS%use_MARBL_tracers = .true.
if (.not. use_MARBL_tracers) then
CS%use_MARBL_tracers = .false.
return
endif

Expand Down Expand Up @@ -229,7 +229,7 @@ subroutine convert_driver_fields_to_forcings(atm_fine_dust_flux, atm_coarse_dust
real :: ndep_conversion !< Factor to convert nitrogen deposition from kg m-2 s-1 -> mmol m-3 (m s-1)
!! [s m2 kg-1 conc Z T-1 ~> mmol kg-1]

if (.not. CS%use_marbl_tracers) return
if (.not. CS%use_MARBL_tracers) return

is = G%isc ; ie = G%iec ; js = G%jsc ; je = G%jec
ndep_conversion = (1.e6/14.) * (US%m_to_Z * US%T_to_s)
Expand Down
Loading
Loading