mom-ocean · dougiesquire · Jul 12, 2024
diff --git a/src/mom5/ocean_bgc/ocean_generic_tracer.F90 b/src/mom5/ocean_bgc/ocean_generic_tracer.F90
@@ -27,10 +27,10 @@ module ocean_generic_mod
 
   use field_manager_mod, only: fm_get_index,fm_string_len, fm_new_value
   use generic_tracer, only: generic_tracer_init, generic_tracer_source, generic_tracer_update_from_bottom
-  use generic_tracer, only: generic_tracer_coupler_get, generic_tracer_coupler_set, generic_tracer_register_diag
+  use generic_tracer, only: generic_tracer_coupler_accumulate, generic_tracer_coupler_set, generic_tracer_register_diag
   use generic_tracer, only: generic_tracer_end, generic_tracer_get_list, do_generic_tracer, generic_tracer_register
   use generic_tracer, only: generic_tracer_coupler_zero, generic_tracer_vertdiff_G, generic_tracer_vertdiff_M
-  use generic_tracer, only: generic_tracer_diag
+  use generic_tracer, only: generic_tracer_diag, generic_tracer_update_from_coupler
 
   use g_tracer_utils,   only: g_tracer_get_name,g_tracer_get_alias,g_tracer_set_values,g_tracer_get_common
   use g_tracer_utils,   only: g_tracer_get_next,g_tracer_type,g_tracer_is_prog,g_tracer_flux_init
@@ -51,6 +51,7 @@ module ocean_generic_mod
   public ocean_generic_sum_sfc
   public ocean_generic_zero_sfc
   public ocean_generic_sbc
+  public ocean_generic_sbc_adjust
   public ocean_generic_init
   public ocean_generic_flux_init
   public ocean_generic_column_physics
