landfast ice parameterisation with lateral drag (MITgcm#741)

* cleaned up version of Yuqing Liu's lateral drag parameterisation

- Liu et al 2021, JGR, A new parameterization of coastal drag to
simulate landfast ice in deep marginal seas in the Arctic

* update SEAICE_OPTIONS.h and compile SEAICE_ALLOW_SIDEDRAG code

* add new file required for new parameterization

* add diagnostics of lateral drag terms

* add minimal documenation for basal and lateral drag

by listing the runtime parameters and add references

* make lateral drag work with EVP, JFNK, and Krylov solvers

* fix a typo

* add minimal documentation for some new fields and parameters

* check for more erroneously set parameters

* make it work with TAF, more store directives and include code

-compile new code in offline_exf_seaice.ad without using it

* more consistent names and some extra explanations

* redefine the position of u/vCoastLineFile to be at U/V-points

- add some more explanation in SEAICE_PARAMS.h
- remove confusing averaging but keep normalisation

* add missing variable declaration

* remove unnecessary extra error message

* update changed variable names in taf-store directive

spotted by jm-c (thanks!)

* more consistent upper/lower case

* tweak some comments and description

* document adding landfast ice parameterisation

---------

Co-authored-by: Jean-Michel Campin <[email protected]>
Co-authored-by: Jean-Michel Campin <[email protected]>

doc/manual_references.bib

@@ -1298,6 +1298,19 @@ @Article{lemieux:12
   pages =	 {5926--5944},
   doi =		 {10.1016/j.jcp.2012.05.024}}
 
+@Article{lemieux:15,
+  author =	 {Jean-Fran{\c{c}}ois Lemieux and L. Bruno Tremblay,
+                  Fr{\'e}d{\'e}ric Dupont and Mathieu Plante and
+                  Gregory C. Smith and Dany Dumont},
+  title =	 {A Basal Stress Parameterization for Modeling
+                  Land-Fast Ice},
+  journal =	 jgr,
+  year =	 2015,
+  volume =	 120,
+  pages =	 {3157--3173},
+  doi =		 {10.1002/2014JC010678}
+}
+
 @Article{leppaeranta:83,
   author =	 {Lepp{\"a}ranta, M.},
   title =	 {A Growth Model for Black Ice, Snow Ican and Snow Thickness in Subarctic Basins},
@@ -1338,6 +1351,20 @@ @Article{lhans:74
  pages        = {118-133}
 }
 
+@Article{liu:22,
+  author =	 {Yuqing Liu and Martin Losch and Nils Hutter and
+                  Longjiang Mu},
+  title =	 {A New Parameterization of Coastal Drag to Simulate
+                  Landfast Ice in Deep Marginal Seas in the {A}rctic.},
+  journal =	 jgr,
+  year =	 2022,
+  volume =	 127,
+  number =	 6,
+  pages =	 {e2022JC018413},
+  url =		 {https://doi.org/10.1029/2022JC018413},
+  doi =		 {10.1029/2022JC018413}
+}
+
 @Article{lipscomb:01,
   author =   {Lipscomb, W. H.},
   title =  {Remapping the thickness distribution in sea ice

doc/phys_pkgs/seaice.rst

@@ -73,6 +73,7 @@ automatically undefines more recent features, see :filelink:`SEAICE_OPTIONS.h
    :varlink:`SEAICE_ZETA_SMOOTHREG`, #define, use differentiable regularization for viscosities
    :varlink:`SEAICE_DELTA_SMOOTHREG`, #undef, use differentiable regularization :math:`\Delta_{\mathrm{reg}}=\sqrt{\Delta^2+\Delta_{\min}}` instead of :math:`\max`-function for :math:`1/\Delta_{\mathrm{reg}}`
    :varlink:`SEAICE_ALLOW_BOTTOMDRAG`, #undef, enable grounding parameterization for improved fastice in shallow seas
+   :varlink:`SEAICE_ALLOW_SIDEDRAG`, #undef, enable lateral drag parameterization for improved fastice along coastlines and islands
    :varlink:`SEAICE_BGRID_DYNAMICS`, #undef, use sea ice dynamics code on legacy B-grid; most of the previous flags are not available with B-grid
    :varlink:`SEAICE_BICE_STRESS`, #undef, B-grid only for backward compatiblity: turn on ice-stress on ocean; defined by default if :varlink:`SEAICE_BGRID_DYNAMICS` is defined
    :varlink:`EXPLICIT_SSH_SLOPE`, #undef, B-grid only for backward compatiblity: use ETAN for tilt computations rather than geostrophic velocities; defined by default if :varlink:`SEAICE_BGRID_DYNAMICS` is defined
@@ -331,6 +332,20 @@ General flags and parameters
   +------------------------------------+------------------------------+-------------------------------------------------------------------------+
   | :varlink:`SEAICE_useMultDimSnow`   | TRUE                         | use same fixed pdf for snow as for multi-thickness-category ice         |
   +------------------------------------+------------------------------+-------------------------------------------------------------------------+
+  | :varlink:`SEAICEbasalDragK1`       | 8.0                          | basal drag parameter K\ :sub:`1` :cite:`lemieux:15`                     |
+  +------------------------------------+------------------------------+-------------------------------------------------------------------------+
+  | :varlink:`SEAICEbasalDragK2`       | 0.0                          | basal drag parameter K\ :sub:`2`                                        |
+  +------------------------------------+------------------------------+-------------------------------------------------------------------------+
+  | :varlink:`SEAICE_cBasalStar`       | :varlink:`SEAICE_cStar` value| basal drag parameter (no units)                                         |
+  +------------------------------------+------------------------------+-------------------------------------------------------------------------+
+  | :varlink:`SEAICEbasalDragU0`       | 5.E-5                        | basal drag parameter (m/s)                                              |
+  +------------------------------------+------------------------------+-------------------------------------------------------------------------+
+  | :varlink:`SEAICESideDrag`          | 0.0                          | lateral drag coefficient :cite:`liu:22`                                 |
+  +------------------------------------+------------------------------+-------------------------------------------------------------------------+
+  | :varlink:`uCoastLineFile`          | unset                        | filename for coastline length for u-equation                            |
+  +------------------------------------+------------------------------+-------------------------------------------------------------------------+
+  | :varlink:`vCoastLineFile`          | unset                        | filename for coastline length for v-equation                            |
+  +------------------------------------+------------------------------+-------------------------------------------------------------------------+
 
 
 The following dynamical ice thickness distribution and ridging parameters in

doc/tag-index

@@ -1,6 +1,8 @@
     Notes on tags used in MITgcmUV
     ==============================
 
+o pkg/seaice: 
+  - add new lateral drag parameterisation for landfast ice (Liu etal 2022, JGR)
 o verification/*/results:
   - update 5 adm + 3 tlm ref. output affected by TAF new version 6.4.5 & 6.5.0
 o pkg/ecco:

pkg/seaice/SEAICE.h

@@ -188,6 +188,20 @@ C     CbobC :: (linear) bottom drag coefficient for basals stress param.
       _RL CbotC      (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
 # endif /* SEAICE_ALLOW_BOTTOMDRAG */
 
+# ifdef SEAICE_ALLOW_SIDEDRAG
+C     coastRoughU/V :: coast line roughness (w/out units) in U and V direction,
+C                      computed from the coastline length at corner points,
+C                      interpolated to U/V points, and scaled by the grid cell
+C                      width
+C     sideDragU/V   :: drag coefficients for lateral drag a parameterisation
+      COMMON/SEAICE_SIDEDRAG/ sideDragU, sideDragV,
+     &     coastRoughU, coastRoughV
+      _RL sideDragU  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
+      _RL sideDragV  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
+      _RL coastRoughU(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
+      _RL coastRoughV(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
+# endif /* SEAICE_ALLOW_SIDEDRAG */
+
 # if ( defined SEAICE_ALLOW_JFNK ||  defined SEAICE_ALLOW_KRYLOV )
 C     diagnostics for the JFNK and Krylov solver
       INTEGER totalNewtonIters

pkg/seaice/SEAICE_OPTIONS.h

@@ -188,6 +188,10 @@ C     iteration depends on the residual. Good for when a non-linear
 C     convergence criterion must be satified
 # undef SEAICE_ALLOW_LSR_FLEX
 
+C     Use parameterisation of explicit lateral drag for a better
+C     representation of fastice along coast lines and islands
+# undef SEAICE_ALLOW_SIDEDRAG
+
 #endif /* SEAICE_CGRID */
 
 #ifdef SEAICE_BGRID_DYNAMICS

pkg/seaice/SEAICE_PARAMS.h

@@ -318,12 +318,26 @@ C
      &     SEAICE_debugPointJ
 
 C--   COMMON /SEAICE_PARM_C/ Character valued sea ice model parameters.
-C     AreaFile          :: File containing initial sea-ice concentration
-C     HsnowFile         :: File containing initial snow thickness
-C     HsaltFile         :: File containing initial sea ice salt content
-C     HeffFile          :: File containing initial sea-ice thickness
-C     uIceFile          :: File containing initial sea-ice U comp. velocity
-C     vIceFile          :: File containing initial sea-ice V comp. velocity
+C     AreaFile       :: File containing initial sea-ice concentration
+C     HsnowFile      :: File containing initial snow thickness
+C     HsaltFile      :: File containing initial sea ice salt content
+C     HeffFile       :: File containing initial sea-ice thickness
+C     uIceFile       :: File containing initial sea-ice U comp. velocity
+C     vIceFile       :: File containing initial sea-ice V comp. velocity
+C     uCoastLineFile :: File containing the some measure of coastline
+C                       roughness length (in m) at the U-points in the
+C                       X-direction (i.e. for the U-equation).
+C     vCoastLineFile :: Files containing the some measure of coastline
+C                       roughness length (in m) at the V-points in the
+C                       Y-direction (i.e. for the V-equation).
+C
+C                       This roughness length can be the subgrid
+C                       scale length of the coastline in a grid cell
+C                       projected in the direction normal to the u/v-
+C                       direction as in Liu et al. (2022), but it can
+C                       also be anything that is a good proxy of coast
+C                       line roughness.
+C
 C        !!! NOTE !!! Initial sea-ice thickness can also be set using
 C        SEAICE_initialHEFF below.  But a constant initial condition
 C        can mean large artificial fluxes of heat and freshwater in
@@ -335,9 +349,11 @@ C
       CHARACTER*(MAX_LEN_FNAM) HeffFile
       CHARACTER*(MAX_LEN_FNAM) uIceFile
       CHARACTER*(MAX_LEN_FNAM) vIceFile
+      CHARACTER*(MAX_LEN_FNAM) uCoastLineFile
+      CHARACTER*(MAX_LEN_FNAM) vCoastLineFile
       COMMON /SEAICE_PARM_C/
      &   AreaFile, HsnowFile, HsaltFile, HeffFile,
-     &   uIceFile, vIceFile
+     &   uIceFile, vIceFile, uCoastLineFile, vCoastLineFile
 
 C--   COMMON /SEAICE_PARM_RL/ Real valued parameters of sea ice model.
 C     SEAICE_deltaTtherm :: Seaice timestep for thermodynamic equations (s)
@@ -415,6 +431,11 @@ C     SEAICEbasalDragU0 (default = 5e-5)
 C     SEAICEbasalDragK1 (default = 8)
 C     SEAICEbasalDragK2  :: if > 0, turns on basal drag
 C                           (default = 0, Lemieux suggests 15)
+C     SEAICEsideDrag     :: if > 0, turns on lateral static drag
+C                           if < 0, turns on lateral quadratic drag
+C                           both are different landfast ice parameterisations
+C                           (Liu et al 2022 use 2e-4,
+C                            the default = 0 turns off the parameterisations)
 C
 C     SEAICE_wetAlbTemp  :: Temp (deg.C) above which wet-albedo values are used
 C     SEAICE_waterAlbedo :: water albedo
@@ -508,6 +529,7 @@ C
       _RL SEAICE_drySnowAlb_south, SEAICE_wetSnowAlb_south, HO_south
       _RL SEAICE_cBasalStar, SEAICEbasalDragU0
       _RL SEAICEbasalDragK1, SEAICEbasalDragK2
+      _RL SEAICEsideDrag
       _RL SEAICE_wetAlbTemp, SEAICE_waterAlbedo
       _RL SEAICE_strength, SEAICE_cStar, SEAICEpressReplFac
       _RL SEAICE_tensilFac, SEAICE_tensilDepth
@@ -566,6 +588,7 @@ C
      &    SEAICE_drySnowAlb_south, SEAICE_wetSnowAlb_south, HO_south,
      &    SEAICE_cBasalStar, SEAICEbasalDragU0,
      &    SEAICEbasalDragK1, SEAICEbasalDragK2,
+     &    SEAICEsideDrag,
      &    SEAICE_wetAlbTemp, SEAICE_waterAlbedo,
      &    SEAICE_strength, SEAICE_cStar, SEAICE_eccen, SEAICE_eccfr,
      &    SEAICEtdMU, SEAICEmcMu,

pkg/seaice/seaice_ad_check_lev1_dir.h

@@ -30,6 +30,10 @@ CADJ STORE dwatn       = comlev1, key=ikey_dynamics, kind=isbyte
 # ifdef SEAICE_ALLOW_BOTTOMDRAG
 CADJ STORE cbotc       = comlev1, key=ikey_dynamics, kind=isbyte
 # endif
+# ifdef SEAICE_ALLOW_SIDEDRAG
+CADJ STORE sideDragU   = comlev1, key=ikey_dynamics, kind=isbyte
+CADJ STORE sideDragV   = comlev1, key=ikey_dynamics, kind=isbyte
+# endif
 CADJ STORE seaicemasku = comlev1, key=ikey_dynamics, kind=isbyte
 CADJ STORE seaicemaskv = comlev1, key=ikey_dynamics, kind=isbyte
 # ifdef SEAICE_ALLOW_EVP

pkg/seaice/seaice_ad_check_lev2_dir.h

@@ -15,6 +15,10 @@ CADJ STORE VICE              = tapelev2, key = ilev_2
 CADJ STORE TICES             = tapelev2, key = ilev_2
 # ifdef SEAICE_CGRID
 CADJ STORE dwatn             = tapelev2, key = ilev_2
+# ifdef SEAICE_ALLOW_SIDEDRAG
+CADJ STORE sideDragU         = tapelev2, key = ilev_2
+CADJ STORE sideDragV         = tapelev2, key = ilev_2
+# endif
 CADJ STORE stressDivergenceX = tapelev2, key = ilev_2
 CADJ STORE stressDivergenceY = tapelev2, key = ilev_2
 #  ifdef SEAICE_ALLOW_EVP

pkg/seaice/seaice_ad_check_lev3_dir.h

@@ -15,6 +15,10 @@ CADJ STORE VICE              = tapelev3, key = ilev_3
 CADJ STORE TICES             = tapelev3, key = ilev_3
 # ifdef SEAICE_CGRID
 CADJ STORE dwatn             = tapelev3, key = ilev_3
+# ifdef SEAICE_ALLOW_SIDEDRAG
+CADJ STORE sideDragU         = tapelev3, key = ilev_3
+CADJ STORE sideDragV         = tapelev3, key = ilev_3
+# endif
 CADJ STORE stressDivergenceX = tapelev3, key = ilev_3
 CADJ STORE stressDivergenceY = tapelev3, key = ilev_3
 #  ifdef SEAICE_ALLOW_EVP

pkg/seaice/seaice_ad_check_lev4_dir.h

@@ -15,6 +15,10 @@ CADJ STORE VICE              = tapelev4, key = ilev_4
 CADJ STORE TICES             = tapelev4, key = ilev_4
 # ifdef SEAICE_CGRID
 CADJ STORE dwatn             = tapelev4, key = ilev_4
+# ifdef SEAICE_ALLOW_SIDEDRAG
+CADJ STORE sideDragU         = tapelev4, key = ilev_4
+CADJ STORE sideDragV         = tapelev4, key = ilev_4
+# endif
 CADJ STORE stressDivergenceX = tapelev4, key = ilev_4
 CADJ STORE stressDivergenceY = tapelev4, key = ilev_4
 #  ifdef SEAICE_ALLOW_EVP

pkg/seaice/seaice_ad_diff.list

@@ -35,3 +35,4 @@ seaice_calc_ice_strength.f
 seaice_reg_ridge.f
 seaice_bottomdrag_coeffs.f
 seaice_calc_stressdiv.f
+seaice_sidedrag_stress.f

pkg/seaice/seaice_calc_lhs.F

@@ -176,6 +176,14 @@ SUBROUTINE SEAICE_CALC_LHS(
      &         +  DWATN(I,J-1,bi,bj) * 0.5 _d 0 *
      &         (uIceLoc(I,J-1,bi,bj)+uIceLoc(I+1,J-1,bi,bj))
      &         ) ) * areaS(I,J)
+#ifdef SEAICE_ALLOW_SIDEDRAG
+C     lateral drag terms: + sideDragU*uIce
+          uIceLHS(I,J,bi,bj) = uIceLHS(I,J,bi,bj)
+     &         + sideDragU(I,J,bi,bj)*uIceLoc(I,J,bi,bj)
+C                         + sideDragV*vIce
+          vIceLHS(I,J,bi,bj) = vIceLHS(I,J,bi,bj)
+     &         + sideDragV(I,J,bi,bj)*vIceLoc(I,J,bi,bj)
+#endif /* SEAICE_ALLOW_SIDEDRAG */
 C     apply masks for interior (important when we have open boundaries)
           uIceLHS(I,J,bi,bj) = uIceLHS(I,J,bi,bj)*maskinW(I,J,bi,bj)
           vIceLHS(I,J,bi,bj) = vIceLHS(I,J,bi,bj)*maskinS(I,J,bi,bj)

pkg/seaice/seaice_calc_residual.F

@@ -25,9 +25,9 @@ SUBROUTINE SEAICE_CALC_RESIDUAL(
 #include "SIZE.h"
 #include "EEPARAMS.h"
 #include "SEAICE_SIZE.h"
-#ifdef SEAICE_ALLOW_BOTTOMDRAG
-#include "SEAICE_PARAMS.h"
-#endif /* SEAICE_ALLOW_BOTTOMDRAG */
+#if ( defined SEAICE_ALLOW_BOTTOMDRAG || defined SEAICE_ALLOW_SIDEDRAG )
+# include "SEAICE_PARAMS.h"
+#endif
 #include "SEAICE.h"
 
 C     !INPUT/OUTPUT PARAMETERS:
@@ -91,6 +91,12 @@ SUBROUTINE SEAICE_CALC_RESIDUAL(
      O     CbotC,
      I     krylovIter, myTime, myIter, myThid )
 #endif /* SEAICE_ALLOW_BOTTOMDRAG */
+#ifdef SEAICE_ALLOW_SIDEDRAG
+        IF ( SEAICEsideDrag .NE. 0. _d 0 ) CALL SEAICE_SIDEDRAG_STRESS(
+     I       uIceLoc, vIceLoc, coastRoughU, coastRoughV, AREA,
+     O       sideDragU, sideDragV,
+     I       krylovIter, myTime, myIter, myThid )
+#endif /* SEAICE_ALLOW_SIDEDRAG */
       CALL SEAICE_CALC_STRAINRATES(
      I     uIceLoc, vIceLoc,
      O     e11, e22, e12,

pkg/seaice/seaice_check.F

@@ -1169,7 +1169,7 @@ SUBROUTINE SEAICE_CHECK( myThid )
 
       IF ( SEAICEpressReplFac .LT. 0. _d 0 .OR.
      &     SEAICEpressReplFac .GT. 1. _d 0 ) THEN
-       WRITE(msgBuf,'(A,I2)')
+       WRITE(msgBuf,'(A,F5.2)')
      &      'SEAICE_CHECK: SEAICEpressReplFac = ', SEAICEpressReplFac
        CALL PRINT_ERROR( msgBuf, myThid )
        WRITE(msgBuf,'(A)') 'SEAICE_CHECK: cannot be < 0 or > 1'
@@ -1212,7 +1212,7 @@ SUBROUTINE SEAICE_CHECK( myThid )
 
 #ifndef SEAICE_ALLOW_BOTTOMDRAG
       IF ( SEAICEbasalDragK2 .GT. 0. _d 0 ) THEN
-       WRITE(msgBuf,'(A,I2)')
+       WRITE(msgBuf,'(A,E12.5)')
      &      'SEAICE_CHECK: SEAICEbasalDragK2 = ', SEAICEbasalDragK2
        CALL PRINT_ERROR( msgBuf, myThid )
        WRITE(msgBuf,'(2A)') 'SEAICE_CHECK: is greater than 0, ',
@@ -1222,6 +1222,32 @@ SUBROUTINE SEAICE_CHECK( myThid )
       ENDIF
 #endif /* SEAICE_ALLOW_BOTTOMDRAG */
 
+#ifndef SEAICE_ALLOW_SIDEDRAG
+      IF ( SEAICEsideDrag .NE. 0. _d 0 ) THEN
+       WRITE(msgBuf,'(A,E12.5,A)')
+     &      'SEAICE_CHECK: SEAICEsideDrag = ', SEAICEsideDrag,
+     &      ' is set, but SEAICE_ALLOW_SIDEDRAG is not defined'
+       CALL PRINT_ERROR( msgBuf, myThid )
+       errCount = errCount + 1
+      ENDIF
+      IF ( uCoastLineFile .NE. ' ' ) THEN
+       i = ILNBLNK(uCoastLineFile)
+       WRITE(msgBuf,'(3A)')
+     &      'SEAICE_CHECK: uCoastLineFile = ', uCoastLineFile(1:i),
+     &      ' is set, but SEAICE_ALLOW_SIDEDRAG is not defined'
+       CALL PRINT_ERROR( msgBuf, myThid )
+       errCount = errCount + 1
+      ENDIF
+      IF( vCoastLineFile .NE. ' ' ) THEN
+       i = ILNBLNK(vCoastLineFile)
+       WRITE(msgBuf,'(3A)')
+     &      'SEAICE_CHECK: vCoastLineFile = ', vCoastLineFile(1:i),
+     &      ' is set, but SEAICE_ALLOW_SIDEDRAG is not defined'
+       CALL PRINT_ERROR( msgBuf, myThid )
+       errCount = errCount + 1
+      ENDIF
+#endif
+
 #ifdef SEAICE_ITD
 C     The ice thickness distribution (ITD) module can only be used with
 C     the zero-layer thermodynamics of S/R SEAICE_GROWTH and the

pkg/seaice/seaice_diagnostics_init.F

@@ -343,6 +343,24 @@ SUBROUTINE SEAICE_DIAGNOSTICS_INIT( myThid )
 c     CALL DIAGNOSTICS_ADDTOLIST( diagNum,
 c    I     diagName, diagCode, diagUnits, diagTitle, diagMate, myThid )
 
+# ifdef SEAICE_ALLOW_SIDEDRAG
+      diagName  = 'SIlatDgU'
+      diagTitle = 'SEAICE lateral coastal drag on U momentum'
+      diagUnits = 'N/m^2           '
+      diagCode  = 'UU      M1      '
+      diagMate  = diagNum + 2
+      CALL DIAGNOSTICS_ADDTOLIST( diagNum,
+     I     diagName, diagCode, diagUnits, diagTitle, 0, myThid )
+
+      diagName  = 'SIlatDgV'
+      diagTitle = 'SEAICE lateral coastal drag on V momentum'
+      diagUnits = 'N/m^2           '
+      diagCode  = 'VV      M1      '
+      diagMate  = diagNum
+      CALL DIAGNOSTICS_ADDTOLIST( diagNum,
+     I     diagName, diagCode, diagUnits, diagTitle, 0, myThid )
+#endif
+
 C     SIqnet, Qnet, and QNETtave are identical.
 C     With #undef NONLIN_FRSURF SIqnet is identical to -(TFLUX-TRELAX).
 C     Except over land and under sea ice, SIqnet is also identical to

pkg/seaice/seaice_dynsolver.F

@@ -415,6 +415,36 @@ SUBROUTINE SEAICE_DYNSOLVER( myTime, myIter, myThid )
        CALL DIAGNOSTICS_FILL(eta    ,'SIeta   ',0,1,0,1,1,myThid)
        CALL DIAGNOSTICS_FILL(press  ,'SIpress ',0,1,0,1,1,myThid)
        CALL DIAGNOSTICS_FILL(deltaC ,'SIdelta ',0,1,0,1,1,myThid)
+# ifdef SEAICE_ALLOW_SIDEDRAG
+C     recompute lateral coast drag terms
+       IF ( DIAGNOSTICS_IS_ON('SIlatDgU',myThid) ) THEN
+        DO bj = myByLo(myThid), myByHi(myThid)
+         DO bi = myBxLo(myThid), myBxHi(myThid)
+C     use sig1 as a temporary field
+          DO j=1,sNy
+           DO i=1,sNx+1
+            sig1(i,j) = sideDragU(i,j,bi,bj)*uIce(i,j,bi,bj)
+           ENDDO
+          ENDDO
+          CALL DIAGNOSTICS_FILL(sig1,'SIlatDgU',0,1,2,bi,bj,myThid)
+         ENDDO
+        ENDDO
+       ENDIF
+       IF ( DIAGNOSTICS_IS_ON('SIlatDgV',myThid) ) THEN
+        DO bj = myByLo(myThid), myByHi(myThid)
+         DO bi = myBxLo(myThid), myBxHi(myThid)
+C     use sig1 as a temporary field
+          DO j=1,sNy+1
+           DO i=1,sNx
+            sig1(i,j) = sideDragV(i,j,bi,bj)*vIce(i,j,bi,bj)
+           ENDDO
+          ENDDO
+          CALL DIAGNOSTICS_FILL(sig1,'SIlatDgV',0,1,2,bi,bj,myThid)
+         ENDDO
+        ENDDO
+       ENDIF
+# endif /* SEAICE_ALLOW_SIDEDRAG */
+
        IF ( DIAGNOSTICS_IS_ON('SItensil',myThid) ) THEN
         DO bj = myByLo(myThid), myByHi(myThid)
          DO bi = myBxLo(myThid), myBxHi(myThid)