Updates to Charnock parameter not set in ST6 (ticket NOAA-EMC#1307)

Previous commit 3780f17 Changes ------- [x] Add capped Charnock coefficient. SIN6 namelist parameters `CHKCAP`, `CHKINF`, `CHKSIG`. Cap is active for CHKCAP > 0.1 [x] Add regression test ww3_ts1 -i input_capcha
stefanzieger · Nov 25, 2024 · dbca9aa · dbca9aa
1 parent 65e3151
commit dbca9aa
Show file tree

Hide file tree

Showing 7 changed files with 86 additions and 32 deletions.
diff --git a/model/bin/comp.tmpl b/model/bin/comp.tmpl
@@ -73,7 +73,7 @@
   # oasis coupler include dir
   if [ "$oasis_mod" = 'yes' ]
   then
-    opt="$opt -I$OASISDIR/build/lib/psmile.MPI1"
+    opt="$opt -I$OASISDIR/include -I$OASISDIR/build/lib/psmile.MPI1"
   fi
 
   # netcdf include dir

diff --git a/model/src/w3gdatmd.F90 b/model/src/w3gdatmd.F90
@@ -915,9 +915,9 @@ MODULE W3GDATMD
     !
 #ifdef W3_ST6
     REAL                  :: SIN6A0, SDS6A1, SDS6A2, SWL6B1, &
-         SIN6WS, SIN6FC, SIN6AHAT
+         SIN6WS, SIN6FC, SIN6CHKMIN, SIN6CHKINF, SIN6CHKCAP, SIN6CHKSIG
     INTEGER               :: SDS6P1, SDS6P2
-    LOGICAL               :: SDS6ET, SWL6S6, SWL6CSTB1
+    LOGICAL               :: SDS6ET, SWL6S6, SWL6CSTB1, SIN6FLCAP
 #endif
   END TYPE SRCP
   !
@@ -1333,9 +1333,9 @@ MODULE W3GDATMD
 #endif
 #ifdef W3_ST6
   REAL, POINTER           :: SIN6A0, SDS6A1, SDS6A2, SWL6B1, &
-       SIN6WS, SIN6FC, SIN6AHAT
+       SIN6WS, SIN6FC, SIN6CHKMIN, SIN6CHKINF, SIN6CHKCAP, SIN6CHKSIG
   INTEGER, POINTER        :: SDS6P1, SDS6P2
-  LOGICAL, POINTER        :: SDS6ET, SWL6S6, SWL6CSTB1
+  LOGICAL, POINTER        :: SDS6ET, SWL6S6, SWL6CSTB1, SIN6FLCAP
 #endif
   REAL, POINTER           :: WWNMEANPTAIL, SSTXFTFTAIL
   !/
@@ -2710,7 +2710,11 @@ SUBROUTINE W3SETG ( IMOD, NDSE, NDST )
     SIN6A0 => MPARS(IMOD)%SRCPS%SIN6A0
     SIN6WS => MPARS(IMOD)%SRCPS%SIN6WS
     SIN6FC => MPARS(IMOD)%SRCPS%SIN6FC
-    SIN6AHAT => MPARS(IMOD)%SRCPS%SIN6AHAT
+    SIN6CHKMIN => MPARS(IMOD)%SRCPS%SIN6CHKMIN
+    SIN6CHKINF => MPARS(IMOD)%SRCPS%SIN6CHKINF
+    SIN6CHKCAP => MPARS(IMOD)%SRCPS%SIN6CHKCAP
+    SIN6CHKSIG => MPARS(IMOD)%SRCPS%SIN6CHKSIG
+    SIN6FLCAP => MPARS(IMOD)%SRCPS%SIN6FLCAP
     SDS6ET => MPARS(IMOD)%SRCPS%SDS6ET
     SDS6A1 => MPARS(IMOD)%SRCPS%SDS6A1
     SDS6P1 => MPARS(IMOD)%SRCPS%SDS6P1

diff --git a/model/src/w3gridmd.F90 b/model/src/w3gridmd.F90
@@ -869,8 +869,8 @@ MODULE W3GRIDMD
 #endif
   !
 #ifdef W3_ST6
-  REAL                    :: SINA0, SINWS, SINFC, SINAHAT,        &
-       SDSA1, SDSA2, SWLB1
+  REAL                    :: SINA0, SINWS, SINFC, CHKMIN,        &
+       CHKINF, CHKCAP, CHKSIG, SDSA1, SDSA2, SWLB1
   INTEGER                 :: SDSP1, SDSP2
   LOGICAL                 :: SDSET, CSTB1
 #endif
@@ -1056,7 +1056,7 @@ MODULE W3GRIDMD
 #endif
 
 #ifdef W3_ST6
-  NAMELIST /SIN6/ SINA0, SINWS, SINFC, SINAHAT
+  NAMELIST /SIN6/ SINA0, SINWS, SINFC, CHKMIN, CHKINF, CHKCAP, CHKSIG
   NAMELIST /SDS6/ SDSET, SDSA1, SDSA2, SDSP1, SDSP2
   NAMELIST /SWL6/ SWLB1, CSTB1
 #endif
@@ -1747,7 +1747,10 @@ SUBROUTINE W3GRID()
     SINA0  = 0.09
     SINWS  = 32.0
     SINFC  = 6.0
-    SINAHAT  = 0.0095
+    CHKMIN  = 0.0095
+    CHKINF  = 1.000E-4
+    CHKCAP  = 0.0000
+    CHKSIG  = 5.0000
 #endif
     !
 #ifdef W3_ST1
@@ -1843,10 +1846,20 @@ SUBROUTINE W3GRID()
     SIN6A0 = SINA0
     SIN6WS = SINWS
     SIN6FC = SINFC
-    SIN6AHAT = SINAHAT
+    SIN6CHKMIN = CHKMIN
+    SIN6CHKINF = CHKINF
+    SIN6CHKCAP = CHKCAP
+    SIN6CHKSIG = CHKSIG
+    SIN6FLCAP = .TRUE.
     J = 1
-    IF ( SIN6A0.LE.0. ) J = 2
-    WRITE (NDSO,921) YESXNO(J), SIN6A0, SIN6WS, SIN6FC, SIN6AHAT
+    JJ = 1
+    IF ( SIN6A0.LE.0.0 ) J = 2
+    IF ( SIN6CHKCAP.LT.1.0E-1 ) THEN
+        JJ = 2
+        SIN6FLCAP = .FALSE.
+    END IF
+    WRITE (NDSO,921) YESXNO(J), SIN6A0, SIN6WS, SIN6FC, SIN6CHKMIN, YESXNO(JJ)
+    IF ( JJ.EQ.1 ) WRITE (NDSO,9210) SIN6CHKINF, SIN6CHKCAP, SIN6CHKSIG
 #endif
     !
     ! 6.e Define Snl.
@@ -3251,7 +3264,7 @@ SUBROUTINE W3GRID()
            CAPCHA, CHAMIN, CHA0, UCAP, SIGMAUCAP
 #endif
 #ifdef W3_ST6
-      WRITE (NDSO,2920) SINA0, SINWS, SINFC, SINAHAT
+      WRITE (NDSO,2920) SINA0, SINWS, SINFC, CHKMIN, CHKINF, CHKCAP, CHKSIG
 #endif
 #ifdef W3_NL1
       WRITE (NDSO,2922) LAMBDA, NLPROP, KDCONV, KDMIN,       &
@@ -6313,9 +6326,17 @@ SUBROUTINE W3GRID()
          '  attenuation factor               :  ',F6.2/ &
          '  wind speed scaling factor        :  ',F6.2/ &
          '  frequency cut-off factor         :  ',F6.2/ &
-         '  minimum Charnock coefficient     :  ',F8.5/ )
+         '  Coupling parameters', /                     &
+         '    minimum Charnock coeficient    :  ',F7.5/ &
+         '    Charnock is capped             :  ',A)
+9210 FORMAT ('      asymptote value              :  ',F7.5/ &
+         '      wind speed threshold         :  ',F5.1/ &
+         '      transition window            :  ',F5.1/ )
+
 2920 FORMAT ( '  &SIN6 SINA0 =', F6.3, ', SINWS =', F6.2, &
-         ', SINFC =', F6.2, ', SINAHAT =', F8.5, ' /')
+         ', SINFC =', F6.2, ', CHKMIN =', F8.5, &
+         ', CHKINF =', F8.5, ', CHKCAP =', F5.1, &
+         ', CHKSIG =', F5.1, ' /')
 #endif
     !
 #ifdef W3_NL0