@@ -336,9 +337,15 @@ subroutine ocean_generic_sbc(Ice_ocean_boundary_fluxes,Disd,Djsd, T_prog, runoff
     character(len=fm_string_len)      :: g_tracer_name
     character(len=fm_string_len), parameter :: sub_name = 'update_generic_tracer_sbc'
 
-    !Extract the tracer surface fields from coupler 
-    call generic_tracer_coupler_get(Ice_ocean_boundary_fluxes)
+    !Extract the tracer surface fields from coupler
+    ! dts: change to use generic_tracer_coupler_accumulate with a weight of 1. for consistency
+    ! with MOM6. Note this means that the generic_?_update_from_coupler method is no longer
+    ! called here. Instead this will be called later in ocean_generic_sbc_adjust when flux
+    ! adjustments are available
+    ! call generic_tracer_coupler_get(Ice_ocean_boundary_fluxes)
+    call generic_tracer_coupler_accumulate(Ice_ocean_boundary_fluxes, 1.)
 
+    ! dts: Is this still needed here, since this is now also done in ocean_generic_sbc_adjust?
     !Update T_prog fields from generic tracer fields
     !
     !Get the tracer list
@@ -394,6 +401,87 @@ subroutine ocean_generic_sbc(Ice_ocean_boundary_fluxes,Disd,Djsd, T_prog, runoff
     enddo
   end subroutine ocean_generic_sbc
 
+  ! <SUBROUTINE NAME="ocean_generic_sbc_adjust">
+  !  <OVERVIEW>
+  !   Adjust tracer surface boundary conditions
+  !  </OVERVIEW>
+  !  <DESCRIPTION>
+  !   This subroutine applies any adjustments to the generic tracers surface boundary conditions,
+  !   for example, adjustment of surface tracer fluxes due to salinity restoring.
+  !  </DESCRIPTION>
+  !  <TEMPLATE>
+  !   call ocean_generic_sbc_adjust(Ice_ocean_boundary_fluxes,Disd,Djsd, T_prog )
+  !  </TEMPLATE>
+  ! </SUBROUTINE>
+  subroutine ocean_generic_sbc_adjust(Disd, Djsd, T_prog, salt_flux_added, runoff)
+   integer, intent(in)                                             :: Disd, Djsd
+   type(ocean_prog_tracer_type), dimension(:), intent(inout)       :: T_prog
+   real, intent(in), dimension(Disd:,Djsd:)                        :: salt_flux_added
+   real, intent(in), dimension(Disd:,Djsd:)                        :: runoff
+
+   type(g_tracer_type), pointer    :: g_tracer_list,g_tracer,g_tracer_next
+   integer                           :: g_tracer_index
+   character(len=fm_string_len)      :: g_tracer_name
+   character(len=fm_string_len), parameter :: sub_name = 'ocean_generic_sbc_adjust'
+
+   ! Adjust tracer fields via the generic_?_update_from_coupler method
+   call generic_tracer_update_from_coupler(Disd, Djsd, salt_flux_added)
+
+   !Update T_prog fields from generic tracer fields
+   !
+   !Get the tracer list
+   call generic_tracer_get_list(g_tracer_list)
+   if(.NOT. associated(g_tracer_list)) call mpp_error(FATAL, trim(sub_name)//&
+        ": No tracer in the list.")
+   !For each tracer name get its T_prog index and get its flux fields
+   g_tracer=>g_tracer_list
+   do
+      if(g_tracer_is_prog(g_tracer)) then
+         call g_tracer_get_alias(g_tracer,g_tracer_name)
+         g_tracer_index = fm_get_index(trim('/ocean_mod/prog_tracers/'//g_tracer_name))
+         if (g_tracer_index .le. 0) &
+              call mpp_error(FATAL,trim(sub_name) // ' Could not get the index for '//g_tracer_name)
+
+         if (_ALLOCATED(g_tracer%stf) )&
+              call g_tracer_get_values(g_tracer,g_tracer_name,'stf',   T_prog(g_tracer_index)%stf,   Disd,Djsd)
+
+         if (_ALLOCATED(g_tracer%btf) )&
+              call g_tracer_get_values(g_tracer,g_tracer_name,'btf',   T_prog(g_tracer_index)%btf,   Disd,Djsd)
+
+         !If the tracer has runoff fill in the T_prog(n)%trunoff and  T_prog(n)%runoff_tracer_flux
+         if(_ALLOCATED(g_tracer%trunoff)) then
+            !Fill in T_prog(n)%trunoff
+
+            call g_tracer_get_values(g_tracer,g_tracer_name,'trunoff',T_prog(g_tracer_index)%trunoff,Disd,Djsd)
+
+            !Fill in T_prog(n)%runoff_tracer_flux
+            T_prog(g_tracer_index)%runoff_tracer_flux = T_prog(g_tracer_index)%trunoff * runoff
+
+            !Set g_tracer%runoff_tracer_flux
+            call g_tracer_set_values(g_tracer,g_tracer_name,'runoff_tracer_flux',T_prog(g_tracer_index)%runoff_tracer_flux,Disd,Djsd)
+            !
+            !Fill in T_prog(n)%triver in MOM
+            !Note: This is done so that MOM can apply the river fluxes through setting either
+            !        the runoff and calving fluxes (when discharge_combine_runoff_calve=.false.)
+            !      or
+            !        the total river concentration (when discharge_combine_runoff_calve=.true.)
+            !
+            !Assume zero calving flux for the generic tracers.
+            !T_prog(g_tracer_index)%tcalving = 0 !MOM default
+            T_prog(g_tracer_index)%triver = T_prog(g_tracer_index)%trunoff
+
+         endif
+
+      endif
+
+      !traverse the linked list till hit NULL
+      call g_tracer_get_next(g_tracer, g_tracer_next)
+      if(.NOT. associated(g_tracer_next)) exit
+      g_tracer=>g_tracer_next
+
+   enddo
+ end subroutine ocean_generic_sbc_adjust
+
   ! <SUBROUTINE NAME="ocean_generic_column_physics">
   !  <OVERVIEW>
   !   Column physics for generic tracers.

diff --git a/src/mom5/ocean_core/ocean_model.F90 b/src/mom5/ocean_core/ocean_model.F90
@@ -1605,7 +1605,7 @@ subroutine update_ocean_model(Ice_ocean_boundary, Ocean_state, Ocean_sfc, &
        ! compute "flux adjustments" (e.g., surface tracer restoring, flux correction)
        call mpp_clock_begin(id_flux_adjust)
        call flux_adjust(Time, T_diag(1:num_diag_tracers), Dens, Ext_mode, &
-                        T_prog(1:num_prog_tracers), Velocity, river, melt, pme)
+                        T_prog(1:num_prog_tracers), Velocity, river, melt, pme, runoff)
 
        call mpp_clock_end(id_flux_adjust)
 

diff --git a/src/mom5/ocean_core/ocean_sbc.F90 b/src/mom5/ocean_core/ocean_sbc.F90
@@ -522,7 +522,7 @@ module ocean_sbc_mod
 use ocean_parameters_mod,     only: MOM_BGRID, MOM_CGRID 
 use ocean_riverspread_mod,    only: spread_river_horz
 use ocean_tempsalt_mod,       only: pottemp_from_contemp
-use ocean_tpm_mod,            only: ocean_tpm_sum_sfc, ocean_tpm_avg_sfc, ocean_tpm_sbc
+use ocean_tpm_mod,            only: ocean_tpm_sum_sfc, ocean_tpm_avg_sfc, ocean_tpm_sbc, ocean_tpm_sbc_adjust
 use ocean_tpm_mod,            only: ocean_tpm_zero_sfc, ocean_tpm_sfc_end
 use ocean_types_mod,          only: ocean_grid_type, ocean_domain_type, ocean_public_type
 use ocean_types_mod,          only: ocean_time_type, ocean_thickness_type 
@@ -4466,7 +4466,7 @@ end subroutine get_ocean_sbc
 ! </DESCRIPTION>
 !
 
-subroutine flux_adjust(Time, T_diag, Dens, Ext_mode, T_prog, Velocity, river, melt, pme)
+subroutine flux_adjust(Time, T_diag, Dens, Ext_mode, T_prog, Velocity, river, melt, pme, runoff)
 #if defined(ACCESS_CM) || defined(ACCESS_OM)
 
   use auscom_ice_parameters_mod, only : use_ioaice, aice_cutoff
@@ -4482,6 +4482,7 @@ subroutine flux_adjust(Time, T_diag, Dens, Ext_mode, T_prog, Velocity, river, me
   real, dimension(isd:,jsd:),     intent(in)    :: river
   real, dimension(isd:,jsd:),     intent(in)    :: melt 
   real, dimension(isd:,jsd:),     intent(inout) :: pme
+  real, dimension(isd:,jsd:),     intent(inout) :: runoff
 
   real, dimension(isd:ied,jsd:jed) :: open_ocean_mask
   real, dimension(isd:ied,jsd:jed) :: pme_restore, flx_restore
@@ -4766,6 +4767,9 @@ subroutine flux_adjust(Time, T_diag, Dens, Ext_mode, T_prog, Velocity, river, me
 
   endif   ! endif for if (id_correction(index_salt) > 0 )
 
+  ! apply adjustments to tracer surface boundary conditions due to salt restoring/correction
+
+  call ocean_tpm_sbc_adjust(Dom, T_prog, flx_restore(:,:)+flx_correct(:,:), runoff)
 
   ! diagnostics for salinity or pme restoring and correction 
 

diff --git a/src/mom5/ocean_tracers/ocean_tpm.F90 b/src/mom5/ocean_tracers/ocean_tpm.F90
@@ -46,6 +46,10 @@ module ocean_tpm_mod  !{
 !               this functionality may be moved into a new, generalized
 !               boundary condition manager.
 !
+!       ocean_tpm_sbc_adjust: Calls specified routines to adjust
+!               surface boundary condition, e.g. due to virtual fluxes
+!               that arise from salinity restoring
+!
 !       ocean_tpm_bbc: Calls specified routines to handle bottom
 !               coundary condition calculations.
 !
@@ -221,6 +225,7 @@ module ocean_tpm_mod  !{
 use ocean_generic_mod, only: ocean_generic_sum_sfc
 use ocean_generic_mod, only: ocean_generic_zero_sfc
 use ocean_generic_mod, only: ocean_generic_sbc
+use ocean_generic_mod, only: ocean_generic_sbc_adjust
 use ocean_generic_mod, only: ocean_generic_init
 use ocean_generic_mod, only: ocean_generic_column_physics
 use ocean_generic_mod, only: ocean_generic_end
@@ -267,6 +272,7 @@ module ocean_tpm_mod  !{
 public ocean_tpm_init
 public ocean_tpm_flux_init
 public ocean_tpm_sbc
+public ocean_tpm_sbc_adjust
 public ocean_tpm_source
 public ocean_tpm_start
 public ocean_tpm_tracer
@@ -1336,6 +1342,54 @@ subroutine ocean_tpm_sbc(Domain, Grid, T_prog, Time, Ice_ocean_boundary_fluxes,
 end subroutine ocean_tpm_sbc  !}
 ! </SUBROUTINE> NAME="ocean_tpm_sbc"
 
+!#######################################################################
+! <SUBROUTINE NAME="ocean_tpm_sbc_adjust">
+!
+! <DESCRIPTION>
+!       call subroutines to adjust surface boundary condition, e.g. due to
+!       virtual fluxes that arise from salinity restoring
+! </DESCRIPTION>
+!
+
+subroutine ocean_tpm_sbc_adjust(Domain, T_prog, salt_flux_added, runoff)
+
+
+implicit none
+
+!
+!-----------------------------------------------------------------------
+!       Arguments
+!-----------------------------------------------------------------------
+!
+
+type(ocean_domain_type), intent(in)                             :: Domain
+type(ocean_prog_tracer_type), dimension(:), intent(inout)       :: T_prog
+real, dimension(Domain%isd:,Domain%jsd:), intent(in)            :: salt_flux_added
+real, dimension(Domain%isd:,Domain%jsd:), intent(in)            :: runoff
+
+!
+!-----------------------------------------------------------------------
+!     local parameters
+!-----------------------------------------------------------------------
+!
+
+!
+!-----------------------------------------------------------------------
+!     local variables
+!-----------------------------------------------------------------------
+!
+
+
+#ifdef USE_OCEAN_BGC
+
+if (do_generic_tracer) call ocean_generic_sbc_adjust(Domain%isd, Domain%jsd, T_prog, salt_flux_added, runoff)
+
+#endif
+
+return
+
+end subroutine ocean_tpm_sbc_adjust  !}
+! </SUBROUTINE> NAME="ocean_tpm_sbc_adjust"
 
 !#######################################################################
 ! <SUBROUTINE NAME="ocean_tpm_init">

diff --git a/src/ocean_shared/generic_tracers/generic_tracer.F90 b/src/ocean_shared/generic_tracers/generic_tracer.F90
@@ -86,7 +86,9 @@ module generic_tracer
   public generic_tracer_source
   public generic_tracer_diag
   public generic_tracer_update_from_bottom
+  public generic_tracer_update_from_coupler
   public generic_tracer_coupler_get
+  public generic_tracer_coupler_accumulate
   public generic_tracer_coupler_set
   public generic_tracer_coupler_zero
   public generic_tracer_end
@@ -304,6 +306,28 @@ subroutine generic_tracer_coupler_get(IOB_struc)
 
   end subroutine generic_tracer_coupler_get
 
+  ! <SUBROUTINE NAME="generic_tracer_coupler_accumulate">
+  !  <OVERVIEW>
+  !   Accumulate the boundary values (%stf and %triver) from coupler fluxes
+  !  </OVERVIEW>
+  !  <DESCRIPTION>
+  !   This subroutine accumulates coupler values for those generic tracers that have flux
+  !   exchange with atmosphere. This routine simply wraps g_tracer_coupler_get for MOM5,
+  !   but the API is provided here for consistency with https://github.com/NOAA-GFDL/ocean_BGC
+  !  </DESCRIPTION>
+  !  <TEMPLATE>
+  !   call generic_tracer_coupler_accumulate(Ice_ocean_boundary_fluxes, weight)
+  !  </TEMPLATE>
+  ! </SUBROUTINE>
+  subroutine  generic_tracer_coupler_accumulate(IOB_struc, weight)
+    type(coupler_2d_bc_type), intent(in)    :: IOB_struc
+    real,                     intent(in)    :: weight
+
+    !All generic tracers
+    !Update tracer boundary values (%stf and %triver) from coupler fluxes foreach tracer in the prog_tracer_list
+    call g_tracer_coupler_get(tracer_list, IOB_struc)
+
+  end subroutine generic_tracer_coupler_accumulate
 
   ! <SUBROUTINE NAME="generic_tracer_diag">
   !  <OVERVIEW>
@@ -481,6 +505,44 @@ subroutine generic_tracer_update_from_bottom(dt, tau, model_time)
 
   end subroutine generic_tracer_update_from_bottom
 
+  ! <SUBROUTINE NAME="generic_tracer_update_from_coupler">
+  !  <OVERVIEW>
+  !   Modify the values obtained from the coupler
+  !  </OVERVIEW>
+  !  <DESCRIPTION>
+  !   Calls the corresponding generic_X_update_from_coupler routine for each package X.
+  !  </DESCRIPTION>
+  !  <TEMPLATE>
+  !   call generic_tracer_update_from_coupler(ilb, jlb, salt_flux_added)
+  !  </TEMPLATE>
+  !  <IN NAME="ilb,jlb" TYPE="integer">
+  !   Lower bounds of x and y extents of input arrays on data domain
+  !  </IN>
+  !  <IN NAME="salt_flux_added" TYPE="real, dimension(ilb:,jlb:), optional">
+  !   Surface salt flux into ocean from restoring or flux adjustment [g/m^2/sec]
+  !  </IN>
+  ! </SUBROUTINE>
+
+  subroutine generic_tracer_update_from_coupler(ilb, jlb, salt_flux_added)
+     integer,                     intent(in) :: ilb, jlb
+     real, dimension(ilb:,jlb:),  intent(in) :: salt_flux_added
+
+    character(len=fm_string_len), parameter :: sub_name = 'generic_tracer_update_from_coupler'
+
+    !Specific tracers
+    !    if(do_generic_CFC)    call generic_CFC_update_from_coupler(tracer_list) !Nothing to do
+
+    if(do_generic_TOPAZ)  call generic_TOPAZ_update_from_coupler(tracer_list)
+
+    if(do_generic_BLING)  call generic_BLING_update_from_coupler(tracer_list)
+
+    if(do_generic_miniBLING)  call generic_miniBLING_update_from_coupler(tracer_list)
+
+    if(do_generic_COBALT)  call generic_COBALT_update_from_coupler(tracer_list)
+
+    return
+
+   end subroutine generic_tracer_update_from_coupler
 
   ! <SUBROUTINE NAME="generic_tracer_vertdiff_G">
   !  <OVERVIEW>