Inflexible calendar year and approximate conversion in option 2 of the freshwater-budget adjustment mechanism

[smueller]

Context

Option 2 of the freshwater-budget adjustment mechanism (nn_fwb = 2, module sbcfwb) uses a fixed calendar year of 365 days rather than that of the selected model calendar to time annual updates of the adjustment flux; this issue has previously been reported in the form of source-code comments (source:/NEMO/releases/r4.0/r4.0-HEAD/src/OCE/SBC/sbcfwb.F90@12578#L129​ and source:/NEMO/releases/r4.0/r4.0-HEAD/src/OCE/SBC/sbcfwb.F90@12578#L135​). Further, the conversion of the combined sea-level, sea-ice, and snow anomaly into an equivalent freshwater flux is approximate and not in line with corresponding conversions used elsewhere in the model.

Analysis

In option 2 of the freshwater-budget adjustment mechanism, year boundaries are identified whenever the reminder of the division of the time step counter by 365 * 86400 / rdt equals zero (source:/NEMO/releases/r4.0/r4.0-HEAD/src/OCE/SBC/sbcfwb.F90@12578#L130​). This i) fixes the adjustment interval to 365 days, ii) sets the date for annual adjustment-flux updates to the date of nit000 - 1 rather than the start of 1 January, and iii) leads to an update of the adjustment flux near the start of the final time step of each year (i.e., before completion of that model year).

Line source:/NEMO/releases/r4.0/r4.0-HEAD/src/OCE/SBC/sbcfwb.F90@12578#L134​ implements the conversion of a combined sea-level, sea ice, and snow anomaly into an equivalent annual-average freshwater flux using a literal value of 1000.0 for the density of sea water and 365 for the number of days in a year rather than variables rau0 and nyear_len(1), respectively.

The impact of these approximations can be demonstrated using reference configuration GYRE_PISCES (which uses calendar years of 360 days) after replacing source:/NEMO/releases/r4.0/r4.0-HEAD/src/OCE/USR/usrdef_sbc.F90@12578#L134​ by zsumemp = 0.0_wp and activation of the freshwater-budget adjustment by setting namelist parameter nn_fwb = 2.

Recommendation

The annual average freshwater adjustment based on the sea-level and sea-ice anomaly at the end of the preceding model year should be computed using daymod variable nyear_len (i.e., the calendar year during which the adjustment will be applied) and using the reference sea-water density rau0 rather than the literal values 365 and 1000.0, respectively. Further, it would be useful for interpretability of the effect of the freshwater-adjustment mechanism to exactly align the annual computation of the freshwater-budget adjustment with the model calendar. Therefore, adjustment-flux updates could occur in the first time step of each new year, which should be identified using variables updated by module daymod such as through the conditional ( ( nday_year == 1 ) .AND. ( nsec_day == NINT( 0.5_wp * rdt ) ) ).

[smueller]

In order to implement the recommendation above, source:/NEMO/releases/r4.0/r4.0-HEAD/src/OCE/SBC/sbcfwb.F90@13098​ could be modified according to:

src/OCE/SBC/sbcfwb.F90

@@ -124 +124 @@
             IF(lwp)WRITE(numout,*)'sbc_fwb : year = ',iyear  , ' freshwater budget correction = ', fwfold
             IF(lwp)WRITE(numout,*)'          year = ',iyear-1, ' freshwater budget read       = ', a_fwb
             IF(lwp)WRITE(numout,*)'          year = ',iyear-2, ' freshwater budget read       = ', a_fwb_b
-         ENDIF   
-         !                                         ! Update fwfold if new year start
-         ikty = 365 * 86400 / rdt                  !!bug  use of 365 days leap year or 360d year !!!!!!!
-         IF( MOD( kt, ikty ) == 0 ) THEN
+         ELSEIF( ( nday_year == 1 ) .AND. ( nsec_day == NINT( 0.5_wp * rdt ) ) ) THEN   ! Update fwfold at the start of 1 January
             a_fwb_b = a_fwb                           ! mean sea level taking into account the ice+snow
                                                       ! sum over the global domain
             a_fwb   = glob_sum( 'sbcfwb', e1e2t(:,:) * ( sshn(:,:) + snwice_mass(:,:) * r1_rau0 ) )
-            a_fwb   = a_fwb * 1.e+3 / ( area * rday * 365. )     ! convert in Kg/m3/s = mm/s
-!!gm        !                                                      !!bug 365d year 
+            a_fwb   = a_fwb * rau0 / ( area * rday * nyear_len(1) )     ! convert into kg/m3/s
             fwfold =  a_fwb                           ! current year freshwater budget correction
             !                                         ! estimate from the previous year budget
          ENDIF

[smueller]

In a discussion with Gurvan and Seb, Gurvan pointed out two further bugs related to this defect:

i) in contrast to model configurations with embedded sea-ice, in the presence of levitating sea-ice/snow (ln_ice_embd = .FALSE., default case) the initial sea-level is not adjusted to compensate for the initial combined sea-ice and snow mass (missing adjustment of initial sshn and sshb values by the global average of the initial snwice_mass field converted to an equivalent sea-level anomaly during initialisation of both SI3 and the CICE coupling interface), which leads to a non-zero global-average sea-level anomaly for an ice/snow-less ocean with identical total mass of water, and, when nn_fwb = 2 to an extra freshwater adjustment fluxes early in the simulation; and

ii) the freshwater-budget adjustment computation does not take into account under-ice-shelf areas (if any) as the computation of the surface area masks the grid-cell areas by tmask(:,:,1), instead of by ssmask(:,:) only (as applied in function glob_sum).

