Fix how fluxes sent to MARBL are treated when DT_THERMO is bigger than coupling interval

config_src/drivers/nuopc_cap/mom_ocean_model_nuopc.F90

@@ -564,9 +564,6 @@ subroutine update_ocean_model(Ice_ocean_boundary, OS, Ocean_sfc, &
 
   if (OS%fluxes%fluxes_used) then
 
-    ! enable_averages() is necessary to post forcing fields to diagnostics
-    call enable_averages(dt_coupling, OS%Time + Ocean_coupling_time_step, OS%diag)
-
     if (do_thermo) &
       call convert_IOB_to_fluxes(Ice_ocean_boundary, OS%fluxes, index_bnds, OS%Time, dt_coupling, &
                                OS%grid, OS%US, OS%forcing_CSp, OS%sfc_state, &

config_src/drivers/nuopc_cap/mom_surface_forcing_nuopc.F90

@@ -1468,7 +1468,8 @@ subroutine surface_forcing_init(Time, G, US, param_file, diag, CS, restore_salt,
 
   call register_forcing_type_diags(Time, diag, US, CS%use_temperature, CS%handles, &
                                    use_berg_fluxes=iceberg_flux_diags, use_waves=use_waves, &
-                                   use_cfcs=CS%use_CFC, use_glc_runoff=glc_runoff_diags)
+                                   use_cfcs=CS%use_CFC, use_MARBL_tracers=CS%use_MARBL_tracers, &
+                                    use_glc_runoff=glc_runoff_diags)
 
   call get_param(param_file, mdl, "ALLOW_FLUX_ADJUSTMENTS", CS%allow_flux_adjustments, &
                  "If true, allows flux adjustments to specified via the "//&

src/core/MOM_forcing_type.F90

@@ -234,7 +234,12 @@ module MOM_forcing_type
     atm_co2 => NULL(),               & !< Atmospheric CO2 Concentration [ppm]
     atm_alt_co2 => NULL(),           & !< Alternate atmospheric CO2 Concentration [ppm]
     dust_flux => NULL(),             & !< Flux of dust into the ocean [R Z T-1 ~> kgN m-2 s-1]
-    iron_flux => NULL()                !< Flux of dust into the ocean [conc Z T-1 ~> conc m s-1]
+    iron_flux => NULL(),             & !< Flux of dust into the ocean [conc Z T-1 ~> conc m s-1]
+    atm_fine_dust_flux => NULL(),    & !< Fine dust flux from atmosphere [R Z T-1 ~> kg m-2 s-1]
+    atm_coarse_dust_flux => NULL(),  & !< Coarse dust flux from atmosphere [R Z T-1 ~> kg m-2 s-1]
+    seaice_dust_flux => NULL(),      & !< Dust flux from seaice [R Z T-1 ~> kg m-2 s-1]
+    atm_bc_flux => NULL(),           & !< Black carbon flux from atmosphere [R Z T-1 ~> kg m-2 s-1]
+    seaice_bc_flux => NULL()           !< Black carbon flux from seaice [R Z T-1 ~> kg m-2 s-1]
 
   real, pointer, dimension(:,:,:) :: &
     fracr_cat   => NULL(),           & !< per-category ice fraction [nondim]
@@ -421,6 +426,11 @@ module MOM_forcing_type
   ! tracer surface flux related diagnostics handles
   integer :: id_ice_fraction = -1
   integer :: id_u10_sqr      = -1
+  integer :: id_atm_fine_dust_flux = -1
+  integer :: id_atm_coarse_dust_flux = -1
+  integer :: id_atm_bc_flux = -1
+  integer :: id_seaice_dust_flux = -1
+  integer :: id_seaice_bc_flux = -1
 
   ! iceberg diagnostic handles
   integer :: id_ustar_berg = -1
@@ -1540,7 +1550,7 @@ end subroutine forcing_SinglePointPrint
 
 !> Register members of the forcing type for diagnostics
 subroutine register_forcing_type_diags(Time, diag, US, use_temperature, handles, use_berg_fluxes, use_waves, &
-                                       use_cfcs, use_glc_runoff)
+                                       use_cfcs, use_MARBL_tracers, use_glc_runoff)
   type(time_type),     intent(in)    :: Time            !< time type
   type(diag_ctrl),     intent(inout) :: diag            !< diagnostic control type
   type(unit_scale_type), intent(in)  :: US              !< A dimensional unit scaling type
@@ -1549,8 +1559,18 @@ subroutine register_forcing_type_diags(Time, diag, US, use_temperature, handles,
   logical, optional,   intent(in)    :: use_berg_fluxes !< If true, allow iceberg flux diagnostics
   logical, optional,   intent(in)    :: use_waves       !< If true, allow wave forcing diagnostics
   logical, optional,   intent(in)    :: use_cfcs        !< If true, allow cfc related diagnostics
+  logical, optional,   intent(in)    :: use_MARBL_tracers  !< If true, allow MARBL related diagnostics
   logical, optional,   intent(in)    :: use_glc_runoff  !< If true, allow separate glacial runoff diagnostics
 
+  logical :: use_cfcs_or_MARBL_tracers
+
+  ! some diagnostics should be registered if either cfc or MARBL tracers are enabled
+  use_cfcs_or_MARBL_tracers = .false.
+  if (present(use_cfcs)) &
+    use_cfcs_or_MARBL_tracers = use_cfcs_or_MARBL_tracers .or. use_cfcs
+  if (present(use_MARBL_tracers)) &
+    use_cfcs_or_MARBL_tracers = use_cfcs_or_MARBL_tracers .or. use_MARBL_tracers
+
   ! Clock for forcing diagnostics
   handles%id_clock_forcing=cpu_clock_id('(Ocean forcing diagnostics)', grain=CLOCK_ROUTINE)
 
@@ -1599,17 +1619,34 @@ subroutine register_forcing_type_diags(Time, diag, US, use_temperature, handles,
     endif
   endif
 
-  ! See:
-  if (present(use_cfcs)) then
-    if (use_cfcs) then
-      handles%id_ice_fraction = register_diag_field('ocean_model', 'ice_fraction', diag%axesT1, Time, &
-          'Fraction of cell area covered by sea ice', 'm2 m-2', conversion=1.0)
+  if (use_cfcs_or_MARBL_tracers) then
+    handles%id_ice_fraction = register_diag_field('ocean_model', 'ice_fraction', diag%axesT1, Time, &
+        'Fraction of cell area covered by sea ice', 'm2 m-2', conversion=1.0)
 
-      handles%id_u10_sqr = register_diag_field('ocean_model', 'u10_sqr', diag%axesT1, Time, &
-          'Wind magnitude at 10m, squared', 'm2 s-2', conversion=US%L_to_m**2*US%s_to_T**2)
-    endif
+    handles%id_u10_sqr = register_diag_field('ocean_model', 'u10_sqr', diag%axesT1, Time, &
+        'Wind magnitude at 10m, squared', 'm2 s-2', conversion=US%L_to_m**2*US%s_to_T**2)
   endif
 
+  if (present(use_MARBL_tracers)) then
+    if (use_MARBL_tracers) then
+      handles%id_atm_fine_dust_flux = register_diag_field('ocean_model', 'ATM_FINE_DUST_FLUX_CPL', &
+          diag%axesT1, Time, 'ATM_FINE_DUST_FLUX from cpl', 'kg m-2 s', &
+          conversion=US%RZ_T_to_kg_m2s)
+      handles%id_atm_coarse_dust_flux = register_diag_field('ocean_model', 'ATM_COARSE_DUST_FLUX_CPL', &
+          diag%axesT1, Time, 'ATM_COARSE_DUST_FLUX from cpl', 'kg m-2 s', &
+          conversion=US%RZ_T_to_kg_m2s)
+      handles%id_atm_bc_flux = register_diag_field('ocean_model', 'ATM_BLACK_CARBON_FLUX_CPL', &
+          diag%axesT1, Time, 'ATM_BLACK_CARBON_FLUX from cpl',  'kg m-2 s', &
+          conversion=US%RZ_T_to_kg_m2s)
+
+      handles%id_seaice_dust_flux = register_diag_field('ocean_model', 'SEAICE_DUST_FLUX_CPL', &
+          diag%axesT1, Time, 'SEAICE_DUST_FLUX from cpl', 'kg m-2 s', &
+          conversion=US%RZ_T_to_kg_m2s)
+      handles%id_seaice_bc_flux = register_diag_field('ocean_model', 'SEAICE_BLACK_CARBON_FLUX_CPL', &
+          diag%axesT1, Time, 'SEAICE_BLACK_CARBON_FLUX from cpl', 'kg m-2 s', &
+          conversion=US%RZ_T_to_kg_m2s)
+    end if
+  end if
   handles%id_psurf = register_diag_field('ocean_model', 'p_surf', diag%axesT1, Time, &
         'Pressure at ice-ocean or atmosphere-ocean interface', &
         'Pa', conversion=US%RL2_T2_to_Pa, cmor_field_name='pso', &
@@ -2338,7 +2375,7 @@ subroutine fluxes_accumulate(flux_tmp, fluxes, G, wt2, forces)
   ! applied based on the time interval stored in flux_tmp.
 
   real :: wt1  ! The relative weight of the previous fluxes [nondim]
-  integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq
+  integer :: i, j, is, ie, js, je, Isq, Ieq, Jsq, Jeq, n
   integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
   is   = G%isc   ; ie   = G%iec    ; js   = G%jsc   ; je   = G%jec
   Isq  = G%IscB  ; Ieq  = G%IecB   ; Jsq  = G%JscB  ; Jeq  = G%JecB
@@ -2499,6 +2536,84 @@ subroutine fluxes_accumulate(flux_tmp, fluxes, G, wt2, forces)
     enddo ; enddo
   endif
 
+  ! Forcings introduced for MARBL
+  ! NOTE: fluxes%salt_flux, %sw, and %p_surf_full are handled above
+  if (associated(fluxes%nhx_dep) .and. associated(flux_tmp%nhx_dep)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%nhx_dep(i,j)  = wt1*fluxes%nhx_dep(i,j) + wt2*flux_tmp%nhx_dep(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%noy_dep) .and. associated(flux_tmp%noy_dep)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%noy_dep(i,j)  = wt1*fluxes%noy_dep(i,j) + wt2*flux_tmp%noy_dep(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%atm_co2) .and. associated(flux_tmp%atm_co2)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%atm_co2(i,j)  = wt1*fluxes%atm_co2(i,j) + wt2*flux_tmp%atm_co2(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%atm_alt_co2) .and. associated(flux_tmp%atm_alt_co2)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%atm_alt_co2(i,j)  = wt1*fluxes%atm_alt_co2(i,j) + wt2*flux_tmp%atm_alt_co2(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%dust_flux) .and. associated(flux_tmp%dust_flux)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%dust_flux(i,j)  = wt1*fluxes%dust_flux(i,j) + wt2*flux_tmp%dust_flux(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%iron_flux) .and. associated(flux_tmp%iron_flux)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%iron_flux(i,j)  = wt1*fluxes%iron_flux(i,j) + wt2*flux_tmp%iron_flux(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%atm_fine_dust_flux) .and. associated(flux_tmp%atm_fine_dust_flux)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%atm_fine_dust_flux(i,j)  = wt1*fluxes%atm_fine_dust_flux(i,j) + wt2*flux_tmp%atm_fine_dust_flux(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%atm_coarse_dust_flux) .and. associated(flux_tmp%atm_coarse_dust_flux)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%atm_coarse_dust_flux(i,j)  = wt1*fluxes%atm_coarse_dust_flux(i,j) + wt2*flux_tmp%atm_coarse_dust_flux(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%atm_bc_flux) .and. associated(flux_tmp%atm_bc_flux)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%atm_bc_flux(i,j)  = wt1*fluxes%atm_bc_flux(i,j) + wt2*flux_tmp%atm_bc_flux(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%seaice_dust_flux) .and. associated(flux_tmp%seaice_dust_flux)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%seaice_dust_flux(i,j)  = wt1*fluxes%seaice_dust_flux(i,j) + wt2*flux_tmp%seaice_dust_flux(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%seaice_bc_flux) .and. associated(flux_tmp%seaice_bc_flux)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%seaice_bc_flux(i,j)  = wt1*fluxes%seaice_bc_flux(i,j) + wt2*flux_tmp%seaice_bc_flux(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%fracr_cat) .and. associated(flux_tmp%fracr_cat)) then
+    do n=1,size(fluxes%fracr_cat,dim=3) ; do j=jsd,jed ; do i=isd,ied
+      fluxes%fracr_cat(i,j,n)  = wt1*fluxes%fracr_cat(i,j,n) + wt2*flux_tmp%fracr_cat(i,j,n)
+    enddo ; enddo ; enddo
+  endif
+  if (associated(fluxes%qsw_cat) .and. associated(flux_tmp%qsw_cat)) then
+    do n=1,size(fluxes%qsw_cat,dim=3) ; do j=jsd,jed ; do i=isd,ied
+      fluxes%qsw_cat(i,j,n)  = wt1*fluxes%qsw_cat(i,j,n) + wt2*flux_tmp%qsw_cat(i,j,n)
+    enddo ; enddo ; enddo
+  endif
+  if (associated(fluxes%ice_fraction) .and. associated(flux_tmp%ice_fraction)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%ice_fraction(i,j)  = wt1*fluxes%ice_fraction(i,j) + wt2*flux_tmp%ice_fraction(i,j)
+    enddo ; enddo
+  endif
+  if (associated(fluxes%u10_sqr) .and. associated(flux_tmp%u10_sqr)) then
+    do j=jsd,jed ; do i=isd,ied
+      fluxes%u10_sqr(i,j)  = wt1*fluxes%u10_sqr(i,j) + wt2*flux_tmp%u10_sqr(i,j)
+    enddo ; enddo
+  endif
+
   if (coupler_type_initialized(fluxes%tr_fluxes) .and. &
       coupler_type_initialized(flux_tmp%tr_fluxes)) &
     call coupler_type_increment_data(flux_tmp%tr_fluxes, fluxes%tr_fluxes, &
@@ -3366,6 +3481,21 @@ subroutine forcing_diagnostics(fluxes_in, sfc_state, G_in, US, time_end, diag, h
     if ((handles%id_u10_sqr > 0) .and. associated(fluxes%u10_sqr)) &
       call post_data(handles%id_u10_sqr, fluxes%u10_sqr, diag)
 
+    if ((handles%id_atm_fine_dust_flux > 0) .and. associated(fluxes%atm_fine_dust_flux)) &
+      call post_data(handles%id_atm_fine_dust_flux, fluxes%atm_fine_dust_flux, diag)
+
+    if ((handles%id_atm_coarse_dust_flux > 0) .and. associated(fluxes%atm_coarse_dust_flux)) &
+      call post_data(handles%id_atm_coarse_dust_flux, fluxes%atm_coarse_dust_flux, diag)
+
+    if ((handles%id_atm_bc_flux > 0) .and. associated(fluxes%atm_bc_flux)) &
+      call post_data(handles%id_atm_bc_flux, fluxes%atm_bc_flux, diag)
+
+    if ((handles%id_seaice_dust_flux > 0) .and. associated(fluxes%seaice_dust_flux)) &
+      call post_data(handles%id_seaice_dust_flux, fluxes%seaice_dust_flux, diag)
+
+    if ((handles%id_seaice_bc_flux > 0) .and. associated(fluxes%seaice_bc_flux)) &
+      call post_data(handles%id_seaice_bc_flux, fluxes%seaice_bc_flux, diag)
+
     ! remaining boundary terms ==================================================
 
     if ((handles%id_psurf > 0) .and. associated(fluxes%p_surf))                      &
@@ -3535,6 +3665,11 @@ subroutine allocate_forcing_by_group(G, fluxes, water, heat, ustar, press, &
   call myAlloc(fluxes%atm_alt_co2,isd,ied,jsd,jed, marbl)
   call myAlloc(fluxes%dust_flux,isd,ied,jsd,jed, marbl)
   call myAlloc(fluxes%iron_flux,isd,ied,jsd,jed, marbl)
+  call myAlloc(fluxes%atm_fine_dust_flux,isd,ied,jsd,jed, marbl)
+  call myAlloc(fluxes%atm_coarse_dust_flux,isd,ied,jsd,jed, marbl)
+  call myAlloc(fluxes%atm_bc_flux,isd,ied,jsd,jed, marbl)
+  call myAlloc(fluxes%seaice_dust_flux,isd,ied,jsd,jed, marbl)
+  call myAlloc(fluxes%seaice_bc_flux,isd,ied,jsd,jed, marbl)
 
   ! These fields should only be allocated when receiving multiple ice categories
   if (present(ice_ncat)) then
@@ -3841,6 +3976,11 @@ subroutine deallocate_forcing_type(fluxes)
   if (associated(fluxes%atm_alt_co2))          deallocate(fluxes%atm_alt_co2)
   if (associated(fluxes%dust_flux))            deallocate(fluxes%dust_flux)
   if (associated(fluxes%iron_flux))            deallocate(fluxes%iron_flux)
+  if (associated(fluxes%atm_fine_dust_flux))   deallocate(fluxes%atm_fine_dust_flux)
+  if (associated(fluxes%atm_coarse_dust_flux)) deallocate(fluxes%atm_coarse_dust_flux)
+  if (associated(fluxes%atm_bc_flux))          deallocate(fluxes%atm_bc_flux)
+  if (associated(fluxes%seaice_dust_flux))     deallocate(fluxes%seaice_dust_flux)
+  if (associated(fluxes%seaice_bc_flux))       deallocate(fluxes%seaice_bc_flux)
   if (associated(fluxes%fracr_cat))            deallocate(fluxes%fracr_cat)
   if (associated(fluxes%qsw_cat))              deallocate(fluxes%qsw_cat)