diff --git a/model/src/w3iogrmd.F90 b/model/src/w3iogrmd.F90
@@ -1573,20 +1573,24 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT &
     IF ( WRITE ) THEN
       WRITE (NDSM) SIN6A0, SDS6ET, SDS6A1, SDS6A2,           &
            SDS6P1, SDS6P2, SWL6S6, SWL6B1, SWL6CSTB1,        &
-           SIN6WS, SIN6FC, SIN6AHAT
+           SIN6WS, SIN6FC, SIN6CHKMIN, SIN6CHKINF,           &
+           SIN6CHKCAP, SIN6CHKSIG, SIN6FLCAP
 #ifdef W3_ASCII
       WRITE (NDSA,*) 'SIN6A0, SDS6ET, SDS6A1, SDS6A2,        &
            SDS6P1, SDS6P2, SWL6S6, SWL6B1, SWL6CSTB1,        &
-           SIN6WS, SIN6FC, AIN6AHAT:',                       &
+           SIN6WS, SIN6FC, SIN6CHKMIN, SIN6CHKINF,           &
+           SIN6CHKCAP, SIN6CHKSIG, SIN6FLCAP:',              &
                    SIN6A0, SDS6ET, SDS6A1, SDS6A2,           &
            SDS6P1, SDS6P2, SWL6S6, SWL6B1, SWL6CSTB1,        &