These two bugs could be fixed by making the modifications

src/ICE/iceistate.F90

@@ -91 +91 @@
       INTEGER, INTENT(in) ::   kt   ! time step 
       !!
       INTEGER  ::   ji, jj, jk, jl         ! dummy loop indices
-      REAL(wp) ::   ztmelts
+      REAL(wp) ::   ztmelts, z1_area, zsshadj
       INTEGER , DIMENSION(4)           ::   itest
       REAL(wp), DIMENSION(jpi,jpj)     ::   z2d
       REAL(wp), DIMENSION(jpi,jpj)     ::   zswitch    ! ice indicator
@@ -424 +424 @@
                gde3w_n(:,:,jk) = gdept_n(:,:,jk  ) - sshn (:,:)
             END DO
          ENDIF
+      ELSE
+         z1_area = 1.0_wp / glob_sum( 'iceistate', e1e2t(:,:) )
+         zsshadj = glob_sum( 'iceistate', e1e2t(:,:) * snwice_mass(:,:) ) * r1_rau0 * z1_area
+         WRITE(ctmp1,'(A34,F10.6,A32)') 'iceistate:   mean ssh adjusted by ', -1.0_wp * zsshadj, ' m to compensate for the initial'
+         CALL ctl_warn( ctmp1,          '             ice+snow mass' )
+         sshn(:,:) = sshn(:,:) - zsshadj
+         sshb(:,:) = sshb(:,:) - zsshadj
       ENDIF
       
       !------------------------------------

src/OCE/SBC/sbcice_cice.F90

@@ -154 +154 @@
       !!---------------------------------------------------------------------
       INTEGER, INTENT( in  ) ::   ksbc                ! surface forcing type
       REAL(wp), DIMENSION(jpi,jpj) :: ztmp1, ztmp2
-      REAL(wp) ::   zcoefu, zcoefv, zcoeff            ! local scalar
+      REAL(wp) ::   zcoefu, zcoefv, zcoeff            ! local scalars
+      REAL(wp) ::   z1_area, zsshadj                  !   "      "
       INTEGER  ::   ji, jj, jl, jk                    ! dummy loop indices
       !!---------------------------------------------------------------------
       !
@@ -265 +266 @@
                   gde3w_n(:,:,jk) = gdept_n(:,:,jk  ) - sshn   (:,:)
                END DO
             ENDIF
+         ELSE
+            z1_area = 1.0_wp / glob_sum( 'sbcice_cice', e1e2t(:,:) )
+            zsshadj = glob_sum( 'sbcice_cice', e1e2t(:,:) * snwice_mass(:,:) ) * r1_rau0 * z1_area
+            WRITE(ctmp1,'(A38,F10.6,A24)') 'sbcicce_cice:   mean ssh adjusted by ', -1.0_wp * zsshadj, ' m to compensate for the'
+            CALL ctl_warn( ctmp1,          '                initial ice+snow mass' )
+            sshn(:,:) = sshn(:,:) - zsshadj
+            sshb(:,:) = sshb(:,:) - zsshadj
          ENDIF
       ENDIF
       !

and

src/OCE/SBC/sbcfwb.F90

@@ -37 +37 @@
    REAL(wp) ::   a_fwb     ! for 2 year before (_b) and before year.
    REAL(wp) ::   fwfold    ! fwfold to be suppressed
    REAL(wp) ::   area      ! global mean ocean surface (interior domain)
+   REAL(wp) ::   area2     ! global ocean surface area (interior domain) incl. ice-shelf cavities
 
    !! * Substitutions
 #  include "vectopt_loop_substitute.h90"
@@ -90 +91 @@
          area = glob_sum( 'sbcfwb', e1e2t(:,:) * tmask(:,:,1))           ! interior global domain surface
          ! isf cavities are excluded because it can feedback to the melting with generation of inhibition of plumes
          ! and in case of no melt, it can generate HSSW.