-           SIN6WS, SIN6FC, SIN6AHAT
+           SIN6WS, SIN6FC, SIN6CHKMIN, SIN6CHKINF,
+           SIN6CHKCAP, SIN6CHKSIG, SIN6FLCAP
 #endif
     ELSE
       READ (NDSM,END=801,ERR=802,IOSTAT=IERR)                &
            SIN6A0, SDS6ET, SDS6A1, SDS6A2,                   &
            SDS6P1, SDS6P2, SWL6S6, SWL6B1, SWL6CSTB1,        &
-           SIN6WS, SIN6FC, SIN6AHAT
+           SIN6WS, SIN6FC, SIN6CHKMIN, SIN6CHKINF,           &
+           SIN6CHKCAP, SIN6CHKSIG, SIN6FLCAP
     END IF
 #endif
     !

diff --git a/model/src/w3src6md.F90 b/model/src/w3src6md.F90
@@ -4,7 +4,7 @@
 !>
 !> @author S. Zieger
 !> @author Q. Liu
-!> @date   26-Jun-2018
+!> @date   11-Oct-2024
 !>
 
 #include "w3macros.h"
@@ -39,7 +39,7 @@ MODULE W3SRC6MD
   !/                                                         (S. Zieger)
   !/    26-Jun-2017 : Recalibration of ST6                ( verison 6.06 )
   !/                                                         (Q. Liu   )
-  !/    11-Oct-2024 : Charnock parameter output added     ( verison 7.14 )
+  !/    11-Oct-2024 : Add Charnock parameter              ( verison 7.14 )
   !/
   !/    Copyright 2009 National Weather Service (NWS),
   !/       National Oceanic and Atmospheric Administration.  All rights
@@ -748,11 +748,16 @@ END SUBROUTINE W3SDS6
   !>
   !> @brief Numerical approximation for the reduction factor.
   !>
-  !> @details Numerical approximation for the reduction factor LFACTOR(f) to
-  !>  reduce energy in the high-frequency part of the resolved part
+  !> @details Numerical approximation for the reduction factor LFACTOR(f)
+  !>  to reduce energy in the high-frequency part of the resolved part
   !>  of the spectrum to meet the constraint on total stress (TAU).
   !>  The constraint is TAU <= TAU_TOT (TAU_TOT = TAU_WAV + TAU_VIS),
   !>  thus the wind input is reduced to match our constraint.
+  !>  This subroutine will compute the charnock parameter
+  !>  CHARN = CHKMIN / SQRT(1-TAU_WAV/TAU_TOT) with CHKMIN being the
+  !>  minimum value. To allow for a reduction of surface drag at high
+  !>  wind speeds a threshold based minimum charnock parameter is
+  !>  optional.
   !>
   !> @param[in]  S      Wind input energy density spectrum.
   !> @param[in]  CINV   Inverse phase speed.
@@ -817,10 +822,16 @@ SUBROUTINE LFACTOR(S, CINV, U10, USTAR, USDIR, SIG, DSII, &
     !                          LFACT(F) = MIN(1,exp((1-U/C(F))*RTAU))
     !        Then alter RTAU and repeat 3) until our constraint is matched.
     !     4) Charnock parameter after equation (3.47) (Komen el al, 1994):