+         area2 = glob_sum( 'sbcfwb', e1e2t(:,:) )                        ! global interior-domain ocean surface incl. ice-shelf cavities
          !
 #if ! defined key_si3 && ! defined key_cice
          snwice_mass_b(:,:) = 0.e0               ! no sea-ice model is being used : no snow+ice mass
@@ -131 +133 @@
             a_fwb_b = a_fwb                           ! mean sea level taking into account the ice+snow
                                                       ! sum over the global domain
             a_fwb   = glob_sum( 'sbcfwb', e1e2t(:,:) * ( sshn(:,:) + snwice_mass(:,:) * r1_rau0 ) )
-            a_fwb   = a_fwb * 1.e+3 / ( area * rday * 365. )     ! convert in Kg/m3/s = mm/s
+            a_fwb   = a_fwb * 1.e+3 / ( area2 * rday * 365. )     ! convert in Kg/m3/s = mm/s
 !!gm        !                                                      !!bug 365d year 
             fwfold =  a_fwb                           ! current year freshwater budget correction
             !                                         ! estimate from the previous year budget

[smueller]

In 13171:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

Gurvan pointed out another bug related to those described in comment:2: in the presence of embedded sea ice (ln_ice_embd = .TRUE.), the adjustment of the vertical scale factors e3t_n and e3t_b following the sea-level adjustment due to the initial sea-ice/snow mass in both SI3 and the CICE coupling interface ignores under-ice-shelf areas (if any). He recommends to make the modifications

src/ICE/iceistate.F90

@@ -387 +387 @@
          !
          IF( .NOT.ln_linssh ) THEN
             !
-            WHERE( ht_0(:,:) > 0 )   ;   z2d(:,:) = 1._wp + sshn(:,:)*tmask(:,:,1) / ht_0(:,:)
-            ELSEWHERE                ;   z2d(:,:) = 1._wp   ;   END WHERE
+            z2d(:,:) = 1._wp + sshn(:,:) * ssmask(:,:) / ( ht_0(:,:) + 1._wp - ssmask(:,:) )
             !
             DO jk = 1,jpkm1                     ! adjust initial vertical scale factors                
                e3t_n(:,:,jk) = e3t_0(:,:,jk) * z2d(:,:)

and

src/OCE/SBC/sbcice_cice.F90

@@ -234 +234 @@
             IF( .NOT.ln_linssh ) THEN
                !
                DO jk = 1,jpkm1                     ! adjust initial vertical scale factors
-                  e3t_n(:,:,jk) = e3t_0(:,:,jk)*( 1._wp + sshn(:,:)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) )
-                  e3t_b(:,:,jk) = e3t_0(:,:,jk)*( 1._wp + sshb(:,:)*tmask(:,:,1)/(ht_0(:,:) + 1.0 - tmask(:,:,1)) )
+                  e3t_n(:,:,jk) = e3t_0(:,:,jk)*( 1._wp + sshn(:,:)*ssmask(:,:)/(ht_0(:,:) + 1._wp - ssmask(:,:)) )
+                  e3t_b(:,:,jk) = e3t_0(:,:,jk)*( 1._wp + sshb(:,:)*ssmask(:,:)/(ht_0(:,:) + 1._wp - ssmask(:,:)) )
                ENDDO
                e3t_a(:,:,:) = e3t_b(:,:,:)
                ! Reconstruction of all vertical scale factors at now and before time-steps

[smueller]

In 13177:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

In 13179:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

In 13180:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

Using source:/NEMO/branches/2020/ticket2487@13180, runs of the model in the AGRIF_DEMO reference configuration terminate with an error during time step 59 of the nested domain 2.

A possible reason for this error could be the substantial differences between the initial sea-level adjustments made in the different domains (about -0.125 m in the global domain, 0.0 m in nested domain 1, about -0.881 m in nested domain 2, and about -1.491 in nested domain 3). For consistency, the initial sea-levels of all domains should be adjusted by the same amount.

[smueller]

In 13294:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

source:/NEMO/branches/2020/ticket2487@13294 passes the AGRIF_DEMO SETTE tests.

In order to reduce the number of reference configurations that are affected by the addition of the sea-level adjustment due to initial sea-ice and snow mass, the sea-level adjustment could be made conditional on a new namelist parameter and be activated in reference configurations with activated freshwater adjustment only (nn_fwb > 0; currently ORCA2_ICE_PISCES).

[smueller]

In 13299:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

In 13300:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

In 13301:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

In 15264:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

In 15410:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

In 15414:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

In 15418:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[smueller]

source:/NEMO/branches/2020/ticket2487@15418 has passed the standard SETTE tests of reference configuration ORCA2_ICE_PISCES.

[clem]

In 15438:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[clem]

In 15439:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found

[clem]

In 15451:

Error: Failed to load processor CommitTicketReference

No macro or processor named 'CommitTicketReference' found