-    !                                         SIN6AHAT
-    !                          CHARN = -----------------------
-    !                                  SQRT( 1.0 - TAU_W/TAU )
-    !
+    !                            CHKMIN
+    !           CHARN = ---------------------------
+    !                   SQRT( 1.0 - TAU_W/TAU_TOT )
+    !
+    !        OPTIONAL: To allow for a reduction of surface drag at high
+    !        wind speeds a threshold based minimum charnock parameter
+    !        CHKMIN is adopted (Breivik et al, 2022, JGR):
+    !                                                            U-UCAP
+    !           CHKMIN = CHKINF + 0.5(CHKMIN - CHKINF)*(1 - TANH ------)
+    !                                                            DELTA
     !  3. Parameters :
     !
     !     Parameter list
@@ -857,7 +868,8 @@ SUBROUTINE LFACTOR(S, CINV, U10, USTAR, USDIR, SIG, DSII, &
     !/
     USE CONSTANTS, ONLY: DAIR, GRAV, TPI
     USE W3GDATMD,  ONLY: NK, NTH, NSPEC, DTH, XFR, ECOS, ESIN
-    USE W3GDATMD,  ONLY: SIN6WS, SIN6AHAT
+    USE W3GDATMD,  ONLY: SIN6WS, SIN6CHKMIN, SIN6FLCAP,             &
+                         SIN6CHKCAP, SIN6CHKINF, SIN6CHKSIG
     USE W3ODATMD,  ONLY: NDST, NDSE, IAPROC, NAPERR
     USE W3TIMEMD,  ONLY: STME21
     USE W3WDATMD,  ONLY: TIME
@@ -893,6 +905,7 @@ SUBROUTINE LFACTOR(S, CINV, U10, USTAR, USDIR, SIG, DSII, &
     REAL              :: TAU_TOT, TAU, TAU_VIS, TAU_WAV
     REAL              :: TAUVX, TAUVY, TAUX, TAUY
     REAL              :: TAU_NND, TAU_INIT(2)
+    REAL              :: CHKMIN
     REAL              :: UPROXY, RTAU, DRTAU, ERR
     LOGICAL           :: OVERSHOT
     CHARACTER(LEN=23) :: IDTIME
@@ -1030,8 +1043,16 @@ SUBROUTINE LFACTOR(S, CINV, U10, USTAR, USDIR, SIG, DSII, &
     !
     LFACT(1:NK) = LF10Hz(1:NK)
     !
-    !/ 4) --- Sea-state depended Charnoc parameter ----------------------- /
-    CHARN = SIN6AHAT / SQRT(1.0 - MIN(TAU_WAV/TAU, 0.999))
+    !/ 4) --- Sea-state depended Charnoc parameter w/ wind speed cap --------- /
+    !
+    CHKMIN = SIN6CHKMIN
+    !
+    IF (SIN6FLCAP) THEN
+      CHKMIN = SIN6CHKINF + 0.5 * (SIN6CHKMIN - SIN6CHKINF) *  &
+             (1.0 - TANH( (U10 - SIN6CHKCAP) / SIN6CHKSIG ))
+    END IF
+    !
+    CHARN = CHKMIN / SQRT(1.0 - MIN(TAU_WAV / TAU_TOT, 0.999))
     !
 #ifdef W3_T6
     WRITE (NDST,273) 'Sin ', IDTIME(1:19), SDENS10Hz*TPI

diff --git a/regtests/ww3_ts1/input_capcha/namelists_CAPCHA.nml b/regtests/ww3_ts1/input_capcha/namelists_CAPCHA.nml
@@ -8,5 +8,8 @@
   TAUWBUG = 1,
   UCAP=28.0
 /
-
+&SIN6
+  CHKMIN=0.0095
+  CHKCAP=28.0
+/
 END OF NAMELISTS
diff --git a/regtests/ww3_ts1/input_capcha/switch_ST6 b/regtests/ww3_ts1/input_capcha/switch_ST6
@@ -0,0 +1 @@
+NOGRB SHRD PR0 FLX4 LN1 ST6 NL1 BT1 DB1 TR0 BS0 IC0 IS0 REF0 WNT1 WNX1 CRT1 CRX1 O0 O1 O2 O3 O4 O5 O6 O7 O10 O11
Original file line number	Diff line number	Diff line change
		@@ -0,0 +1 @@
		NOGRB SHRD PR0 FLX4 LN1 ST6 NL1 BT1 DB1 TR0 BS0 IC0 IS0 REF0 WNT1 WNX1 CRT1 CRX1 O0 O1 O2 O3 O4 O5 O6 O7 O10 O11