From e42793bfa78c53f80579644203c90a4da3757b63 Mon Sep 17 00:00:00 2001 From: ma05cde Date: Tue, 18 Jul 2023 14:40:54 +0000 Subject: [PATCH 01/15] start to integrate ST4 updates --- manual/eqs/ST3.tex | 11 ++++---- manual/eqs/ST4.tex | 20 +++------------ model/src/w3gdatmd.F90 | 8 ++++-- model/src/w3gridmd.F90 | 58 +++++++++++++++++++++++++++--------------- model/src/w3iogrmd.F90 | 17 ++++++++----- 5 files changed, 63 insertions(+), 51 deletions(-) diff --git a/manual/eqs/ST3.tex b/manual/eqs/ST3.tex index 91a7fca10..b287aa2ec 100644 --- a/manual/eqs/ST3.tex +++ b/manual/eqs/ST3.tex @@ -57,16 +57,15 @@ \subsubsection{~$S_{in} + S_{ds}$: \wam\ cycle 4 (ECWAM)} \label{sec:ST3} waves that travel faster than the wind. This accounts for some gustiness in the wind and should possibly be resolution-dependent. For reference, this parameter was not properly set in early versions of the SWAN model, as -discovered by R. Lalbeharry.}. The roughness $z_1$ is defined as, - +discovered by R. Lalbeharry.}. If the friction velocity $u_\star$ is known, +it gives the roughness $z_1$ and the wind speed at altitude $z_u$ (by default $z_u=10$~m), \begin{eqnarray} -U_{10}&=&\frac{u_\star}{\kappa} \log\left(\frac{z_u}{z_1}\right) \\ -z_1&=&\alpha_0 \frac{\tau}{ \sqrt{1-\tau_w/\tau}}, +z_1&=&\alpha_0 \frac{\tau}{ \sqrt{1-\tau_w/u_\star^2}}, \\ +U(z_u)&=&\frac{u_\star}{\kappa} \log\left(\frac{z_u}{z_1}\right) \end{eqnarray} \noindent -where $\tau=u_\star^2$, and $z_u$ is the height at which the wind is -specified. These two equations provide an implicit functional dependence of +In practice these two equations provide an implicit functional dependence of $u_\star$ on $U_{10}$ and $\tau_w/\tau$. This relationship is then tabulated \citep{art:Jan91, rep:Bea07}. diff --git a/manual/eqs/ST4.tex b/manual/eqs/ST4.tex index 624ac6af8..733ec09e7 100644 --- a/manual/eqs/ST4.tex +++ b/manual/eqs/ST4.tex @@ -8,10 +8,10 @@ \subsubsection{~$S_{\mathrm{in}} + S_{\mathrm{ds}}$: Saturation-based dissipatio This family of parameterizations uses a positive part of the wind input taken from WAM cycle 4 with an ad hoc reduction of $u_\star$, implemented in order to allow a balance with a saturation-based dissipation that uses different options for -a cumulative term. There are three main options for defining the saturation and the cumulative term. Chosing one or the other is done with the {\F SDSBCHOICE} parameter, with {\F SDSBCHOICE=1} for \cite{art:Aea10}, {\F SDSBCHOICE=2} for \cite{Filipot&Ardhuin2012}, and {\F SDSBCHOICE=3} for \cite{Romero2019}. That last options uses a saturation that is defined from the local spectral density, and thus gives zero dissipation for directions where the threshold is not reached, leading to much broader directional spectra. Also the stronger bimodality is achieved by having a strong modulation effect as a cumulative term. +a cumulative term. There are three main options for defining the saturation and the cumulative term. Chosing one or the other is done with the {\F SDSBCHOICE} parameter, with {\F SDSBCHOICE=1} for \cite{art:Aea10}, {\F SDSBCHOICE=2} for \cite{Filipot&Ardhuin2012}, and {\F SDSBCHOICE=3} for \cite{Romero2019} and later adjustments including \cite{art:AA23}. That last option uses a saturation that is defined from the local spectral density, and thus gives zero dissipation for directions where the threshold is not reached, leading to much broader directional spectra. Also the stronger bimodality is achieved by having a strong modulation effect as a cumulative term. Many other adjustments can be made by changing the namelist parameters. A few successful combinations -are given by tables \ref{tab:ST4_parSIN} and \ref{tab:ST4_parSDS}, with results described by \citep{art:RA13,art:SAG16}. +are given by tables \ref{tab:ST4_parSIN} and \ref{tab:ST4_parSDS}, with results described by \citep{art:RA13,art:SAG16,art:AA23}. Further calibration to any particular wind field should be done for best performance. Guidance for this is given by \cite{Stopa2018}. %We also note that the particular %set of parameters T400 corresponds to setting IPHYS=1 in the ECWAM code cycle 45R2, with a few differences @@ -216,27 +216,15 @@ \subsubsection{~$S_{\mathrm{in}} + S_{\mathrm{ds}}$: Saturation-based dissipatio direction will typically produce less dissipation than a sea state with all the energy radiated in the same direction. -Based on recent analysis by \cite{Guimaraes2018} and \cite{Peureux&al.2019}, this saturation is enhanced by a factor $M_L$ that represents -the effect of long waves on short waves -\begin{equation} -M_l(k,\theta)=1+M_\theta \sqrt{\mathrm{mss}(k,\theta)} + N_\theta \sqrt{\mathrm{nss}(k,\theta)} \label{defFACSAT}. -\end{equation} -where $M_\theta$ is twice the modulation transfer function for short wave steepness, with -$M_\theta=8$ when following the simplified theory by \cite{art:LHS60} and using the root mean square enhancement of $B$ over a -long wave cycle. $N_\theta$ is an additional straining factor due to the instability of the wave action envelope of short waves -propagating in the direction close to that of the long wave \citep{Peureux&al.2019}. The squared slopes $\mathrm{mss}(k,\theta)$ is -the mean square slope in direction $\theta$, wheras $\mathrm{nss}(k,\theta)$ is a slope of long waves propagating in a narrow window $\pm \delta_\theta$, -around the short wave direction $\theta$. - We finally define our dissipation term as the sum of the saturation-based term and a cumulative breaking term $S_{\mathrm{bk,cu}}$, \begin{eqnarray} \cS_{ds}(k,\theta)& =& \sigma \frac{C_{\mathrm{ds}}^{\mathrm{sat}}}{B^2_r} \left[ \delta_d -\max\left\{ M_l(k,\theta) B\left(k\right) - +\max\left\{ B\left(k\right) - B_r,0\right\}^2 \right. \nonumber \\ - & & + \left(1-\delta_d \right) \left. \max\left\{ M_L(k,\theta) B'\left(k,\theta \right)- B_r + & & + \left(1-\delta_d \right) \left. \max\left\{B'\left(k,\theta \right)- B_r ,0\right\}^2\right]N(k,\theta) \nonumber \\ & & + \cS_{\mathrm{bk,cu}}(k,\theta) + \cS_{\mathrm{turb}}(k,\theta) \label{Sds_all}. \end{eqnarray} diff --git a/model/src/w3gdatmd.F90 b/model/src/w3gdatmd.F90 index 6cd6e91d8..92fc3e040 100644 --- a/model/src/w3gdatmd.F90 +++ b/model/src/w3gdatmd.F90 @@ -1304,7 +1304,7 @@ MODULE W3GDATMD FFXFM, FFXPM, SSDSBRF1, SSDSBRF2, & SSDSBINT, SSDSBCK, SSDSHCK, SSDSABK, & SSDSPBK, SSINBR,SSINTHP,TTAUWSHELTER,& - SSWELLF(:), SSDSC(:), SSDSBR, & + SINTAILPAR(:), SSWELLF(:), SSDSC(:), SSDSBR, & SSDSP, WWNMEANP, SSTXFTF, SSTXFTWN, & SSDSBT, SSDSCOS, SSDSDTH, SSDSBM(:) #endif @@ -2062,7 +2062,10 @@ SUBROUTINE W3DIMS ( IMOD, MK, MTH, NDSE, NDST ) SDSNTH = MTH/2-1 !MIN(NINT(SSDSDTH/(DTH*RADE)),MTH/2-1) ALLOCATE( MPARS(IMOD)%SRCPS%SATINDICES(2*SDSNTH+1,MTH), & MPARS(IMOD)%SRCPS%SATWEIGHTS(2*SDSNTH+1,MTH), & - MPARS(IMOD)%SRCPS%CUMULW(MSPEC,MSPEC), & + STAT=ISTAT ) + CHECK_ALLOC_STATUS ( ISTAT ) + ! IF (SSDSC(3).LT.0.) & + ALLOCATE (MPARS(IMOD)%SRCPS%CUMULW(MSPEC,MSPEC), & STAT=ISTAT ) CHECK_ALLOC_STATUS ( ISTAT ) #endif @@ -2633,6 +2636,7 @@ SUBROUTINE W3SETG ( IMOD, NDSE, NDST ) ZZ0RAT => MPARS(IMOD)%SRCPS%ZZ0RAT ZZALP => MPARS(IMOD)%SRCPS%ZZALP TTAUWSHELTER => MPARS(IMOD)%SRCPS%TTAUWSHELTER + SINTAILPAR => MPARS(IMOD)%SRCPS%SINTAILPAR SSWELLFPAR => MPARS(IMOD)%SRCPS%SSWELLFPAR SSWELLF => MPARS(IMOD)%SRCPS%SSWELLF SSDSC => MPARS(IMOD)%SRCPS%SSDSC diff --git a/model/src/w3gridmd.F90 b/model/src/w3gridmd.F90 index 281ed7a3f..23d779d65 100644 --- a/model/src/w3gridmd.F90 +++ b/model/src/w3gridmd.F90 @@ -835,7 +835,8 @@ MODULE W3GRIDMD #endif ! #ifdef W3_ST4 - INTEGER :: SWELLFPAR, SDSISO, SDSBRFDF + INTEGER :: SWELLFPAR, SDSISO, SDSBRFDF, SINTABLE,& + TAUWBUG, VISCSTRESS REAL :: SDSBCHOICE REAL :: ZWND, ALPHA0, Z0MAX, BETAMAX, SINTHP,& ZALP, Z0RAT, TAUWSHELTER, SWELLF, & @@ -851,7 +852,8 @@ MODULE W3GRIDMD SDSBRF1, & SDSBM0, SDSBM1, SDSBM2, SDSBM3, & SDSBM4, SDSFACMTF, SDSCUMP, SDSNUW, & - SDSL, SDSMWD, SDSMWPOW, SPMSS, SDSNMTF + SDSL, SDSMWD, SDSMWPOW, SPMSS, SDSNMTF, SINTAIL1, SINTAIL2, & + CUMSIGP #endif ! #ifdef W3_ST6 @@ -991,7 +993,7 @@ MODULE W3GRIDMD NAMELIST /SIN4/ ZWND, ALPHA0, Z0MAX, BETAMAX, SINTHP, ZALP, & TAUWSHELTER, SWELLFPAR, SWELLF, & SWELLF2, SWELLF3, SWELLF4, SWELLF5, SWELLF6, & - SWELLF7, Z0RAT, SINBR + SWELLF7, Z0RAT, SINBR, SINTABLE, SINTAIL1, SINTAIL2, TAUWBUG, VISCSTRESS #endif #ifdef W3_NL1 NAMELIST /SNL1/ LAMBDA, NLPROP, KDCONV, KDMIN, & @@ -1031,7 +1033,7 @@ MODULE W3GRIDMD SDSC5, SDSC6, SDSBR, SDSBT, SDSP, SDSISO, & SDSBCK, SDSABK, SDSPBK, SDSBINT, SDSHCK, & SDSDTH, SDSCOS, SDSBRF1, SDSBRFDF, SDSNUW, & - SDSBM0, SDSBM1, SDSBM2, SDSBM3, SDSBM4, & + SDSBM0, SDSBM1, SDSBM2, SDSBM3, SDSBM4, CUMSIGP,& WHITECAPWIDTH, WHITECAPDUR, SDSMWD, SDSMWPOW, SDKOF #endif @@ -1710,6 +1712,12 @@ SUBROUTINE W3GRID() TAUWSHELTER = 0.3 ZALP = 0.006 SINBR = 0. + SINTABLE = 1 + SINTAIL1 = 0. ! TAUWSHELTER FOR TAIL (no table) + SINTAIL2 = 0. ! additional peak in capillary range + TAUWBUG = 1 ! TAUWBUG is 1 is the bug is kept: + ! initializes TAUWX/Y to zero in W3SRCE + VISCSTRESS =0 #endif ! #ifdef W3_ST6 @@ -1793,6 +1801,11 @@ SUBROUTINE W3GRID() SSWELLF(6) = SWELLF6 SSWELLF(7) = SWELLF7 SSWELLFPAR = SWELLFPAR + SINTAILPAR(1) = FLOAT(SINTABLE) + SINTAILPAR(2) = SINTAIL1 + SINTAILPAR(3) = SINTAIL2 + SINTAILPAR(4) = FLOAT(TAUWBUG) + SINTAILPAR(5) = VISCSTRESS #endif ! #ifdef W3_ST6 @@ -2074,8 +2087,8 @@ SUBROUTINE W3GRID() SDSDTH = 80. SDSCOS = 2. SDSISO = 2 - SDSBM0 = 1. - SDSBM1 = 0. + SDSBM0 = 1. ! All these parameters are related to finite depth + SDSBM1 = 0. ! scaling of breaking SDSBM2 = 0. SDSBM3 = 0. SDSBM4 = 0. @@ -2085,8 +2098,9 @@ SUBROUTINE W3GRID() SDSBINT = 0.3 SDSHCK = 1.5 WHITECAPWIDTH = 0.3 - SDSSTRAIN = 0. SDSFACMTF = 400 ! MTF factor for Lambda , Romero (2019) + CUMSIGP = 0. + SDSSTRAIN = 0. SDSSTRAINA = 15. SDSSTRAIN2 = 0. WHITECAPDUR = 0.56 ! breaking duration factor @@ -2097,7 +2111,7 @@ SUBROUTINE W3GRID() ! MTF SPMSS = 0.5 ! cmss^SPMSS SDSNMTF = 1.5 ! MTF power - SDSCUMP = 2. + SDSCUMP = 2. ! 2 for cumulative mss, 1 for cumulative orb. vel. ! MW SDSMWD = .9 ! new AFo SDSMWPOW = 1. ! (k )^pow @@ -2179,9 +2193,9 @@ SUBROUTINE W3GRID() SSDSC(7) = WHITECAPWIDTH SSDSC(8) = SDSSTRAIN ! Straining constant ... SSDSC(9) = SDSL - SSDSC(10) = SDSSTRAINA*NTH/360. ! angle Aor enhanced straining + SSDSC(10) = SDSSTRAINA*NTH/360. ! angle for enhanced straining SSDSC(11) = SDSSTRAIN2 ! straining constant for directional part - SSDSC(12) = SDSBT + SSDSC(12) = CUMSIGP SSDSC(13) = SDSMWD SSDSC(14) = SPMSS SSDSC(15) = SDSMWPOW @@ -3165,7 +3179,7 @@ SUBROUTINE W3GRID() #ifdef W3_ST4 WRITE (NDSO,2920) ZWND, ALPHA0, Z0MAX, BETAMAX, SINTHP, ZALP, & TAUWSHELTER, SWELLFPAR, SWELLF, SWELLF2, SWELLF3, SWELLF4, & - SWELLF5, SWELLF6, SWELLF7, Z0RAT, SINBR + SWELLF5, SWELLF6, SWELLF7, Z0RAT, SINBR, SINTABLE, TAUWBUG, VISCSTRESS, SINTAIL1, SINTAIL2 #endif #ifdef W3_ST6 WRITE (NDSO,2920) SINA0, SINWS, SINFC @@ -3227,7 +3241,7 @@ SUBROUTINE W3GRID() SDSBT, SDSP, SDSISO, SDSCOS, SDSDTH, SDSBRF1, & SDSBRFDF, SDSBM0, SDSBM1, SDSBM2, SDSBM3, SDSBM4, & SPMSS, SDKOF, SDSMWD, SDSFACMTF, SDSNMTF,SDSMWPOW,& - SDSCUMP, SDSNUW, WHITECAPWIDTH, WHITECAPDUR + SDSCUMP, CUMSIGP, SDSNUW, WHITECAPWIDTH, WHITECAPDUR #endif #ifdef W3_ST6 WRITE (NDSO,2924) SDSET, SDSA1, SDSA2, SDSP1, SDSP2 @@ -3270,7 +3284,7 @@ SUBROUTINE W3GRID() JGS_TERMINATE_DIFFERENCE, & JGS_TERMINATE_NORM, & JGS_LIMITER, & - JGS_LIMITER_FUNC, & + JGS_LIMITER_FUNC, & JGS_USE_JACOBI, & JGS_BLOCK_GAUSS_SEIDEL, & JGS_MAXITER, & @@ -3607,7 +3621,7 @@ SUBROUTINE W3GRID() END SELECT IF (FSTOTALIMP .or. FSTOTALEXP) THEN - LPDLIB = .TRUE. + LPDLIB = .TRUE. ENDIF ! IF (SUM(UNSTSCHEMES).GT.1) WRITE(NDSO,1035) @@ -6192,7 +6206,9 @@ SUBROUTINE W3GRID() ' SWELLF =',F8.5,', SWELLF2 =',F8.5, & ', SWELLF3 =',F8.5,', SWELLF4 =',F9.1,','/ & ' SWELLF5 =',F8.5,', SWELLF6 =',F8.5, & - ', SWELLF7 =',F12.2,', Z0RAT =',F8.5,', SINBR =',F8.5,' /') + ', SWELLF7 =',F12.2,', Z0RAT =',F8.5,', SINBR =',F8.5,','/ & + ' SINTABLE =',I2,', TAUWBUG =',I2, & + ', VISCSTRESS =',F8.5,', SINTAIL1 =',F8.5,', SINTAIL2 =',F8.5,' /') #endif ! #ifdef W3_ST6 @@ -6361,7 +6377,7 @@ SUBROUTINE W3GRID() ' SPMSS = ',F5.2, ', SDKOF =',F5.2, & ', SDSMWD =',F5.2,', SDSFACMTF =',F5.1,', '/ & ' SDSMWPOW =',F3.1,', SDSNMTF =', F5.2, & - ', SDSCUMP =', F3.1,', SDSNUW =', E8.3,', '/, & + ', SDSCUMP =', F3.1,', CUMSIGP =', F3.1,', SDSNUW =', E10.3,', '/, & ' WHITECAPWIDTH =',F5.2, ' WHITECAPDUR =',F5.2,' /') #endif ! @@ -6473,12 +6489,12 @@ SUBROUTINE W3GRID() 947 FORMAT (/' Ice scattering ',A,/ & ' --------------------------------------------------') 948 FORMAT (' IS2 Scattering ... '/& - ' scattering coefficient : ',E9.3/ & - ' 0: no back-scattering : ',E9.3/ & + ' scattering coefficient : ',E10.3/ & + ' 0: no back-scattering : ',E10.3/ & ' TRUE: istropic back-scattering : ',L3/ & ' TRUE: update of ICEDMAX : ',L3/ & ' TRUE: keeps updated ICEDMAX : ',L3/ & - ' flexural strength : ',E9.3/ & + ' flexural strength : ',E10.3/ & ' TRUE: uses Robinson-Palmer disp.: ',L3/ & ' attenuation : ',F5.2/ & ' fragility : ',F5.2/ & @@ -6486,7 +6502,7 @@ SUBROUTINE W3GRID() ' pack scattering coef 1 : ',F5.2/ & ' pack scattering coef 2 : ',F5.2/ & ' scaling by concentration : ',F5.2/ & - ' creep B coefficient : ',E9.3/ & + ' creep B coefficient : ',E10.3/ & ' creep C coefficient : ',F5.2/ & ' creep D coefficient : ',F5.2/ & ' creep N power : ',F5.2/ & @@ -6497,7 +6513,7 @@ SUBROUTINE W3GRID() ' energy of activation : ',F5.2/ & ' anelastic coefficient : ',E11.3/ & ' anelastic exponent : ',F5.2) -2948 FORMAT ( ' &SIS2 ISC1 =',E9.3,', IS2BACKSCAT =',E9.3, & +2948 FORMAT ( ' &SIS2 ISC1 =',E10.3,', IS2BACKSCAT =',E10.3, & ', IS2ISOSCAT =',L3,', IS2BREAK =',L3, & ', IS2DUPDATE =',L3,','/ & ' IS2FLEXSTR =',E11.3,', IS2DISP =',L3, & diff --git a/model/src/w3iogrmd.F90 b/model/src/w3iogrmd.F90 index 3aa2688ab..d954d41c4 100644 --- a/model/src/w3iogrmd.F90 +++ b/model/src/w3iogrmd.F90 @@ -1242,10 +1242,11 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT ) SSTXFTFTAIL, SSTXFTWN, SSTXFTF, SSTXFTWN, & SSDSBRF1, SSDSBRF2, SSDSBRFDF,SSDSBCK, SSDSABK, & SSDSPBK, SSDSBINT, FFXPM, FFXFM, FFXFA, & - SSDSHCK, DELUST, DELTAIL, DELTAUW, & - DELU, DELALP, TAUT, TAUHFT, TAUHFT2, & + SSDSHCK, & IKTAB, DCKI, QBI, SATINDICES, SATWEIGHTS, & - DIKCUMUL, CUMULW + DIKCUMUL, CUMULW, SINTAILPAR + IF (SINTAILPAR(1).GT.0.5) WRITE (NDSM) DELUST, DELTAIL,& + DELTAUW, DELU, DELALP, TAUT, TAUHFT, TAUHFT2 ELSE READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & ZZWND, AALPHA, ZZ0MAX, BBETA, SSINTHP, ZZALP, & @@ -1256,10 +1257,14 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT ) SSTXFTFTAIL, SSTXFTWN, SSTXFTF, SSTXFTWN, & SSDSBRF1, SSDSBRF2, SSDSBRFDF,SSDSBCK, SSDSABK, & SSDSPBK, SSDSBINT, FFXPM, FFXFM, FFXFA, & - SSDSHCK, DELUST, DELTAIL, DELTAUW, & - DELU, DELALP, TAUT, TAUHFT, TAUHFT2, & + SSDSHCK, & IKTAB, DCKI, QBI, SATINDICES, SATWEIGHTS, & - DIKCUMUL, CUMULW + DIKCUMUL, CUMULW, SINTAILPAR + IF (SINTAILPAR(1).GT.0.5) THEN + CALL INSIN4(.FALSE.) + READ (NDSM) DELUST, DELTAIL, & + DELTAUW, DELU, DELALP, TAUT, TAUHFT, TAUHFT2 + END IF END IF #endif ! From c9876661c97635bc0f8339fc5c198009a1bc4284 Mon Sep 17 00:00:00 2001 From: ma05cde Date: Tue, 18 Jul 2023 15:15:25 +0000 Subject: [PATCH 02/15] import ST4 changes --- model/src/w3iorsmd.F90 | 6 +- model/src/w3src4md.F90 | 489 ++++++++++-------- model/src/w3srcemd.F90 | 21 +- model/src/ww3_ounp.F90 | 2 +- regtests/ww3_ts1/input/namelists_ST4_T700.nml | 2 +- 5 files changed, 282 insertions(+), 238 deletions(-) diff --git a/model/src/w3iorsmd.F90 b/model/src/w3iorsmd.F90 index 3bd2aa4ea..05f7e9163 100644 --- a/model/src/w3iorsmd.F90 +++ b/model/src/w3iorsmd.F90 @@ -628,7 +628,7 @@ SUBROUTINE W3IORS ( INXOUT, NDSR, DUMFPI, IMOD, FLRSTRT ) ! Original non-server version writing of spectra ! IF ( .NOT.IOSFLG .OR. (NAPROC.EQ.1.AND.NAPRST.EQ.1) ) THEN -#ifdef W3_MPI +#ifdef W3_MPI DO JSEA=1, NSEAL CALL INIT_GET_ISEA(ISEA, JSEA) NREC = ISEA + 2 @@ -637,7 +637,7 @@ SUBROUTINE W3IORS ( INXOUT, NDSR, DUMFPI, IMOD, FLRSTRT ) WRITEBUFF(1:NSPEC) = VA(1:NSPEC,JSEA) WRITE (NDSR,POS=RPOS,ERR=803,IOSTAT=IERR) WRITEBUFF END DO -#else +#else DO JSEA=1, NSEA ISEA = JSEA NREC = ISEA + 2 @@ -646,7 +646,7 @@ SUBROUTINE W3IORS ( INXOUT, NDSR, DUMFPI, IMOD, FLRSTRT ) WRITEBUFF(1:NSPEC) = VA(1:NSPEC,JSEA) WRITE (NDSR,POS=RPOS,ERR=803,IOSTAT=IERR) WRITEBUFF END DO -#endif +#endif ! ! I/O server version writing of spectra ( !/MPI ) ! diff --git a/model/src/w3src4md.F90 b/model/src/w3src4md.F90 index e2bf12c9a..2a7ebd1a9 100644 --- a/model/src/w3src4md.F90 +++ b/model/src/w3src4md.F90 @@ -39,6 +39,7 @@ MODULE W3SRC4MD !/ 02-Sep-2011 : Clean up and time optimization ( version 4.04 ) !/ 04-Sep-2011 : Estimation of whitecap stats. ( version 4.04 ) !/ 13-Nov-2013 : Reduced frequency range with IG ( version 4.13 ) + !/ 01-Mar-2023 : Clean up of SDS4 ( version 7.14 ) !/ ! 1. Purpose : ! @@ -90,11 +91,9 @@ MODULE W3SRC4MD !air kinematic viscosity (used in WAM) INTEGER, PARAMETER :: ITAUMAX=200,JUMAX=200 INTEGER, PARAMETER :: IUSTAR=100,IALPHA=200, ILEVTAIL=50 - REAL :: TAUT(0:ITAUMAX,0:JUMAX), DELTAUW, DELU - ! Table for H.F. stress as a function of 2 variables - REAL :: TAUHFT(0:IUSTAR,0:IALPHA), DELUST, DELALP - ! Table for H.F. stress as a function of 3 variables - REAL :: TAUHFT2(0:IUSTAR,0:IALPHA,0:ILEVTAIL) + ! Tables for total stress and H.F. stress as a function of 2 or 3 variables + REAL, ALLOCATABLE :: TAUT(:,:),TAUHFT(:,:),TAUHFT2(:,:,:) + REAL :: DELUST, DELALP,DELTAUW, DELU ! Table for swell damping REAL :: DELTAIL REAL, PARAMETER :: UMAX = 50. @@ -231,7 +230,7 @@ SUBROUTINE W3SPR4 (A, CG, WN, EMEAN, FMEAN, FMEAN1, WNMEAN, & USE W3GDATMD, ONLY: NK, NTH, NSPEC, SIG, DTH, DDEN, WWNMEANP, & WWNMEANPTAIL, FTE, FTF, SSTXFTF, SSTXFTWN,& SSTXFTFTAIL, SSWELLF, ESIN, ECOS, AAIRCMIN, & - AAIRGB, AALPHA, ZZWND + AAIRGB, AALPHA, ZZWND, SSDSC #ifdef W3_S USE W3SERVMD, ONLY: STRACE #endif @@ -267,7 +266,7 @@ SUBROUTINE W3SPR4 (A, CG, WN, EMEAN, FMEAN, FMEAN1, WNMEAN, & #endif REAL :: TAUW, EBAND, EMEANWS,UNZ, & - EB(NK),EB2(NK),ELCS, ELSN + EB(NK),EB2(NK),ELCS, ELSN, SIGFAC !/ !/ ------------------------------------------------------------------- / !/ @@ -294,17 +293,18 @@ SUBROUTINE W3SPR4 (A, CG, WN, EMEAN, FMEAN, FMEAN1, WNMEAN, & DO IK=1, NK EB(IK) = 0. EB2(IK) = 0. + SIGFAC=SIG(IK)**SSDSC(12) * DDEN(IK) / CG(IK) DO ITH=1, NTH IS=ITH+(IK-1)*NTH EB(IK) = EB(IK) + A(ITH,IK) - ELCS = ELCS + A(ITH,IK)*ECOS(IS)*DDEN(IK) / CG(IK) - ELSN = ELSN + A(ITH,IK)*ESIN(IS)*DDEN(IK) / CG(IK) + ELCS = ELCS + A(ITH,IK)*ECOS(IS)*SIGFAC + ELSN = ELSN + A(ITH,IK)*ESIN(IS)*SIGFAC IF (LLWS(IS)) EB2(IK) = EB2(IK) + A(ITH,IK) AMAX = MAX ( AMAX , A(ITH,IK) ) END DO END DO - - DLWMEAN=ATAN2(ELSN,ELCS); + ! + DLWMEAN=ATAN2(ELSN,ELCS) ! ! 2. Integrate over directions -------------------------------------- * ! @@ -358,7 +358,6 @@ SUBROUTINE W3SPR4 (A, CG, WN, EMEAN, FMEAN, FMEAN1, WNMEAN, & CALL W3FLX5 ( ZZWND, U, UDIR, TAUA, TAUADIR, DAIR, & USTAR, USDIR, Z0, CD, CHARN ) #else - Z0=0. CALL CALC_USTAR(U,TAUW,USTAR,Z0,CHARN) UNZ = MAX ( 0.01 , U ) CD = (USTAR/UNZ)**2 @@ -510,7 +509,7 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & USE W3GDATMD, ONLY: NK, NTH, NSPEC, DDEN, SIG, SIG2, TH, & ESIN, ECOS, EC2, ZZWND, AALPHA, BBETA, ZZALP,& TTAUWSHELTER, SSWELLF, DDEN2, DTH, SSINTHP, & - ZZ0RAT, SSINBR + ZZ0RAT, SSINBR, SINTAILPAR #ifdef W3_S USE W3SERVMD, ONLY: STRACE #endif @@ -561,13 +560,21 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & REAL XI,DELI1,DELI2 REAL XJ,DELJ1,DELJ2 REAL XK,DELK1,DELK2 - REAL :: CONST, CONST0, CONST2, TAU1 + REAL :: CONST, CONST0, CONST2, TAU1, TAU1NT, ZINF, TENSK REAL X,ZARG,ZLOG,UST REAL :: COSWIND, XSTRESS, YSTRESS, TAUHF REAL TEMP, TEMP2 INTEGER IND,J,I,ISTAB REAL DSTAB(3,NSPEC), DVISC, DTURB REAL STRESSSTAB(3,2),STRESSSTABN(3,2) + ! + INTEGER, PARAMETER :: JTOT=50 + REAL , PARAMETER :: KM=363.,CMM=0.2325 ! K and C at phase speed minimum in rad/m + REAL :: OMEGACC, OMEGA, ZZ0, ZX, ZBETA, USTR, TAUR, & + CONST1, LEVTAIL0, X0, Y, DELY, YC, ZMU, & + LEVTAIL, CGTAIL, ALPHAM, FM, ALPHAT, FMEAN + + REAL, ALLOCATABLE :: W(:) #ifdef W3_T0 REAL :: DOUT(NK,NTH) #endif @@ -591,6 +598,11 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & STRESSSTAB =0. STRESSSTABN =0. ! + ! Coupling coefficient times density ratio DRAT + ! + CONST1=BBETA/KAPPA**2 ! needed for the tail + CONST0=CONST1*DRAT ! needed for the resolved spectrum + ! ! 1.a estimation of surface roughness parameters ! Z0VISC = 0.1*nu_air/MAX(USTAR,0.0001) @@ -615,9 +627,9 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & ! At this point UORB and AORB are the variances of the orbital velocity and surface elevation ! UORB = UORB + EB *SIG(IK)**2 * DDEN(IK) / CG(IK) - AORB = AORB + EB * DDEN(IK) / CG(IK) !deep water only + AORB = AORB + EB * DDEN(IK) / CG(IK) !correct for deep water only END DO - + ! FMEAN = SQRT((UORB+1E-6)/(AORB+1E-6)) UORB = 2*SQRT(UORB) ! significant orbital amplitude AORB1 = 2*AORB**(1-0.5*SSWELLF(6)) ! half the significant wave height ... if SWELLF(6)=1 RE = 4*UORB*AORB1 / NU_AIR ! Reynolds number @@ -695,10 +707,6 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & STRESSSTAB(ISTAB,:)=0. STRESSSTABN(ISTAB,:)=0. ! - ! Coupling coefficient times density ratio DRAT - ! - CONST0=BBETA*DRAT/(kappa**2) - ! DO IK=1, NK TAUPX=TAUX-ABS(TTAUWSHELTER)*STRESSSTAB(ISTAB,1) TAUPY=TAUY-ABS(TTAUWSHELTER)*STRESSSTAB(ISTAB,2) @@ -813,13 +821,22 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & DOUT(IK,ITH) = D(ITH+(IK-1)*NTH) END DO END DO - CALL PRT2DS (NDST, NK, NK, NTH, DOUT, SIG(1:NK), ' ', 1., & + CALL PRT2DS (NDST, NK, NK, NTH, DOUT, SIG(1), ' ', 1., & 0.0, 0.001, 'Diag Sin', ' ', 'NONAME') #endif ! #ifdef W3_T1 CALL OUTMAT (NDST, D, NTH, NTH, NK, 'diag Sin') #endif + ! + TAUPX=TAUX-ABS(TTAUWSHELTER)*XSTRESS + TAUPY=TAUY-ABS(TTAUWSHELTER)*YSTRESS + USTP=(TAUPX**2+TAUPY**2)**0.25 + USDIRP=ATAN2(TAUPY,TAUPX) + + UST=USTP + ! + ! Computes HF tail ! ! Computes the high-frequency contribution ! the difference in spectal density (kx,ky) to (f,theta) @@ -832,36 +849,112 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & COSWIND=(ECOS(IS)*COSU+ESIN(IS)*SINU) TEMP=TEMP+A(IS)*(MAX(COSWIND,0.))**3 END DO + ! + LEVTAIL0= CONST0*TEMP ! LEVTAIL is sum of A(k,theta)*cos^3(theta-wind)*DTH*SIG^5/(g^2*2pi)*2*pi*SIG/CG + ! ! which is the same as sum of E(f,theta)*cos^3(theta-wind)*DTH*SIG^5/(g^2*2pi) + ! reminder: sum of E(f,theta)*DTH*SIG^5/(g^2*2pi) is 2*k^3*E(k) + IF (SINTAILPAR(1).LT.0.5) THEN + ALLOCATE(W(JTOT)) + W(2:JTOT-1)=1. + W(1)=0.5 + W(JTOT)=0.5 + X0 = 0.05 + ! + USTR= UST + ZZ0=Z0 + OMEGACC = MAX(SIG(NK),X0*GRAV/UST) + YC = OMEGACC*SQRT(ZZ0/GRAV) - TAUPX=TAUX-ABS(TTAUWSHELTER)*XSTRESS - TAUPY=TAUY-ABS(TTAUWSHELTER)*YSTRESS - USTP=(TAUPX**2+TAUPY**2)**0.25 - USDIRP=ATAN2(TAUPY,TAUPX) + ! DELY = MAX((1.-YC)/REAL(JTOT),0.) + ! Changed integration variable from Y to LOG(Y) and to log(K) + !ZINF = LOG(YC) + !DELY = MAX((1.-ZINF)/REAL(JTOT),0.) + ZINF = LOG(SIG(NK)**2/GRAV) + DELY = (LOG(TPI/0.005)-ZINF)/REAL(JTOT) - UST=USTP - ! finds the values in the tabulated stress TAUHFT - XI=UST/DELUST - IND = MAX(1,MIN (IUSTAR-1, INT(XI))) - DELI1= MAX(MIN (1. ,XI-FLOAT(IND)),0.) - DELI2= 1. - DELI1 - XJ=MAX(0.,(GRAV*Z0/MAX(UST,0.00001)**2-AALPHA) / DELALP) - J = MAX(1 ,MIN (IALPHA-1, INT(XJ))) - DELJ1= MAX(0.,MIN (1. , XJ-FLOAT(J))) - DELJ2=1. - DELJ1 - IF (TTAUWSHELTER.GT.0) THEN - XK = CONST0*TEMP / DELTAIL - I = MIN (ILEVTAIL-1, INT(XK)) - DELK1= MIN (1. ,XK-FLOAT(I)) - DELK2=1. - DELK1 - TAU1 =((TAUHFT2(IND,J,I)*DELI2+TAUHFT2(IND+1,J,I)*DELI1 )*DELJ2 & - +(TAUHFT2(IND,J+1,I)*DELI2+TAUHFT2(IND+1,J+1,I)*DELI1)*DELJ1)*DELK2 & - +((TAUHFT2(IND,J,I+1)*DELI2+TAUHFT2(IND+1,J,I+1)*DELI1 )*DELJ2 & - +(TAUHFT2(IND,J+1,I+1)*DELI2+TAUHFT2(IND+1,J+1,I+1)*DELI1)*DELJ1)*DELK1 + TAUR=UST**2 + TAU1=0. + + ! Integration loop over the tail wavenumbers or frequencies ... + DO J=1,JTOT + !Y = YC+REAL(J-1)*DELY + !OMEGA = Y*SQRT(GRAV/ZZ0) + !OMEGA = SQRT(GRAV*Y) + ! This is the deep water phase speed... No surface tension !! + !CM = GRAV/OMEGA + ! With this form, Y is the wavenumber in the tail; + Y= EXP(ZINF+REAL(J-1)*DELY) + TENSK =1+(Y/KM)**2 + OMEGA = SQRT(GRAV*Y*TENSK) + CM = SQRT(GRAV*TENSK/Y) + CGTAIL = 0.5*(3*(Y/KM)**2+1)*SQRT(GRAV/(Y*TENSK)) + !this is the inverse wave age, shifted by ZZALP (tuning) + ZX = USTR/CM +ZZALP + ZARG = MIN(KAPPA/ZX,20.) + ! ZMU corresponds to EXP(ZCN) + ZMU = MIN(GRAV*ZZ0/CM**2*EXP(ZARG),1.) + ZLOG = MIN(ALOG(ZMU),0.) + ZBETA = CONST1*ZMU*ZLOG**4 + ! + ! Optional addition of capillary wave peak + ! + IF (SINTAILPAR(3).GT.0) THEN + IF (USTR.LT.CM) THEN + ALPHAM=MAX(0.,0.01*(1.+ALOG(USTR/CM))) + ELSE + ALPHAM=0.01*(1+3.*ALOG(USTR/CM)) + END IF + FM=EXP(-0.25*(Y/KM-1)**2) + + ALPHAT=ALPHAM*(CMM/CM)*FM ! equivalent to 2*Bh in Elfouhaily et al. + LEVTAIL=LEVTAIL0*0.5*(1-tanh((Y-20)/5))+SINTAILPAR(3)*0.5*(1+TANH((Y-20)/5))*ALPHAT + ELSE + LEVTAIL=LEVTAIL0 + END IF + ! WRITE(991,*) 'TAIL??',SINTAILPAR(3),LEVTAIL0,LEVTAIL,ALPHAT,Y,Y/KM,OMEGA/(TPI) + + !TAU1=TAU1+W(J)*ZBETA*(USTR/UST)**2/Y*DELY ! integration over LOG(Y) + TAU1=TAU1+W(J)*ZBETA*USTR**2*LEVTAIL*DELY*CGTAIL/CM ! integration over LOG(K) + + ! NB: the factor ABS(TTAUWSHELTER) was forgotten in the TAUHFT2 table + !TAUR=TAUR-W(J)*ABS(TTAUWSHELTER)*USTR**2*ZBETA*LEVTAIL/Y*DELY + !TAUR=TAUR-W(J)*USTR**2*ZBETA*LEVTAIL*DELY ! integration over LOG(Y) + TAUR=TAUR-W(J)*SINTAILPAR(2)*USTR**2*ZBETA*LEVTAIL*DELY*CGTAIL/CM ! DK/K*CG/C = D OMEGA / OMEGA + USTR=SQRT(MAX(TAUR,0.)) + END DO + DEALLOCATE(W) + TAU1NT=TAU1 + TAUHF = TAU1 + ! + ! In this case, uses tables for high frequency contribution to TAUW. + ! ELSE - TAU1 =(TAUHFT(IND,J)*DELI2+TAUHFT(IND+1,J)*DELI1 )*DELJ2 & - +(TAUHFT(IND,J+1)*DELI2+TAUHFT(IND+1,J+1)*DELI1)*DELJ1 - END IF - TAUHF = CONST0*TEMP*UST**2*TAU1 + ! finds the values in the tabulated stress TAUHFT + XI=UST/DELUST + IND = MAX(1,MIN (IUSTAR-1, INT(XI))) + DELI1= MAX(MIN (1. ,XI-FLOAT(IND)),0.) + DELI2= 1. - DELI1 + XJ=MAX(0.,(GRAV*Z0/MAX(UST,0.00001)**2-AALPHA) / DELALP) + J = MAX(1 ,MIN (IALPHA-1, INT(XJ))) + DELJ1= MAX(0.,MIN (1. , XJ-FLOAT(J))) + DELJ2=1. - DELJ1 + IF (TTAUWSHELTER.GT.0) THEN + XK = LEVTAIL / DELTAIL + I = MIN (ILEVTAIL-1, INT(XK)) + DELK1= MIN (1. ,XK-FLOAT(I)) + DELK2=1. - DELK1 + TAU1 =((TAUHFT2(IND,J,I)*DELI2+TAUHFT2(IND+1,J,I)*DELI1 )*DELJ2 & + +(TAUHFT2(IND,J+1,I)*DELI2+TAUHFT2(IND+1,J+1,I)*DELI1)*DELJ1)*DELK2 & + +((TAUHFT2(IND,J,I+1)*DELI2+TAUHFT2(IND+1,J,I+1)*DELI1 )*DELJ2 & + +(TAUHFT2(IND,J+1,I+1)*DELI2+TAUHFT2(IND+1,J+1,I+1)*DELI1)*DELJ1)*DELK1 + ELSE + TAU1 =(TAUHFT(IND,J)*DELI2+TAUHFT(IND+1,J)*DELI1 )*DELJ2 & + +(TAUHFT(IND,J+1)*DELI2+TAUHFT(IND+1,J+1)*DELI1)*DELJ1 + END IF + ! + TAUHF = LEVTAIL0*UST**2*TAU1 + END IF ! End of test on use of table + TAUWX = XSTRESS+TAUHF*COS(USDIRP) TAUWY = YSTRESS+TAUHF*SIN(USDIRP) ! @@ -975,7 +1068,8 @@ SUBROUTINE INSIN4(FLTABS) SSDSDTH, SSDSCOS, TH, DTH, XFR, ECOS, ESIN, & SSDSC, SSDSBRF1, SSDSBCK, SSDSBINT, SSDSPBK, & SSDSABK, SSDSHCK, IKTAB, DCKI, SATINDICES, & - SATWEIGHTS, CUMULW, NKHS, NKD, NDTAB, QBI + SATWEIGHTS, CUMULW, NKHS, NKD, NDTAB, QBI, & + SINTAILPAR #ifdef W3_S USE W3SERVMD, ONLY: STRACE #endif @@ -1016,11 +1110,16 @@ SUBROUTINE INSIN4(FLTABS) ! ! These precomputed tables are written in mod_def.ww3 ! - IF (FLTABS) THEN - CALL TABU_STRESS - CALL TABU_TAUHF(SIG(NK) ) !tabulate high-frequency stress: 2D table + IF (SINTAILPAR(1).GT.0.5) THEN + IF (.NOT. ALLOCATED(TAUT)) ALLOCATE(TAUT(0:ITAUMAX,0:JUMAX)) + IF (.NOT. ALLOCATED(TAUHFT)) ALLOCATE(TAUHFT(0:IUSTAR,0:IALPHA)) + IF (FLTABS) THEN + CALL TABU_STRESS + CALL TABU_TAUHF(SIG(NK) ) !tabulate high-frequency stress: 2D table + END IF IF (TTAUWSHELTER.GT.0) THEN - CALL TABU_TAUHF2(SIG(NK) ) !tabulate high-frequency stress: 3D table + IF (.NOT. ALLOCATED(TAUHFT2)) ALLOCATE(TAUHFT2(0:IUSTAR,0:IALPHA,0:ILEVTAIL)) + IF (FLTABS) CALL TABU_TAUHF2(SIG(NK) ) !tabulate high-frequency stress: 3D table END IF END IF ! @@ -1146,7 +1245,7 @@ SUBROUTINE INSIN4(FLTABS) ! Precomputes the weights for the cumulative effect (TEST 441 and 500) ! DIKCUMUL = 0 - IF (SSDSC(3).NE.0) THEN + IF (SSDSC(3).LT.0.) THEN ! DIKCUMUL is the integer difference in frequency bands ! between the "large breakers" and short "wiped-out waves" DIKCUMUL = NINT(SSDSBRF1/(XFR-1.)) @@ -1264,7 +1363,7 @@ SUBROUTINE TABU_STRESS ! ---------------------------------------------------------------------- INTEGER I,J,ITER REAL ZTAUW,UTOP,CDRAG,WCD,USTOLD,TAUOLD - REAL X,UST,ZZ0,ZNU,F,DELF,ZZ00 + REAL X,UST,ZZ0,F,DELF,ZZ00 ! ! DELU = UMAX/FLOAT(JUMAX) @@ -1755,6 +1854,7 @@ SUBROUTINE CALC_USTAR(WINDSPEED,TAUW,USTAR,Z0,CHARN) ! 2. Method : ! ! Computation of u* based on Quasi-linear theory + ! uses Charnock relation with modified roughness Z1=Z0/SQRT(1-TAUW/TAU) ! ! 3. Parameters : ! @@ -1791,8 +1891,8 @@ SUBROUTINE CALC_USTAR(WINDSPEED,TAUW,USTAR,Z0,CHARN) ! ! 10. Source code : !-----------------------------------------------------------------------------! - USE CONSTANTS, ONLY: GRAV, KAPPA - USE W3GDATMD, ONLY: ZZWND, AALPHA + USE CONSTANTS, ONLY: GRAV, KAPPA, NU_AIR + USE W3GDATMD, ONLY: ZZWND, AALPHA, ZZ0MAX, SINTAILPAR #ifdef W3_T USE W3ODATMD, ONLY: NDST #endif @@ -1800,22 +1900,60 @@ SUBROUTINE CALC_USTAR(WINDSPEED,TAUW,USTAR,Z0,CHARN) REAL, intent(in) :: WINDSPEED,TAUW REAL, intent(out) :: USTAR, Z0, CHARN ! local variables - REAL SQRTCDM1 - REAL XI,DELI1,DELI2,XJ,delj1,delj2 - REAL TAUW_LOCAL - INTEGER IND,J - ! - TAUW_LOCAL=MAX(MIN(TAUW,TAUWMAX),0.) - XI = SQRT(TAUW_LOCAL)/DELTAUW - IND = MIN ( ITAUMAX-1, INT(XI)) ! index for stress table - DELI1 = MIN(1.,XI - REAL(IND)) !interpolation coefficient for stress table - DELI2 = 1. - DELI1 - XJ = WINDSPEED/DELU - J = MIN ( JUMAX-1, INT(XJ) ) - DELJ1 = MIN(1.,XJ - REAL(J)) - DELJ2 = 1. - DELJ1 - USTAR=(TAUT(IND,J)*DELI2+TAUT(IND+1,J )*DELI1)*DELJ2 & - + (TAUT(IND,J+1)*DELI2+TAUT(IND+1,J+1)*DELI1)*DELJ1 + REAL :: SQRTCDM1 + REAL :: XI,DELI1,DELI2,XJ,delj1,delj2 ! used for table version + INTEGER :: IND,J + REAL :: TAUW_LOCAL + REAL :: TAUOLD,CDRAG,WCD,USTOLD,X,UST,ZZ0,ZNU,ZZ00,F,DELF + INTEGER, PARAMETER :: NITER=10 + REAL , PARAMETER :: XM=0.50, EPS1=0.00001 + INTEGER :: ITER + ! VARIABLE. TYPE. PURPOSE. + ! *XM* REAL POWER OF TAUW/TAU IN ROUGHNESS LENGTH. + ! *XNU* REAL KINEMATIC VISCOSITY OF AIR. + ! *NITER* INTEGER NUMBER OF ITERATIONS TO OBTAIN TOTAL STRESS + ! *EPS1* REAL SMALL NUMBER TO MAKE SURE THAT A SOLUTION + ! IS OBTAINED IN ITERATION WITH TAU>TAUW. + + ! + IF (SINTAILPAR(1).GT.0.5) THEN + TAUW_LOCAL=MAX(MIN(TAUW,TAUWMAX),0.) + XI = SQRT(TAUW_LOCAL)/DELTAUW + IND = MIN ( ITAUMAX-1, INT(XI)) ! index for stress table + DELI1 = MIN(1.,XI - REAL(IND)) !interpolation coefficient for stress table + DELI2 = 1. - DELI1 + XJ = WINDSPEED/DELU + J = MIN ( JUMAX-1, INT(XJ) ) + DELJ1 = MIN(1.,XJ - REAL(J)) + DELJ2 = 1. - DELJ1 + USTAR=(TAUT(IND,J)*DELI2+TAUT(IND+1,J )*DELI1)*DELJ2 & + + (TAUT(IND,J+1)*DELI2+TAUT(IND+1,J+1)*DELI1)*DELJ1 + ELSE + ! This max is for comparison ... to be removed later + ! TAUW_LOCAL=MAX(MIN(TAUW,TAUWMAX),0.) + TAUW_LOCAL=TAUW + CDRAG = 0.0012875 + WCD = SQRT(CDRAG) + USTOLD = WINDSPEED*WCD + TAUOLD = MAX(USTOLD**2, TAUW_LOCAL+EPS1) + ! Newton method to solve for ustar in U=ustar*log(Z/Z0) + DO ITER=1,NITER + X = TAUW_LOCAL/TAUOLD + UST = SQRT(TAUOLD) + ZZ00=AALPHA*TAUOLD/GRAV + IF (ZZ0MAX.NE.0) ZZ00=MIN(ZZ00,ZZ0MAX) + ! Corrects roughness ZZ00 for quasi-linear effect + ZZ0 = ZZ00/(1.-X)**XM + ZNU = 0.11*nu_air/MAX(UST,1E-6) + ZZ0 = SINTAILPAR(5)*ZNU+ZZ0 + F = UST-KAPPA*WINDSPEED/(ALOG(ZZWND/ZZ0)) + DELF= 1.-KAPPA*WINDSPEED/(ALOG(ZZWND/ZZ0))**2*2./UST & + *(1.-(XM+1)*X)/(1.-X) + UST = UST-F/DELF + TAUOLD= MAX(UST**2., TAUW_LOCAL+EPS1) + END DO + USTAR=UST + END IF ! ! Determines roughness length ! @@ -1832,6 +1970,7 @@ SUBROUTINE CALC_USTAR(WINDSPEED,TAUW,USTAR,Z0,CHARN) END IF CHARN = AALPHA END IF + ! WRITE(6,*) 'CALC_USTAR:',WINDSPEED,TAUW,AALPHA,CHARN,Z0,USTAR ! RETURN END SUBROUTINE CALC_USTAR @@ -1882,10 +2021,11 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & !/ 06-Jun-2018 : Add optional DEBUGSRC ( version 6.04 ) !/ 22-Feb-2020 : Option to use Romero (GRL 2019) ( version 7.06 ) !/ 13-Aug-2021 : Consider DAIR a variable ( version 7.14 ) + !/ 01-Mar-2023 : Clean up of SDS4 ( version 7.xx ) !/ ! 1. Purpose : ! - ! Calculate whitecapping source term and diagonal term of derivative. + ! Calculate wave dissipation source term and diagonal term of derivative. ! ! 2. Method : ! @@ -1949,7 +2089,7 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & SSDSISO, SSDSDTH, SSDSBM, AAIRCMIN, & SSDSBRFDF, SSDSBCK, IKTAB, DCKI, & SATINDICES, SATWEIGHTS, CUMULW, NKHS, NKD, & - NDTAB, QBI + NDTAB, QBI, DSIP, SSDSBRF1,XFR #ifdef W3_IG1 USE W3GDATMD, ONLY: IGPARS #endif @@ -1989,7 +2129,7 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & INTEGER :: IK, IK1, ITH, IK2, JTH, ITH2, & IKHS, IKD, SDSNTH, IT, IKM, NKM INTEGER :: NSMOOTH(NK) - REAL :: C, COSWIND, ASUM, SDIAGISO + REAL :: C, C2, CUMULWISO, COSWIND, ASUM, SDIAGISO REAL :: COEF1, COEF2, COEF4(NK), & COEF5(NK) @@ -2004,19 +2144,15 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & REAL :: FACSAT, DKHS, FACSTRAINB, FACSTRAINL REAL :: BTH0(NK) !saturation spectrum REAL :: BTH(NSPEC) !saturation spectrum - REAL :: BTH0S(NK) !smoothed saturation spectrum - REAL :: BTHS(NSPEC) !smoothed saturation spectrum - INTEGER :: IMSSMAX(NK), NTHSUM - REAL :: MSSSUM(NK,5), WTHSUM(NTH), FACHF - REAL :: MSSSUM2(NK,NTH) - REAL :: MSSLONG(NK,NTH) + REAL :: MSSSUM(NK,5), FACHF + REAL :: MSSLONG REAL :: MSSPCS, MSSPC2, MSSPS2, MSSP, MSSD, MSSTH REAL :: MICHE, X, KLOC #ifdef W3_T0 REAL :: DOUT(NK,NTH) #endif REAL :: QB(NK), S2(NK) - REAL :: TSTR, TMAX, DT, T, MFT + REAL :: TSTR, TMAX, DT, T, MFT, DIRFORCUM REAL :: PB(NSPEC), PB2(NSPEC), BRM12(NK), BTOVER REAL :: KO, LMODULATION(NTH) !/ @@ -2034,10 +2170,10 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & ! found in certain compilers NSMOOTH=0 S1=0.; E1=0. - NTIMES=0;IKSUP=0;IMSSMAX=0 + NTIMES=0;IKSUP=0 DK=0.; HS=0.; KBAR=0.; DCK=0.; EFDF=0. - BTH0=0.; BTH=0.; BTH0S=0.; DDIAG=0.; SRHS=0.; PB=0. - BTHS=0.; MSSSUM(:,:)=0. + BTH0=0.; BTH=0.; DDIAG=0.; SRHS=0.; PB=0. + MSSSUM(:,:)=0. #ifdef W3_T0 DOUT=0. #endif @@ -2047,50 +2183,33 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & ! 1. Initialization and numerical factors ! FACTURB=SSDSC(5)*USTAR**2/GRAV*DAIR/DWAT + DIKCUMUL = NINT(SSDSBRF1/(XFR-1.)) BREAKFRACTION=0. RENEWALFREQ=0. IK1=1 #ifdef W3_IG1 IK1=NINT(IGPARS(5))+1 #endif - NTHSUM=MIN(FLOOR(SSDSC(10)+0.5),NTH-1) ! number of angular bins for enhanced modulation - IF (NTHSUM.GT.0) THEN - WTHSUM(1:NTHSUM)=1 - WTHSUM(NTHSUM+1)=SSDSC(10)+0.5-NTHSUM - ELSE - WTHSUM(1)=2*SSDSC(10) - END IF ! - ! 1.b MSS parameters used for Modulation factors for B or lambda + ! 1.b MSS parameters used for Modulation factors for lambda (Romero ) ! IF (SSDSC(8).GT.0.OR.SSDSC(11).GT.0.OR.SSDSC(18).GT.0) THEN - MSSSUM2(:,:)=0. DO IK=1,NK - IMSSMAX (IK) = 1 MSSP = 0. MSSPC2 = 0. MSSPS2 = 0. MSSPCS = 0. ! - ! Sums the contributions to the directional MSS for all ITH + ! Sums the contributions to the directional MSS for all angles ! DO ITH=1,NTH IS=ITH+(IK-1)*NTH - MSSLONG(IK,ITH) = K(IK)**SSDSC(20) * A(IS) * DDEN(IK) / CG(IK) ! contribution to MSS - END DO - DO ITH=1,NTH - DO JTH=-NTHSUM,NTHSUM - ITH2 = 1+MOD(ITH-1+JTH+NTH,NTH) - MSSSUM2(IK,ITH) = MSSSUM2(IK,ITH)+MSSLONG(IK,ITH2)*WTHSUM(ABS(JTH)+1) - END DO - MSSPC2 = MSSPC2 +MSSLONG(IK,ITH)*EC2(ITH) - MSSPS2 = MSSPS2 +MSSLONG(IK,ITH)*ES2(ITH) - MSSPCS = MSSPCS +MSSLONG(IK,ITH)*ESC(ITH) - MSSP = MSSP +MSSLONG(IK,ITH) + MSSLONG = K(IK)**SSDSC(20) * A(IS) * DDEN(IK) / CG(IK) ! contribution to MSS + MSSPC2 = MSSPC2 +MSSLONG*EC2(ITH) + MSSPS2 = MSSPS2 +MSSLONG*ES2(ITH) + MSSPCS = MSSPCS +MSSLONG*ESC(ITH) + MSSP = MSSP +MSSLONG END DO - ! - ! Now sums over IK - ! MSSSUM (IK:NK,1) = MSSSUM (IK:NK,1) +MSSP MSSSUM (IK:NK,3) = MSSSUM (IK:NK,3) +MSSPC2 MSSSUM (IK:NK,4) = MSSSUM (IK:NK,4) +MSSPS2 @@ -2100,18 +2219,13 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & ! MSSD=0.5*(ATAN2(2*MSSSUM(IK,5),MSSSUM(IK,3)-MSSSUM(IK,4))) IF (MSSD.LT.0) MSSD = MSSD + PI - IMSSMAX (IK)=1+NINT(MSSD *NTH/TPI) - ! - ! mss along perpendicular direction - ! - MSSSUM (IK,2) = MAX(0.,MSSSUM(IK,4)*COS(MSSD)**2 & - -2*MSSSUM(IK,5)*SIN(MSSD)*COS(MSSD)+ & - MSSSUM(IK,3)*SIN(MSSD)**2 ) + MSSSUM (IK,2) = MSSD END DO END IF ! SSDSC(8).GT.0) THEN ! ! 2. Estimation of spontaneous breaking from local saturation ! + !############################################################################################" SELECT CASE (NINT(SSDSC(1))) CASE (1) ! @@ -2121,10 +2235,6 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & ! ! 2.a.1 Computes saturation ! - SDSNTH = MIN(NINT(SSDSDTH/(DTH*RADE)),NTH/2-1) - ! SSDSDIK is the integer difference in frequency bands - ! between the "large breakers" and short "wiped-out waves" - ! BTH(:) = 0. DO IK=IK1, NK @@ -2134,99 +2244,20 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & BTH(IS0+1)=0. ASUM = SUM(A(IS0+1:IS0+NTH)) BTH0(IK)=ASUM*FACSAT - IKC = MAX(1,IK-DIKCUMUL) - KLOC=K(IK)**(2-SSDSC(20)) ! local wavenumber factor, if mss not used. - + ! IF (SSDSDTH.GE.180) THEN ! integrates around full circle BTH(IS0+1:IS0+NTH)=BTH0(IK) ELSE DO ITH=1,NTH ! partial integration IS=ITH+(IK-1)*NTH - - ! straining effect of long waves on short waves - ! extended from Longuet-Higgins and Stewart (JFM 1960, eq. 2.27) the amplitude modulation - ! in deep water is equal to the long wave slope k*a cos(theta1-theta2) - ! Here we assume that the saturation is modulated as (1 + SSDSC(8) * sqrt(mss) ) - ! where mss_theta is the mss in direction ITH. - ! - ! Note: SSDSC(8) is sqrt(2)*times the mss MTF: equal to 4*sqrt(2) according to Longuet-Higgins and Stewart - ! - IF (SSDSC(8).GT.0.OR.SSDSC(11).GT.0) THEN - ! - MSSTH=(MSSSUM(IKC,1)-MSSSUM(IKC,2))*EC2(1+ABS(ITH-IMSSMAX (IKC))) & - +MSSSUM(IKC,2)*ES2(1+ABS(ITH-IMSSMAX (IKC)))*KLOC - ! - FACSTRAINB=1+SSDSC(8)*SQRT(MSSTH)+SSDSC(11)*SQRT(MSSSUM2(IKC,ITH)*KLOC) - ELSE - FACSTRAINB=1 - END IF - ! BTH(IS)=DOT_PRODUCT(SATWEIGHTS(:,ITH), A(IS0+SATINDICES(:,ITH)) ) & - *FACSAT*FACSTRAINB + *FACSAT END DO - IF (SSDSISO.NE.1) THEN - BTH0(IK)=MAXVAL(BTH(IS0+1:IS0+NTH)) - END IF + BTH0(IK)=MAXVAL(BTH(IS0+1:IS0+NTH)) END IF ! - END DO !NK END - ! - ! Optional smoothing of B and B0 over frequencies - ! - IF (SSDSBRFDF.GT.0.AND.SSDSBRFDF.LT.NK/2) THEN - BTH0S(:)=BTH0(:) - BTHS(:)=BTH(:) - NSMOOTH(:)=1 - DO IK=1, SSDSBRFDF - BTH0S(1+SSDSBRFDF)=BTH0S(1+SSDSBRFDF)+BTH0(IK) - NSMOOTH(1+SSDSBRFDF)=NSMOOTH(1+SSDSBRFDF)+1 - DO ITH=1,NTH - IS=ITH+(IK-1)*NTH - BTHS(ITH+SSDSBRFDF*NTH)=BTHS(ITH+SSDSBRFDF*NTH)+BTH(IS) - END DO - END DO - DO IK=IK1+1+SSDSBRFDF,1+2*SSDSBRFDF - BTH0S(1+SSDSBRFDF)=BTH0S(1+SSDSBRFDF)+BTH0(IK) - NSMOOTH(1+SSDSBRFDF)=NSMOOTH(1+SSDSBRFDF)+1 - DO ITH=1,NTH - IS=ITH+(IK-1)*NTH - BTHS(ITH+SSDSBRFDF*NTH)=BTHS(ITH+SSDSBRFDF*NTH)+BTH(IS) - END DO - END DO - DO IK=SSDSBRFDF,IK1,-1 - BTH0S(IK)=BTH0S(IK+1)-BTH0(IK+SSDSBRFDF+1) - NSMOOTH(IK)=NSMOOTH(IK+1)-1 - DO ITH=1,NTH - IS=ITH+(IK-1)*NTH - BTHS(IS)=BTHS(IS+NTH)-BTH(IS+(SSDSBRFDF+1)*NTH) - END DO - END DO - ! - DO IK=IK1+1+SSDSBRFDF,NK-SSDSBRFDF - BTH0S(IK)=BTH0S(IK-1)-BTH0(IK-SSDSBRFDF-1)+BTH0(IK+SSDSBRFDF) - NSMOOTH(IK)=NSMOOTH(IK-1) - DO ITH=1,NTH - IS=ITH+(IK-1)*NTH - BTHS(IS)=BTHS(IS-NTH)-BTH(IS-(SSDSBRFDF+1)*NTH)+BTH(IS+(SSDSBRFDF)*NTH) - END DO - END DO - ! - DO IK=NK-SSDSBRFDF+1,NK - BTH0S(IK)=BTH0S(IK-1)-BTH0(IK-SSDSBRFDF) - NSMOOTH(IK)=NSMOOTH(IK-1)-1 - DO ITH=1,NTH - IS=ITH+(IK-1)*NTH - BTHS(IS)=BTHS(IS-NTH)-BTH(IS-(SSDSBRFDF+1)*NTH) - END DO - END DO - ! division by NSMOOTH - BTH0(:)=MAX(0.,BTH0S(:)/NSMOOTH(:)) - DO IK=IK1,NK - IS0=(IK-1)*NTH - BTH(IS0+1:IS0+NTH)=MAX(0.,BTHS(IS0+1:IS0+NTH)/NSMOOTH(IK)) - END DO - END IF ! end of optional smoothing + END DO !IK=NK ! ! 2.a.2 Computes spontaneous breaking dissipation rate ! @@ -2238,7 +2269,8 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & MICHE=1. ELSE X=TANH(MIN(K(IK)*DEPTH,10.)) - MICHE=(X*(SSDSBM(1)+X*(SSDSBM(2)+X*(SSDSBM(3)+X*SSDSBM(4)))))**2 ! Correction of saturation level for shallow-water kinematics + ! Correction of saturation threshold for shallow-water kinematics + MICHE=(X*(SSDSBM(1)+X*(SSDSBM(2)+X*(SSDSBM(3)+X*SSDSBM(4)))))**2 END IF COEF1=(SSDSBR*MICHE) ! @@ -2268,7 +2300,7 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & BRLAMBDA = PB / (2.*PI**2.) SRHS = DDIAG * A - ! + !############################################################################################" CASE(2) ! ! 2.b Computes spontaneous breaking for T500 (Filipot et al. JGR 2010) @@ -2412,7 +2444,7 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & ! Compute Lambda = PB* l(k,th) ! with l(k,th)=1/(2*pi²)= the breaking crest density BRLAMBDA = PB / (2.*PI**2.) - ! + !############################################################################################" CASE(3) ! ! 2c Romero (GRL 2019) @@ -2425,27 +2457,15 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & KLOC=K(IK)**(2-SSDSC(20)) ! local wavenumber factor, if mss not used. BTH(1:NTH)=MAX(A(IS0+1:IS0+NTH)*SIG(IK)*K(IK)**3,.00000000000001) ! - IF (SSDSC(8).GT.0) THEN ! Applies modulation factor on B - DO ITH=1,NTH - MSSTH=(MSSSUM(IK,1)-MSSSUM(IK,2))*EC2(1+ABS(ITH-IMSSMAX (IK))) & - +MSSSUM(IK,2)*ES2(1+ABS(ITH-IMSSMAX (IK)))*KLOC - FACSTRAINB=(1.+SSDSC(8)*SQRT(MSSTH)+SSDSC(11)*SQRT(MSSSUM2(IK,ITH))*KLOC) - BTH(ITH)=BTH(ITH)*FACSTRAINB - END DO - END IF - ! + DIRFORCUM=DLWMEAN + IF (SSDSC(11).GT.0) DIRFORCUM=MSSSUM(IK,2) + C=SIG(IK)/K(IK) BTH0(IK)=sum(BTH(1:NTH)*DTH) IF (SSDSC(18).GT.0) THEN ! Applies modulation factor on Lambda DO ITH=1,NTH - IF (SSDSC(11).GT.0) THEN - MSSTH=(MSSSUM(IK,1)-MSSSUM(IK,2))*EC2(1+ABS(ITH-IMSSMAX (IK))) & - +MSSSUM(IK,2)*ES2(1+ABS(ITH-IMSSMAX (IK)))*KLOC - FACSTRAINL=1.+SSDSC(18)*SQRT(MSSTH)+SSDSC(11)*SQRT(MSSSUM2(IK,ITH)*KLOC) - ELSE - FACSTRAINL=1.+SSDSC(18)*((MSSSUM(IK,1)*KLOC)**SSDSC(14) * & ! Romero - (ECOS(ITH)*COS(DLWMEAN)+ESIN(ITH)*SIN(DLWMEAN))**2) - ENDIF + FACSTRAINL=1.+SSDSC(18)*((MSSSUM(IK,1)*KLOC)**SSDSC(14) * & ! Romero + (ECOS(ITH)*COS(DIRFORCUM)+ESIN(ITH)*SIN(DIRFORCUM))**2) LMODULATION(ITH)= FACSTRAINL**SSDSC(19) END DO ELSE @@ -2470,7 +2490,7 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & PB = BRLAMBDA *C ! END SELECT - ! + !############################################################################################" ! ! !/ ------------------------------------------------------------------- / @@ -2482,16 +2502,29 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & ! IF ( (SSDSC(3).NE.0.) .OR. (SSDSC(5).NE.0.) .OR. (SSDSC(21).NE.0.) ) THEN DO IK=IK1, NK + RENEWALFREQ = 0. FACTURB2=-2.*SIG(IK)*K(IK)*FACTURB DVISC=-4.*SSDSC(21)*K(IK)*K(IK) + C = SIG(IK)/K(IK) ! phase speed ! + IF (SSDSC(3).GT.0 .AND. IK.GT.DIKCUMUL) THEN + ! this is the cheap isotropic version + DO IK2=IK1,IK-DIKCUMUL + C2 = SIG(IK2)/K(IK2) + IS2=(IK2-1)*NTH + CUMULWISO=ABS(C2-C)*DSIP(IK2)/(0.5*C2) * DTH + RENEWALFREQ=RENEWALFREQ-CUMULWISO*SUM(BRLAMBDA(IS2+1:IS2+NTH)) + END DO + END IF + DO ITH=1,NTH IS=ITH+(IK-1)*NTH ! ! Computes cumulative effect from Breaking probability ! - RENEWALFREQ = 0. - IF (SSDSC(3).NE.0 .AND. IK.GT.DIKCUMUL) THEN + IF (SSDSC(3).LT.0 .AND. IK.GT.DIKCUMUL) THEN + RENEWALFREQ = 0. + ! this is the expensive and largely useless version DO IK2=IK1,IK-DIKCUMUL IF (BTH0(IK2).GT.SSDSBR) THEN IS2=(IK2-1)*NTH diff --git a/model/src/w3srcemd.F90 b/model/src/w3srcemd.F90 index 6aa708bb8..3500c0842 100644 --- a/model/src/w3srcemd.F90 +++ b/model/src/w3srcemd.F90 @@ -555,7 +555,7 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, & #endif #ifdef W3_ST4 USE W3SRC4MD, ONLY : W3SPR4, W3SIN4, W3SDS4 - USE W3GDATMD, ONLY : ZZWND, FFXFM, FFXPM, FFXFA + USE W3GDATMD, ONLY : ZZWND, FFXFM, FFXPM, FFXFA, SINTAILPAR #endif #ifdef W3_ST6 USE W3SRC6MD @@ -1033,10 +1033,14 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, & TWS = 1./FMEANWS #endif #ifdef W3_ST4 - TAUWX=0. - TAUWY=0. - IF ( IT .eq. 0 ) THEN + IF (SINTAILPAR(4).GT.0.5) THEN ! this is designed to keep the bug as an option + TAUWX=0. + TAUWY=0. + END IF + IF ( IT .EQ. 0 ) THEN LLWS(:) = .TRUE. + TAUWX=0. + TAUWY=0. USTAR=0. USTDIR=0. ELSE @@ -1060,7 +1064,7 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, & #endif #ifdef W3_ST4 - CALL W3SIN4 ( SPEC, CG1, WN2, U10ABS, USTAR, DRAT, AS, & + IF (SINTAILPAR(4).GT.0.5) CALL W3SIN4 ( SPEC, CG1, WN2, U10ABS, USTAR, DRAT, AS, & U10DIR, Z0, CD, TAUWX, TAUWY, TAUWAX, TAUWAY, & VSIN, VDIN, LLWS, IX, IY, BRLAMBDA ) END IF @@ -1902,6 +1906,13 @@ SUBROUTINE W3SRCE ( srce_call, IT, ISEA, JSEA, IX, IY, IMOD, & CALL W3SIN4 ( SPEC, CG1, WN2, U10ABS, USTAR, DRAT, AS, & U10DIR, Z0, CD, TAUWX, TAUWY, TAUWAX, TAUWAY, & VSIN, VDIN, LLWS, IX, IY, BRLAMBDA ) + IF (SINTAILPAR(4).LT.0.5) CALL W3SPR4 (SPEC, CG1, WN1, EMEAN, FMEAN, FMEAN1, WNMEAN,& + AMAX, U10ABS, U10DIR, & +#ifdef W3_FLX5 + TAUA, TAUADIR, DAIR, & +#endif + USTAR, USTDIR, & + TAUWX, TAUWY, CD, Z0, CHARN, LLWS, FMEANWS, DLWMEAN) #endif ! diff --git a/model/src/ww3_ounp.F90 b/model/src/ww3_ounp.F90 index 409888da1..dffc82ef2 100644 --- a/model/src/ww3_ounp.F90 +++ b/model/src/ww3_ounp.F90 @@ -2157,7 +2157,7 @@ SUBROUTINE W3EXNC(I,NCID,NREQ,INDREQ,ORDER) RHOAIR, USTAR, USTD, Z0, CD, CHARN ) #endif ! - DO ITT=1, 3 + DO ITT=1, 4 #ifdef W3_ST2 CALL W3SIN2 (A, CG, WN2, UABS, UDIRR, CD, Z0, & FPI, XIN, DIA ) diff --git a/regtests/ww3_ts1/input/namelists_ST4_T700.nml b/regtests/ww3_ts1/input/namelists_ST4_T700.nml index aa6ecdf70..b47cc70d0 100644 --- a/regtests/ww3_ts1/input/namelists_ST4_T700.nml +++ b/regtests/ww3_ts1/input/namelists_ST4_T700.nml @@ -1,4 +1,4 @@ &SDS4 SDSBCHOICE=3, SDSC2 = -3.8, SDSBR = 0.005, - SDSSTRAIN =0., SDSSTRAIN2 = 0., FXFM3 = 20., SDSFACMTF = 400., + FXFM3 = 20., SDSFACMTF = 400., SDSCUM=0., SDSC5 =0. / END OF NAMELISTS From 283b9895ea878f131f797d1964a3db1b328ba46f Mon Sep 17 00:00:00 2001 From: ma05cde Date: Tue, 18 Jul 2023 15:23:36 +0000 Subject: [PATCH 03/15] add regtest --- regtests/ww3_ts1/input/namelists_ST4_T475.nml | 7 + regtests/ww3_ts1/input/namelists_ST4_T702.nml | 12 + regtests/ww3_ts1/input/ww3_grid_ST4_T475.nml | 225 ++++++++++++++++++ regtests/ww3_ts1/input/ww3_grid_ST4_T702.nml | 225 ++++++++++++++++++ .../ww3_ts1/input_10ms/namelists_Romero.nml | 20 ++ .../ww3_ts1/input_10ms/namelists_ST4_T471.nml | 2 + .../ww3_ts1/input_10ms/namelists_ST4_T475.nml | 7 + .../ww3_ts1/input_10ms/namelists_ST4_T701.nml | 20 ++ .../ww3_ts1/input_10ms/namelists_ST4_T702.nml | 14 ++ regtests/ww3_ts1/input_10ms/points.list | 1 + regtests/ww3_ts1/input_10ms/switch | 1 + regtests/ww3_ts1/input_10ms/switch_ST4 | 1 + .../ww3_ts1/input_10ms/ww3_grid_Romero.nml | 225 ++++++++++++++++++ .../ww3_ts1/input_10ms/ww3_grid_ST4_T471.nml | 225 ++++++++++++++++++ .../ww3_ts1/input_10ms/ww3_grid_ST4_T475.nml | 225 ++++++++++++++++++ .../ww3_ts1/input_10ms/ww3_grid_ST4_T701.nml | 225 ++++++++++++++++++ .../ww3_ts1/input_10ms/ww3_grid_ST4_T702.nml | 225 ++++++++++++++++++ regtests/ww3_ts1/input_10ms/ww3_ounf.nml | 29 +++ regtests/ww3_ts1/input_10ms/ww3_ounp_spec.nml | 48 ++++ regtests/ww3_ts1/input_10ms/ww3_shel.nml | 54 +++++ 20 files changed, 1791 insertions(+) create mode 100644 regtests/ww3_ts1/input/namelists_ST4_T475.nml create mode 100644 regtests/ww3_ts1/input/namelists_ST4_T702.nml create mode 100644 regtests/ww3_ts1/input/ww3_grid_ST4_T475.nml create mode 100644 regtests/ww3_ts1/input/ww3_grid_ST4_T702.nml create mode 100644 regtests/ww3_ts1/input_10ms/namelists_Romero.nml create mode 100644 regtests/ww3_ts1/input_10ms/namelists_ST4_T471.nml create mode 100644 regtests/ww3_ts1/input_10ms/namelists_ST4_T475.nml create mode 100644 regtests/ww3_ts1/input_10ms/namelists_ST4_T701.nml create mode 100644 regtests/ww3_ts1/input_10ms/namelists_ST4_T702.nml create mode 100644 regtests/ww3_ts1/input_10ms/points.list create mode 100644 regtests/ww3_ts1/input_10ms/switch create mode 100644 regtests/ww3_ts1/input_10ms/switch_ST4 create mode 100644 regtests/ww3_ts1/input_10ms/ww3_grid_Romero.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T471.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T475.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T701.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T702.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_ounf.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_ounp_spec.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_shel.nml diff --git a/regtests/ww3_ts1/input/namelists_ST4_T475.nml b/regtests/ww3_ts1/input/namelists_ST4_T475.nml new file mode 100644 index 000000000..e104247aa --- /dev/null +++ b/regtests/ww3_ts1/input/namelists_ST4_T475.nml @@ -0,0 +1,7 @@ +&SIN4 BETAMAX = 1.75, SWELLF = 0.66, TAUWSHELTER = 0.3, + SWELLF3 = 0.022, SWELLF4 = 115000.0, SWELLF7 = 432000.00 / +&SDS4 FXFM3 = 2.5 / +&SIC2 IC2ROUGH = 0.001000, IC2VISC = 2.000, IC2DMAX =0.300 / +&SIS2 ISC1 =0.200E+00, IS2BREAK = T, IS2DUPDATE = F, IS2CREEPB = 0.200E+08 / + +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input/namelists_ST4_T702.nml b/regtests/ww3_ts1/input/namelists_ST4_T702.nml new file mode 100644 index 000000000..444c02e29 --- /dev/null +++ b/regtests/ww3_ts1/input/namelists_ST4_T702.nml @@ -0,0 +1,12 @@ +&SIN4 BETAMAX = 1.7, SWELLF = 0.60, TAUWSHELTER = 0.2, + SWELLF3 = 0.022, SWELLF4 = 115000.0, SWELLF7 = 432000.00 / +&SDS4 SDSBCHOICE = 3, SDSC2 = -3.80, FXFM3 = 20.00, + SDSBR = 0.005, SDSBT = 0.0011, SDSCUM = 0.300, SDSC5 = 1.0, + SDSMWD = 0.00, SDSFACMTF = 400 / +&SNL1 NLPROP = 25000000.0 / +&SIC2 IC2ROUGH = 0.001000, IC2VISC = 2.000, IC2DMAX =0.300 / +&SIS2 ISC1 =0.200E+00, IS2BREAK = T, IS2DUPDATE = F, IS2CREEPB = 0.200E+08 / +&MISC ICEHINIT = 0.5, ICEHMIN = 0.1, CICE0 = 0.25, NOSW =6, + CICEN = 2.00, LICE = 40000., FLAGTR = 4, FACBERG = 0.2 , + WCOR1=21., WCOR2=0.5 / +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input/ww3_grid_ST4_T475.nml b/regtests/ww3_ts1/input/ww3_grid_ST4_T475.nml new file mode 100644 index 000000000..7987e9528 --- /dev/null +++ b/regtests/ww3_ts1/input/ww3_grid_ST4_T475.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.0485 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 24 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 5s and 60s. +! A common value is 10s. +! DTMIN ~= 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input/namelists_ST4_T475.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input/ww3_grid_ST4_T702.nml b/regtests/ww3_ts1/input/ww3_grid_ST4_T702.nml new file mode 100644 index 000000000..48135e1d9 --- /dev/null +++ b/regtests/ww3_ts1/input/ww3_grid_ST4_T702.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.0485 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 24 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 5s and 60s. +! A common value is 10s. +! DTMIN ~= 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input/namelists_ST4_T702.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/namelists_Romero.nml b/regtests/ww3_ts1/input_10ms/namelists_Romero.nml new file mode 100644 index 000000000..919c786d9 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/namelists_Romero.nml @@ -0,0 +1,20 @@ +&SIN4 BETAMAX = 1.43, SWELLF = 0.66, TAUWSHELTER = 0.3, Z0MAX = 0.0008, + SWELLF3 = 0.022, SWELLF4 = 150000.0, SWELLF7 = 360000.00, ZALP = 0.006 / +&SDS4 SDSBCHOICE = 3, SDSC2 = -3.80, FXFM3 = 20.00, WNMEANP = 1.0 , + SDSSTRAINA = 0.00, SDSSTRAIN = 0.00, SDSSTRAIN2 = 0.00, + SDSBR = 0.005, SDSBT = 0.0011, SDSCUM = 0.300, SDSC5 = 1.0, + SDSMWD = 0.90, SDSFACMTF = 400 / +&SNL1 NLPROP = 25000000.0 / + +&OUTS P2SF = 1, E3D = 1, I1P2SF = 1, I2P2SF = 36 / +&PRO3 WDTHCG = 1.50, WDTHTH = 1.50 / +&REF1 REFCOAST = 0.1, REFCOSP_STRAIGHT = 4, REFFREQ = 0., REFICEBERG = 0.4, + REFMAP = 0., REFSLOPE = 0., REFSUBGRID = 0.2, REFRMAX = 0.5 / +&SIC2 IC2DISPER = F, IC2TURB = 0.5 , IC2ROUGH = 0.0001, + IC2REYNOLDS = 150000, IC2SMOOTH = 200000., IC2VISC = 1.0 / +&SIS2 ISC1 = 1., IS2C2 = 0.000000, IS2C3 = 0. , IS2BACKSCAT = 1. , + IS2BREAK = T, IS2DUPDATE = F , IS2CREEPB = 5E8 , IS2CREEPD = 0.3 / +&MISC ICEHINIT = 0.5, ICEHMIN = 0.1, CICE0 = 0.25, NOSW =6, + CICEN = 2.00, LICE = 40000., FACBERG = 0.2 , + WCOR1=21., WCOR2=0.5 / +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input_10ms/namelists_ST4_T471.nml b/regtests/ww3_ts1/input_10ms/namelists_ST4_T471.nml new file mode 100644 index 000000000..5ab7abb58 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/namelists_ST4_T471.nml @@ -0,0 +1,2 @@ +! T471 corresponds to the default parameter values for ST4. +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input_10ms/namelists_ST4_T475.nml b/regtests/ww3_ts1/input_10ms/namelists_ST4_T475.nml new file mode 100644 index 000000000..e104247aa --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/namelists_ST4_T475.nml @@ -0,0 +1,7 @@ +&SIN4 BETAMAX = 1.75, SWELLF = 0.66, TAUWSHELTER = 0.3, + SWELLF3 = 0.022, SWELLF4 = 115000.0, SWELLF7 = 432000.00 / +&SDS4 FXFM3 = 2.5 / +&SIC2 IC2ROUGH = 0.001000, IC2VISC = 2.000, IC2DMAX =0.300 / +&SIS2 ISC1 =0.200E+00, IS2BREAK = T, IS2DUPDATE = F, IS2CREEPB = 0.200E+08 / + +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input_10ms/namelists_ST4_T701.nml b/regtests/ww3_ts1/input_10ms/namelists_ST4_T701.nml new file mode 100644 index 000000000..b107fef73 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/namelists_ST4_T701.nml @@ -0,0 +1,20 @@ +&SIN4 BETAMAX = 1.7, SWELLF = 0.60, TAUWSHELTER = 0.3, + SWELLF3 = 0.022, SWELLF4 = 115000.0, SWELLF7 = 432000.00 / +&SDS4 SDSBCHOICE = 3, SDSC2 = -3.80, FXFM3 = 20.00, CUMSIGP = 2, + SDSBR = 0.005, SDSBT = 0.0011, SDSCUM = 0.300, SDSC5 = 1.0, + SDSMWD = 2.00, SDSFACMTF = 400 / +&SNL1 NLPROP = 25000000.0 / + +&OUTS P2SF = 1, E3D = 1, I1P2SF = 1, I2P2SF = 36 / +&PRO3 WDTHCG = 1.50, WDTHTH = 1.50 / +&REF1 REFCOAST = 0.1, REFCOSP_STRAIGHT = 4, REFFREQ = 0., REFICEBERG = 0.4, + REFMAP = 0., REFSLOPE = 0., REFSUBGRID = 0.2, REFRMAX = 0.5 / + +&SIC2 IC2ROUGH = 0.001000, IC2VISC = 2.000, IC2DMAX =0.300 / +&SIS2 ISC1 =0.200E+00, IS2BREAK = T, IS2DUPDATE = F, IS2CREEPB = 0.200E+08, + IS2CREEPD = 0.50 / + +&MISC ICEHINIT = 0.5, ICEHMIN = 0.1, CICE0 = 0.25, NOSW =6, + CICEN = 2.00, LICE = 40000., FACBERG = 0.2 , + WCOR1=21., WCOR2=0.5 / +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input_10ms/namelists_ST4_T702.nml b/regtests/ww3_ts1/input_10ms/namelists_ST4_T702.nml new file mode 100644 index 000000000..1b673567a --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/namelists_ST4_T702.nml @@ -0,0 +1,14 @@ +&SIN4 BETAMAX = 1.7, SWELLF = 0.60, TAUWSHELTER = 0.2, + SWELLF3 = 0.022, SWELLF4 = 115000.0, SWELLF7 = 432000.00 / +&SDS4 SDSBCHOICE = 3, SDSC2 = -3.80, FXFM3 = 20.00, CUMSIGP = 2, + SDSSTRAINA = 0.00, SDSSTRAIN = 0.00, SDSSTRAIN2 = 0.00, + SDSBR = 0.005, SDSBT = 0.0011, SDSCUM = 0.300, SDSC5 = 1.0, + SDSMWD = 0.00, SDSFACMTF = 400 / +&SNL1 NLPROP = 25000000.0 / +&SIC2 IC2ROUGH = 0.001000, IC2VISC = 2.000, IC2DMAX =0.300 / +&SIS2 ISC1 =0.200E+00, IS2BREAK = T, IS2DUPDATE = F, IS2CREEPB = 0.200E+08 / +! DO NOT FORGET TO ADD FLAGTR = 4 for real life runs ... +&MISC ICEHINIT = 0.5, ICEHMIN = 0.1, CICE0 = 0.25, NOSW =6, + CICEN = 2.00, LICE = 40000., FACBERG = 0.2 , + WCOR1=21., WCOR2=0.5 / +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input_10ms/points.list b/regtests/ww3_ts1/input_10ms/points.list new file mode 100644 index 000000000..5ad8fde50 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/points.list @@ -0,0 +1 @@ +0.0 0.0 'The_point' diff --git a/regtests/ww3_ts1/input_10ms/switch b/regtests/ww3_ts1/input_10ms/switch new file mode 100644 index 000000000..c3b8938ee --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/switch @@ -0,0 +1 @@ +NOGRB SHRD PR0 FLX0 LN1 ST4 NL1 BT1 DB1 TR0 BS0 IC0 IS0 REF0 WNT1 WNX1 CRT1 CRX1 O0 O1 O2 O3 O4 O5 O6 O7 O10 O11 diff --git a/regtests/ww3_ts1/input_10ms/switch_ST4 b/regtests/ww3_ts1/input_10ms/switch_ST4 new file mode 100644 index 000000000..c3b8938ee --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/switch_ST4 @@ -0,0 +1 @@ +NOGRB SHRD PR0 FLX0 LN1 ST4 NL1 BT1 DB1 TR0 BS0 IC0 IS0 REF0 WNT1 WNX1 CRT1 CRX1 O0 O1 O2 O3 O4 O5 O6 O7 O10 O11 diff --git a/regtests/ww3_ts1/input_10ms/ww3_grid_Romero.nml b/regtests/ww3_ts1/input_10ms/ww3_grid_Romero.nml new file mode 100644 index 000000000..c510784c1 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_grid_Romero.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.034 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 36 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 1s and 60s. +! A common value is 10s. +! DTMIN = 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input_10ms/namelists_Romero.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T471.nml b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T471.nml new file mode 100644 index 000000000..80069e4a1 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T471.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.034 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 36 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 1s and 60s. +! A common value is 10s. +! DTMIN = 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input_10ms/namelists_ST4_T471.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T475.nml b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T475.nml new file mode 100644 index 000000000..7d41e0b7d --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T475.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.0485 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 24 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 5s and 60s. +! A common value is 10s. +! DTMIN ~= 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input_10ms/namelists_ST4_T475.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T701.nml b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T701.nml new file mode 100644 index 000000000..d7aca045c --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T701.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.034 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 36 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 1s and 60s. +! A common value is 10s. +! DTMIN = 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input_10ms/namelists_ST4_T701.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T702.nml b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T702.nml new file mode 100644 index 000000000..7669b24c6 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T702.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.034 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 36 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 5s and 60s. +! A common value is 10s. +! DTMIN = 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input_10ms/namelists_ST4_T702.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/ww3_ounf.nml b/regtests/ww3_ts1/input_10ms/ww3_ounf.nml new file mode 100644 index 000000000..cb00b5dbf --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_ounf.nml @@ -0,0 +1,29 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III ww3_ounf.nml - Grid output post-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the output fields to postprocess via FIELD_NML namelist +! -------------------------------------------------------------------- ! +&FIELD_NML + FIELD%TIMESTART = '20000101 000000' + FIELD%TIMESTRIDE = '10' + FIELD%TIMECOUNT = '8000' + FIELD%LIST = 'DPT QP WND ICE HS MSS MSD FAW WCC WCF WCH WCM FOC TAW CHA UST' + FIELD%PARTITION = '0 1 2' + FIELD%TYPE = 4 +/ + +! -------------------------------------------------------------------- ! +! Define the content of the output file via FILE_NML namelist +! -------------------------------------------------------------------- ! +&FILE_NML + FILE%IX0 = 2 + FILE%IXN = 2 + FILE%IY0 = 2 + FILE%IYN = 2 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/ww3_ounp_spec.nml b/regtests/ww3_ts1/input_10ms/ww3_ounp_spec.nml new file mode 100644 index 000000000..34bac9764 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_ounp_spec.nml @@ -0,0 +1,48 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III ww3_ounp.nml - Point output post-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the output fields to postprocess via POINT_NML namelist +! -------------------------------------------------------------------- ! +&POINT_NML + POINT%TIMESTART = '20000101 000000' + POINT%TIMESTRIDE = '1800.' + POINT%TIMECOUNT = '1000' + POINT%TIMESPLIT = 4 + POINT%BUFFER = 100 + POINT%TYPE = 3 +/ + +! -------------------------------------------------------------------- ! +! Define the content of the output file via FILE_NML namelist +! -------------------------------------------------------------------- ! +&FILE_NML +/ + +! -------------------------------------------------------------------- ! +! Define the type 0, inventory of file +! -------------------------------------------------------------------- ! + + +! -------------------------------------------------------------------- ! +! Define the type 1, spectra via SPECTRA_NML namelist +! -------------------------------------------------------------------- ! +&SPECTRA_NML +/ + +! -------------------------------------------------------------------- ! +! Define the type 2, mean parameter via PARAM_NML namelist +! -------------------------------------------------------------------- ! +&PARAM_NML +/ + +! -------------------------------------------------------------------- ! +! Define the type 3, source terms via SOURCE_NML namelist +! -------------------------------------------------------------------- ! +&SOURCE_NML +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/ww3_shel.nml b/regtests/ww3_ts1/input_10ms/ww3_shel.nml new file mode 100644 index 000000000..81b6371d5 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_shel.nml @@ -0,0 +1,54 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III ww3_shel.nml - single-grid model ! +! -------------------------------------------------------------------- ! + + +! -------------------------------------------------------------------- ! +! Define top-level model parameters via DOMAIN_NML namelist +! -------------------------------------------------------------------- ! +&DOMAIN_NML + DOMAIN%START = '20000101 000000' + DOMAIN%STOP = '20000106 000000' +/ + +! -------------------------------------------------------------------- ! +! Define each forcing via the INPUT_NML namelist +! -------------------------------------------------------------------- ! +&INPUT_NML + INPUT%FORCING%WINDS = 'H' +/ + +! -------------------------------------------------------------------- ! +! Define the output types point parameters via OUTPUT_TYPE_NML namelist +! -------------------------------------------------------------------- ! +&OUTPUT_TYPE_NML + TYPE%FIELD%LIST = 'DPT QP WND ICE HS MSS MSD FAW WCC WCF WCH WCM FOC TAW CHA UST' + TYPE%POINT%FILE = '../input_10ms/points.list' +/ + +! -------------------------------------------------------------------- ! +! Define output dates via OUTPUT_DATE_NML namelist +! -------------------------------------------------------------------- ! +&OUTPUT_DATE_NML + DATE%FIELD = '19680606 000000' '1800' '20230618 000000' + DATE%POINT = '19680606 000000' '1800' '20230618 000000' +/ + +! -------------------------------------------------------------------- ! +! Define homogeneous input via HOMOG_COUNT_NML and HOMOG_INPUT_NML namelist +! -------------------------------------------------------------------- ! +&HOMOG_COUNT_NML + HOMOG_COUNT%N_CUR = 0 + HOMOG_COUNT%N_WND = 1 +/ + +&HOMOG_INPUT_NML + HOMOG_INPUT(1)%NAME = 'WND' + HOMOG_INPUT(1)%VALUE1 = 10. + HOMOG_INPUT(1)%VALUE2 = 270. + HOMOG_INPUT(1)%VALUE3 = 0. +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! From ba991c5c4338b67b9c64047bd59e1feded247172 Mon Sep 17 00:00:00 2001 From: ma05cde Date: Wed, 20 Sep 2023 10:43:10 +0000 Subject: [PATCH 04/15] add SINTAILPAR init --- model/src/w3gdatmd.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/model/src/w3gdatmd.F90 b/model/src/w3gdatmd.F90 index 92fc3e040..04c2afbe2 100644 --- a/model/src/w3gdatmd.F90 +++ b/model/src/w3gdatmd.F90 @@ -886,7 +886,7 @@ MODULE W3GDATMD REAL, POINTER :: DCKI(:,:), SATWEIGHTS(:,:),CUMULW(:,:),QBI(:,:) REAL :: AALPHA, BBETA, ZZ0MAX, ZZ0RAT, ZZALP,& SSINTHP, TTAUWSHELTER, SSWELLF(1:7), & - SSDSC(1:21), SSDSBR, & + SSDSC(1:21), SSDSBR, SINTAILPAR(1:5),& SSDSP, WWNMEANP, SSTXFTF, SSTXFTWN, & FFXPM, FFXFM, FFXFA, & SSDSBRF1, SSDSBRF2, SSDSBINT,SSDSBCK,& From fcf9fa63ae853c80d36ab140945e480ef99321e2 Mon Sep 17 00:00:00 2001 From: ma05cde Date: Wed, 27 Sep 2023 07:45:45 +0000 Subject: [PATCH 05/15] add namelists --- regtests/ww3_ts1/input/namelists_ST4_T707.nml | 13 + regtests/ww3_ts1/input/ww3_grid_ST4_T707.nml | 225 ++++++++++++++++++ .../ww3_ts1/input_10ms/namelists_ST4_T707.nml | 14 ++ .../ww3_ts1/input_10ms/namelists_ST4_T713.nml | 15 ++ .../ww3_ts1/input_10ms/ww3_grid_ST4_T707.nml | 225 ++++++++++++++++++ .../ww3_ts1/input_10ms/ww3_grid_ST4_T713.nml | 225 ++++++++++++++++++ 6 files changed, 717 insertions(+) create mode 100644 regtests/ww3_ts1/input/namelists_ST4_T707.nml create mode 100644 regtests/ww3_ts1/input/ww3_grid_ST4_T707.nml create mode 100644 regtests/ww3_ts1/input_10ms/namelists_ST4_T707.nml create mode 100644 regtests/ww3_ts1/input_10ms/namelists_ST4_T713.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T707.nml create mode 100644 regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T713.nml diff --git a/regtests/ww3_ts1/input/namelists_ST4_T707.nml b/regtests/ww3_ts1/input/namelists_ST4_T707.nml new file mode 100644 index 000000000..b1d0a9727 --- /dev/null +++ b/regtests/ww3_ts1/input/namelists_ST4_T707.nml @@ -0,0 +1,13 @@ + &SNL1 IQTYPE = -2, GQMNF1 = 11, GQMNT1 = 6, GQMNQ_OM2 = 6, + TAILNL=-5.0, GQMTHRSAT=5E-5, GQMTHRCOU = 0.05, GQAMP1=1., + GQAMP2=0.0022, GQAMP3=1., GQAMP4=1.0 / + &SIN4 BETAMAX = 1.6, TAUWSHELTER = 0.0 / + &SDS4 SDSBCHOICE=3, SDSC2 = -2.3, SDSBR = 0.005, CUMSIGP =2.0, + FXFM3 = 20, SDSFACMTF = 400., + SDSMWD = 2., SDSCUM = 0.35, SDSNUW =0, SDSC5=1., SDSBRF1=0.5 / +&SIC2 IC2ROUGH = 0.001000, IC2VISC = 2.000, IC2DMAX =0.300 / +&SIS2 ISC1 =0.200E+00, IS2BREAK = T, IS2DUPDATE = F, IS2CREEPB = 0.200E+08 / +&MISC ICEHINIT = 0.5, ICEHMIN = 0.1, CICE0 = 0.25, NOSW =6, + CICEN = 2.00, LICE = 40000., FLAGTR = 4, FACBERG = 0.2 , + WCOR1=21., WCOR2=0.5 / +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input/ww3_grid_ST4_T707.nml b/regtests/ww3_ts1/input/ww3_grid_ST4_T707.nml new file mode 100644 index 000000000..e6ef84a56 --- /dev/null +++ b/regtests/ww3_ts1/input/ww3_grid_ST4_T707.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.0485 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 24 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 5s and 60s. +! A common value is 10s. +! DTMIN ~= 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input/namelists_ST4_T707.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/namelists_ST4_T707.nml b/regtests/ww3_ts1/input_10ms/namelists_ST4_T707.nml new file mode 100644 index 000000000..8d19dd444 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/namelists_ST4_T707.nml @@ -0,0 +1,14 @@ + &SNL1 IQTYPE = -2, GQMNF1 = 11, GQMNT1 = 6, GQMNQ_OM2 = 6, + TAILNL=-5.0, GQMTHRSAT=5E-5, GQMTHRCOU = 0.05, GQAMP1=1., + GQAMP2=0.0022, GQAMP3=1., GQAMP4=1.0 / + &SIN4 BETAMAX = 1.6, TAUWSHELTER = 0.0 / + &SDS4 SDSBCHOICE=3, SDSC2 = -2.3, SDSBR = 0.005, CUMSIGP =2.0, + FXFM3 = 20, SDSFACMTF = 400., + SDSMWD = 2., SDSCUM = 0.35, SDSNUW =0, SDSC5=1., SDSBRF1=0.5 / +&SIC2 IC2ROUGH = 0.001000, IC2VISC = 2.000, IC2DMAX =0.300 / +&SIS2 ISC1 =0.200E+00, IS2BREAK = T, IS2DUPDATE = F, IS2CREEPB = 0.200E+08 / +! DO NOT FORGET TO ADD FLAGTR = 4 for real life runs ... +&MISC ICEHINIT = 0.5, ICEHMIN = 0.1, CICE0 = 0.25, NOSW =6, + CICEN = 2.00, LICE = 40000., FACBERG = 0.2 , + WCOR1=21., WCOR2=0.5 / +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input_10ms/namelists_ST4_T713.nml b/regtests/ww3_ts1/input_10ms/namelists_ST4_T713.nml new file mode 100644 index 000000000..fa4a7eb78 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/namelists_ST4_T713.nml @@ -0,0 +1,15 @@ + &SNL1 IQTYPE = -2, GQMNF1 = 11, GQMNT1 = 6, GQMNQ_OM2 = 6, + TAILNL=-5.0, GQMTHRSAT=5E-5, GQMTHRCOU = 0.05, GQAMP1=1., + GQAMP2=0.0022, GQAMP3=2. / +&SIN4 BETAMAX = 1.1, TAUWSHELTER = 0.0, TAUWBUG = 0, + VISCSTRESS =1., SINTABLE=0 / +&SDS4 SDSBCHOICE=3, SDSC2 = -2.5, SDSBR = 0.005, CUMSIGP =2.0, + SDSSTRAIN2 =1.,SDSCUMP=1., FXFM3 = 20, SDSFACMTF = 200., + SDSMWD = 0.9, SDSCUM = 0.3, SDSNUW =0, SDSC5=0.5, SDSBRF1=0.5 / +&SIC2 IC2ROUGH = 0.001000, IC2VISC = 2.000, IC2DMAX =0.300 / +&SIS2 ISC1 =0.200E+00, IS2BREAK = T, IS2DUPDATE = F, IS2CREEPB = 0.200E+08 / +! DO NOT FORGET TO ADD FLAGTR = 4 for real life runs ... +&MISC ICEHINIT = 0.5, ICEHMIN = 0.1, CICE0 = 0.25, NOSW =6, + CICEN = 2.00, LICE = 40000., FACBERG = 0.2 , + WCOR1=21., WCOR2=0.5 / +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T707.nml b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T707.nml new file mode 100644 index 000000000..5378ebec3 --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T707.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.034 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 36 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 5s and 60s. +! A common value is 10s. +! DTMIN = 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 5. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input_10ms/namelists_ST4_T707.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! diff --git a/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T713.nml b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T713.nml new file mode 100644 index 000000000..3efd65adf --- /dev/null +++ b/regtests/ww3_ts1/input_10ms/ww3_grid_ST4_T713.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.034 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 36 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 5s and 60s. +! A common value is 10s. +! DTMIN = 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 5. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input_10ms/namelists_ST4_T713.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! From f2a218b87999a6966fb82d88fc0497862f5f340c Mon Sep 17 00:00:00 2001 From: ma05cde Date: Tue, 3 Oct 2023 15:07:31 +0000 Subject: [PATCH 06/15] address issue #1085 address issue #1074 add regtest --- model/src/w3src4md.F90 | 3 +- regtests/bin/matrix.base | 5 + regtests/ww3_ts1/input/namelists_ST4_T500.nml | 3 + regtests/ww3_ts1/input/ww3_grid_ST4_T500.nml | 225 ++++++++++++++++++ 4 files changed, 235 insertions(+), 1 deletion(-) create mode 100644 regtests/ww3_ts1/input/namelists_ST4_T500.nml create mode 100644 regtests/ww3_ts1/input/ww3_grid_ST4_T500.nml diff --git a/model/src/w3src4md.F90 b/model/src/w3src4md.F90 index 2a7ebd1a9..07142f7ac 100644 --- a/model/src/w3src4md.F90 +++ b/model/src/w3src4md.F90 @@ -2444,6 +2444,7 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & ! Compute Lambda = PB* l(k,th) ! with l(k,th)=1/(2*pi²)= the breaking crest density BRLAMBDA = PB / (2.*PI**2.) + SRHS = DDIAG * A !############################################################################################" CASE(3) ! @@ -2553,7 +2554,7 @@ SUBROUTINE W3SDS4 (A, K, CG, USTAR, USDIR, DEPTH, DAIR, SRHS, & RETURN END IF ! - WHITECAP(1:2) = 0. + WHITECAP(1:4) = 0. ! ! precomputes integration of Lambda over direction ! times wavelength times a (a=5 in Reul&Chapron JGR 2003) times dk diff --git a/regtests/bin/matrix.base b/regtests/bin/matrix.base index 88e7ee835..cac401bff 100755 --- a/regtests/bin/matrix.base +++ b/regtests/bin/matrix.base @@ -905,12 +905,17 @@ echo "$rtst -s ST2 -w work_ST2 $ww3 ww3_ts1" >> matrix.body echo "$rtst -s ST3 -w work_ST3 $ww3 ww3_ts1" >> matrix.body echo "$rtst -s ST4 -w work_ST4 $ww3 ww3_ts1" >> matrix.body + echo "$rtst -s ST4 -w work_ST4_T500 -g ST4_T500 -N $ww3 ww3_ts1" >> matrix.body echo "$rtst -s ST4 -w work_ST4_T700 -g ST4_T700 -N $ww3 ww3_ts1" >> matrix.body echo "$rtst -s ST4_WRT -w work_ST4_WRT $ww3 ww3_ts1" >> matrix.body echo "$rtst -s ST4_GMD -w work_ST4_GMD $ww3 ww3_ts1" >> matrix.body echo "$rtst -s ST4_TSA -w work_ST4_TSA $ww3 ww3_ts1" >> matrix.body echo "$rtst -s ST6 -w work_ST6 $ww3 ww3_ts1" >> matrix.body echo "$rtst -w work_NL5 -i input_nl5_matrix $ww3 ww3_ts1" >> matrix.body + echo "$rtst -g Romero -w work_Romero -i input_10ms -N $ww3 ww3_ts1" >> matrix.body + echo "$rtst -g ST4_T701 -w work_T701 -i input_10ms -N $ww3 ww3_ts1" >> matrix.body + echo "$rtst -g ST4_T702 -w work_T702 -i input_10ms -N $ww3 ww3_ts1" >> matrix.body + fi # fetch limited growth, no switch sharing here diff --git a/regtests/ww3_ts1/input/namelists_ST4_T500.nml b/regtests/ww3_ts1/input/namelists_ST4_T500.nml new file mode 100644 index 000000000..317705f93 --- /dev/null +++ b/regtests/ww3_ts1/input/namelists_ST4_T500.nml @@ -0,0 +1,3 @@ + &SDS4 SDSBCHOICE=2, SDSC2 = 0.0, SDSBR = 0.005, + FXFM3 = 9., SDSBCK = 0.185, SDSHCK = 1.5/ +END OF NAMELISTS diff --git a/regtests/ww3_ts1/input/ww3_grid_ST4_T500.nml b/regtests/ww3_ts1/input/ww3_grid_ST4_T500.nml new file mode 100644 index 000000000..bef18d975 --- /dev/null +++ b/regtests/ww3_ts1/input/ww3_grid_ST4_T500.nml @@ -0,0 +1,225 @@ +! -------------------------------------------------------------------- ! +! WAVEWATCH III - ww3_grid.nml - Grid pre-processing ! +! -------------------------------------------------------------------- ! + +! -------------------------------------------------------------------- ! +! Define the spectrum parameterization via SPECTRUM_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! SPECTRUM%XFR = 0. ! frequency increment +! SPECTRUM%FREQ1 = 0. ! first frequency (Hz) +! SPECTRUM%NK = 0 ! number of frequencies (wavenumbers) +! SPECTRUM%NTH = 0 ! number of direction bins +! SPECTRUM%THOFF = 0. ! relative offset of first direction [-0.5,0.5] +! -------------------------------------------------------------------- ! +&SPECTRUM_NML + SPECTRUM%XFR = 1.10 + SPECTRUM%FREQ1 = 0.0485 + SPECTRUM%NK = 36 + SPECTRUM%NTH = 24 +/ + +! -------------------------------------------------------------------- ! +! Define the run parameterization via RUN_NML namelist +! +! * namelist must be terminated with / +! * definitions & defaults: +! RUN%FLDRY = F ! dry run (I/O only, no calculation) +! RUN%FLCX = F ! x-component of propagation +! RUN%FLCY = F ! y-component of propagation +! RUN%FLCTH = F ! direction shift +! RUN%FLCK = F ! wavenumber shift +! RUN%FLSOU = F ! source terms +! -------------------------------------------------------------------- ! +&RUN_NML + RUN%FLSOU = T +/ + +! -------------------------------------------------------------------- ! +! Define the timesteps parameterization via TIMESTEPS_NML namelist +! +! * It is highly recommended to set up time steps which are multiple +! between them. +! +! * The first time step to calculate is the maximum CFL time step +! which depend on the lowest frequency FREQ1 previously set up and the +! lowest spatial grid resolution in meters DXY. +! reminder : 1 degree=60minutes // 1minute=1mile // 1mile=1.852km +! The formula for the CFL time is : +! Tcfl = DXY / (G / (FREQ1*4*Pi) ) with the constants Pi=3,14 and G=9.8m/s²; +! DTXY ~= 90% Tcfl +! DTMAX ~= 3 * DTXY (maximum global time step limit) +! +! * The refraction time step depends on how strong can be the current velocities +! on your grid : +! DTKTH ~= DTMAX / 2 ! in case of no or light current velocities +! DTKTH ~= DTMAX / 10 ! in case of strong current velocities +! +! * The source terms time step is usually defined between 5s and 60s. +! A common value is 10s. +! DTMIN ~= 10 +! +! * namelist must be terminated with / +! * definitions & defaults: +! TIMESTEPS%DTMAX = 0. ! maximum global time step (s) +! TIMESTEPS%DTXY = 0. ! maximum CFL time step for x-y (s) +! TIMESTEPS%DTKTH = 0. ! maximum CFL time step for k-th (s) +! TIMESTEPS%DTMIN = 0. ! minimum source term time step (s) +! -------------------------------------------------------------------- ! +&TIMESTEPS_NML + TIMESTEPS%DTMAX = 900. + TIMESTEPS%DTXY = 900. + TIMESTEPS%DTKTH = 900. + TIMESTEPS%DTMIN = 15. +/ + +! -------------------------------------------------------------------- ! +! Define the grid to preprocess via GRID_NML namelist +! +! * the tunable parameters for source terms, propagation schemes, and +! numerics are read using namelists. +! * Any namelist found in the folowing sections is temporarily written +! to param.scratch, and read from there if necessary. +! * The order of the namelists is immaterial. +! * Namelists not needed for the given switch settings will be skipped +! automatically +! +! * grid type can be : +! 'RECT' : rectilinear +! 'CURV' : curvilinear +! 'UNST' : unstructured (triangle-based) +! +! * coordinate system can be : +! 'SPHE' : Spherical (degrees) +! 'CART' : Cartesian (meters) +! +! * grid closure can only be applied in spherical coordinates +! +! * grid closure can be : +! 'NONE' : No closure is applied +! 'SMPL' : Simple grid closure. Grid is periodic in the +! : i-index and wraps at i=NX+1. In other words, +! : (NX+1,J) => (1,J). A grid with simple closure +! : may be rectilinear or curvilinear. +! 'TRPL' : Tripole grid closure : Grid is periodic in the +! : i-index and wraps at i=NX+1 and has closure at +! : j=NY+1. In other words, (NX+1,J<=NY) => (1,J) +! : and (I,NY+1) => (NX-I+1,NY). Tripole +! : grid closure requires that NX be even. A grid +! : with tripole closure must be curvilinear. +! +! * The coastline limit depth is the value which distinguish the sea +! points to the land points. All the points with depth values (ZBIN) +! greater than this limit (ZLIM) will be considered as excluded points +! and will never be wet points, even if the water level grows over. +! It can only overwrite the status of a sea point to a land point. +! The value must have a negative value under the mean sea level +! +! * The minimum water depth allowed to compute the model is the absolute +! depth value (DMIN) used in the model if the input depth is lower to +! avoid the model to blow up. +! +! * namelist must be terminated with / +! * definitions & defaults: +! GRID%NAME = 'unset' ! grid name (30 char) +! GRID%NML = 'namelists.nml' ! namelists filename +! GRID%TYPE = 'unset' ! grid type +! GRID%COORD = 'unset' ! coordinate system +! GRID%CLOS = 'unset' ! grid closure +! +! GRID%ZLIM = 0. ! coastline limit depth (m) +! GRID%DMIN = 0. ! abs. minimum water depth (m) +! -------------------------------------------------------------------- ! +&GRID_NML + GRID%NAME = 'HOMOGENEOUS SOURCE TERM TEST' + GRID%NML = '../input/namelists_ST4_T500.nml' + GRID%TYPE = 'RECT' + GRID%COORD = 'SPHE' + GRID%CLOS = 'NONE' + GRID%ZLIM = -5. + GRID%DMIN = 5.75 +/ + +! -------------------------------------------------------------------- ! +! Define the rectilinear grid type via RECT_NML namelist +! - only for RECT grids - +! +! * The minimum grid size is 3x3. +! +! * If the grid increments SX and SY are given in minutes of arc, the scaling +! factor SF must be set to 60. to provide an increment factor in degree. +! +! * If CSTRG='SMPL', then SX is forced to 360/NX. +! +! * value <= value_read / scale_fac +! +! * namelist must be terminated with / +! * definitions & defaults: +! RECT%NX = 0 ! number of points along x-axis +! RECT%NY = 0 ! number of points along y-axis +! +! RECT%SX = 0. ! grid increment along x-axis +! RECT%SY = 0. ! grid increment along y-axis +! RECT%SF = 1. ! scaling division factor for x-y axis +! +! RECT%X0 = 0. ! x-coordinate of lower-left corner (deg) +! RECT%Y0 = 0. ! y-coordinate of lower-left corner (deg) +! RECT%SF0 = 1. ! scaling division factor for x0,y0 coord +! -------------------------------------------------------------------- ! +&RECT_NML + RECT%NX = 3 + RECT%NY = 3 + RECT%SX = 1. + RECT%SY = 1. + RECT%SF = 1.E-2 + RECT%X0 = -1. + RECT%Y0 = -1. + RECT%SF0 = 1.E-2 +/ + +! -------------------------------------------------------------------- ! +! Define the depth to preprocess via DEPTH_NML namelist +! - for RECT and CURV grids - +! +! * if no obstruction subgrid, need to set &MISC FLAGTR = 0 +! +! * The depth value must have negative values under the mean sea level +! +! * value <= value_read * scale_fac +! +! * IDLA : Layout indicator : +! 1 : Read line-by-line bottom to top. (default) +! 2 : Like 1, single read statement. +! 3 : Read line-by-line top to bottom. +! 4 : Like 3, single read statement. +! * IDFM : format indicator : +! 1 : Free format. (default) +! 2 : Fixed format. +! 3 : Unformatted. +! * FORMAT : element format to read : +! '(....)' : auto detected (default) +! '(f10.6)' : float type +! +! * Example : +! IDF SF IDLA IDFM FORMAT FILENAME +! 50 0.001 1 1 '(....)' 'GLOB-30M.bot' +! +! * namelist must be terminated with / +! * definitions & defaults: +! DEPTH%SF = 1. ! scale factor +! DEPTH%FILENAME = 'unset' ! filename +! DEPTH%IDF = 50 ! file unit number +! DEPTH%IDLA = 1 ! layout indicator +! DEPTH%IDFM = 1 ! format indicator +! DEPTH%FORMAT = '(....)' ! formatted read format +! -------------------------------------------------------------------- ! +&DEPTH_NML + DEPTH%SF = -2500. + DEPTH%FILENAME = '../input/HOMOGENEOUS.depth' + DEPTH%IDLA = 3 +/ + +! -------------------------------------------------------------------- ! +! WAVEWATCH III - end of namelist ! +! -------------------------------------------------------------------- ! From 776aba6137b14239c40a917a2ac5139e96146eaf Mon Sep 17 00:00:00 2001 From: ma05cde Date: Tue, 7 Nov 2023 10:21:03 +0000 Subject: [PATCH 07/15] correct typo --- model/src/w3iogrmd.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/model/src/w3iogrmd.F90 b/model/src/w3iogrmd.F90 index 9c01c0949..0f107d3b2 100644 --- a/model/src/w3iogrmd.F90 +++ b/model/src/w3iogrmd.F90 @@ -1560,7 +1560,7 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & SSDSHCK, DELUST, DELTAIL, DELTAUW, & DELU, DELALP, TAUT, TAUHFT, TAUHFT2, & IKTAB, DCKI, QBI, SATINDICES, SATWEIGHTS, & - DIKCUMUL, CUMULW, SINTALPAR + DIKCUMUL, CUMULW, SINTAILPAR IF (SINTAILPAR(1).GT.0.5) THEN CALL INSIN4(.FALSE.) READ (NDSM) DELUST, DELTAIL, DELTAUW, DELU, DELALP, & From e407f715334277983257563d9956def3c16ca18b Mon Sep 17 00:00:00 2001 From: ma05cde Date: Fri, 10 Nov 2023 09:28:16 +0000 Subject: [PATCH 08/15] remove comment --- model/src/w3gdatmd.F90 | 5 +---- 1 file changed, 1 insertion(+), 4 deletions(-) diff --git a/model/src/w3gdatmd.F90 b/model/src/w3gdatmd.F90 index 44c96ff11..7a960ccd6 100644 --- a/model/src/w3gdatmd.F90 +++ b/model/src/w3gdatmd.F90 @@ -2077,10 +2077,7 @@ SUBROUTINE W3DIMS ( IMOD, MK, MTH, NDSE, NDST ) SDSNTH = MTH/2-1 !MIN(NINT(SSDSDTH/(DTH*RADE)),MTH/2-1) ALLOCATE( MPARS(IMOD)%SRCPS%SATINDICES(2*SDSNTH+1,MTH), & MPARS(IMOD)%SRCPS%SATWEIGHTS(2*SDSNTH+1,MTH), & - STAT=ISTAT ) - CHECK_ALLOC_STATUS ( ISTAT ) - ! IF (SSDSC(3).LT.0.) & - ALLOCATE (MPARS(IMOD)%SRCPS%CUMULW(MSPEC,MSPEC), & + MPARS(IMOD)%SRCPS%CUMULW(MSPEC,MSPEC), & STAT=ISTAT ) CHECK_ALLOC_STATUS ( ISTAT ) #endif From e1929a6b3731194f849b92061225401d3cb934da Mon Sep 17 00:00:00 2001 From: ma05cde Date: Fri, 10 Nov 2023 09:29:31 +0000 Subject: [PATCH 09/15] correct typo --- regtests/bin/run_cmake_test | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/regtests/bin/run_cmake_test b/regtests/bin/run_cmake_test index 86248bb4e..844f3e23e 100755 --- a/regtests/bin/run_cmake_test +++ b/regtests/bin/run_cmake_test @@ -109,7 +109,7 @@ EOF # --------------------------------------------------------------------------- # echo ' ' -echo " Running now options: run_test $*" +echo " Running now options: run_cmake_test $*" echo ' ' # 2.a Setup array of command-line arguments @@ -377,7 +377,7 @@ fi if [ $time_count ] then # Add time counter if -T - echo " REGTESTS Time counter: run_test $ARGS" >> time_count.txt + echo " REGTESTS Time counter: run_cmake_test $ARGS" >> time_count.txt Tstart=`date +"%s.%2N"` fi From 7c375bda0f3b6545292137050a1a4a792249d3e6 Mon Sep 17 00:00:00 2001 From: ma05cde Date: Fri, 10 Nov 2023 09:35:38 +0000 Subject: [PATCH 10/15] correct bug in mod_def reading --- model/src/w3iogrmd.F90 | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/model/src/w3iogrmd.F90 b/model/src/w3iogrmd.F90 index 0f107d3b2..c2ef21a7f 100644 --- a/model/src/w3iogrmd.F90 +++ b/model/src/w3iogrmd.F90 @@ -1563,8 +1563,9 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & DIKCUMUL, CUMULW, SINTAILPAR IF (SINTAILPAR(1).GT.0.5) THEN CALL INSIN4(.FALSE.) - READ (NDSM) DELUST, DELTAIL, DELTAUW, DELU, DELALP, & - TAUT, TAUHFT, TAUHFT2 + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & + DELUST, DELTAIL, DELTAUW, DELU, DELALP, & + TAUT, TAUHFT, TAUHFT2 END IF END IF #endif From 3608127770ecdcd2658e01ea5fec87ec607882ad Mon Sep 17 00:00:00 2001 From: ma05cde Date: Fri, 10 Nov 2023 14:27:26 +0000 Subject: [PATCH 11/15] add variable initialization correct duplicated variable read/write --- model/src/w3gdatmd.F90 | 6 ++++++ model/src/w3iogrmd.F90 | 34 +++++++++++++++------------------- 2 files changed, 21 insertions(+), 19 deletions(-) diff --git a/model/src/w3gdatmd.F90 b/model/src/w3gdatmd.F90 index 7a960ccd6..59d3bcddf 100644 --- a/model/src/w3gdatmd.F90 +++ b/model/src/w3gdatmd.F90 @@ -2074,12 +2074,18 @@ SUBROUTINE W3DIMS ( IMOD, MK, MTH, NDSE, NDST ) MPARS(IMOD)%SRCPS%QBI(NKHS,NKD), & STAT=ISTAT ) CHECK_ALLOC_STATUS ( ISTAT ) + MPARS(IMOD)%SRCPS%IKTAB(:,:)=0. + MPARS(IMOD)%SRCPS%DCKI(:,:)=0. + MPARS(IMOD)%SRCPS%QBI(:,:)=0. SDSNTH = MTH/2-1 !MIN(NINT(SSDSDTH/(DTH*RADE)),MTH/2-1) ALLOCATE( MPARS(IMOD)%SRCPS%SATINDICES(2*SDSNTH+1,MTH), & MPARS(IMOD)%SRCPS%SATWEIGHTS(2*SDSNTH+1,MTH), & MPARS(IMOD)%SRCPS%CUMULW(MSPEC,MSPEC), & STAT=ISTAT ) CHECK_ALLOC_STATUS ( ISTAT ) + MPARS(IMOD)%SRCPS%SATINDICES(:,:)=0. + MPARS(IMOD)%SRCPS%SATWEIGHTS(:,:)=0. + MPARS(IMOD)%SRCPS%CUMULW(:,:)=0. #endif ! SGRDS(IMOD)%SINIT = .TRUE. diff --git a/model/src/w3iogrmd.F90 b/model/src/w3iogrmd.F90 index c2ef21a7f..fadd237f6 100644 --- a/model/src/w3iogrmd.F90 +++ b/model/src/w3iogrmd.F90 @@ -1065,12 +1065,12 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & READ(NDSM,END=801,ERR=802,IOSTAT=IERR)GRIDSHIFT #ifdef W3_SEC1 - READ (NDSM) NITERSEC1 + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) NITERSEC1 #endif ! #ifdef W3_RTD !! Read rotated Polat/lon and AnglD from mod_def JGLi12Jun2012 - READ (NDSM) PoLat, PoLon, AnglD, FLAGUNR + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) PoLat, PoLon, AnglD, FLAGUNR #endif ! @@ -1313,35 +1313,35 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & FACP, XREL, XFLT, FXFM, FXPM, XFT, XFC, FACSD, FHMAX, & FFACBERG, DELAB, FWTABLE #ifdef W3_RWND - READ (NDSM) & + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & RWINDC #endif #ifdef W3_WCOR - READ (NDSM) & + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & WWCOR #endif #ifdef W3_REF1 - READ (NDSM) & + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & RREF, REFPARS, REFLC, REFLD #endif #ifdef W3_IG1 - READ (NDSM) & + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & IGPARS(1:12) #endif #ifdef W3_IC2 - READ (NDSM) & + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & IC2PARS(1:8) #endif #ifdef W3_IC3 - READ (NDSM) & + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & IC3PARS #endif #ifdef W3_IC4 - READ (NDSM) & + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & IC4PARS,IC4_KI,IC4_FC #endif #ifdef W3_IC5 - READ (NDSM) & + READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & IC5PARS #endif END IF @@ -1506,8 +1506,7 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & SSTXFTFTAIL, SSTXFTWN, SSTXFTF, SSTXFTWN, & SSDSBRF1, SSDSBRF2, SSDSBRFDF,SSDSBCK, SSDSABK, & SSDSPBK, SSDSBINT, FFXPM, FFXFM, FFXFA, & - SSDSHCK, DELUST, DELTAIL, DELTAUW, & - DELU, DELALP, TAUT, TAUHFT, TAUHFT2, & + SSDSHCK, & IKTAB, DCKI, QBI, SATINDICES, SATWEIGHTS, & DIKCUMUL, CUMULW, SINTAILPAR #ifdef W3_ASCII @@ -1520,10 +1519,9 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & SSTXFTFTAIL, SSTXFTWN, SSTXFTF, SSTXFTWN, & SSDSBRF1, SSDSBRF2, SSDSBRFDF,SSDSBCK, SSDSABK, & SSDSPBK, SSDSBINT, FFXPM, FFXFM, FFXFA, & - SSDSHCK, DELUST, DELTAIL, DELTAUW, & - DELU, DELALP, TAUT, TAUHFT, TAUHFT2, & + SSDSHCK, & IKTAB, DCKI, QBI, SATINDICES, SATWEIGHTS, & - DIKCUMUL, CUMULW:', & + DIKCUMUL, CUMULW, SINTAILPAR:', & ZZWND, AALPHA, ZZ0MAX, BBETA, SSINTHP, ZZALP, & TTAUWSHELTER, SSWELLFPAR, SSWELLF, SSINBR, & ZZ0RAT, SSDSC, & @@ -1532,8 +1530,7 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & SSTXFTFTAIL, SSTXFTWN, SSTXFTF, SSTXFTWN, & SSDSBRF1, SSDSBRF2, SSDSBRFDF,SSDSBCK, SSDSABK, & SSDSPBK, SSDSBINT, FFXPM, FFXFM, FFXFA, & - SSDSHCK, DELUST, DELTAIL, DELTAUW, & - DELU, DELALP, TAUT, TAUHFT, TAUHFT2, & + SSDSHCK, & IKTAB, DCKI, QBI, SATINDICES, SATWEIGHTS, & DIKCUMUL, CUMULW, SINTAILPAR #endif @@ -1557,8 +1554,7 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & SSTXFTFTAIL, SSTXFTWN, SSTXFTF, SSTXFTWN, & SSDSBRF1, SSDSBRF2, SSDSBRFDF,SSDSBCK, SSDSABK, & SSDSPBK, SSDSBINT, FFXPM, FFXFM, FFXFA, & - SSDSHCK, DELUST, DELTAIL, DELTAUW, & - DELU, DELALP, TAUT, TAUHFT, TAUHFT2, & + SSDSHCK, & IKTAB, DCKI, QBI, SATINDICES, SATWEIGHTS, & DIKCUMUL, CUMULW, SINTAILPAR IF (SINTAILPAR(1).GT.0.5) THEN From bd3aa722f55e465e53e7822e672e9e1896190892 Mon Sep 17 00:00:00 2001 From: ma05cde Date: Tue, 21 Nov 2023 15:59:13 +0000 Subject: [PATCH 12/15] manually merge bgufix/gqm --- manual/eqs/NL1.tex | 29 +++++++++++++++-------------- model/src/w3gridmd.F90 | 4 ++-- model/src/w3iogrmd.F90 | 13 ++++++++----- model/src/w3snl1md.F90 | 22 ++++++++++++++-------- 4 files changed, 39 insertions(+), 29 deletions(-) diff --git a/manual/eqs/NL1.tex b/manual/eqs/NL1.tex index d9bc4c521..c45e9294d 100644 --- a/manual/eqs/NL1.tex +++ b/manual/eqs/NL1.tex @@ -55,37 +55,38 @@ \subsubsection{~$S_{nl}$: Discrete Interaction Approximation (\dia)} \label{sec: \sin(\delta_{\theta,3})&=&\sin(\delta_{\theta,2}) (1-\lambda)^2/(1+\lambda)^2. \end{eqnarray} - For these quadruplets, each source term value -$S_{nl}(\bk)$ corresponding to each discrete $(f_r,\theta)$ -we compute the three contributions that correspond to the situation in which $\bk$ takes the role of $\bk$,$\bk_{2,+}$, $\bk_{2,-}$, $\bk_{3,+}$ and $\bk_{3,-}$ in the quadruplet, namely the full source term is +Hence for any $\bk$ one quadruplet selects $\bk_{2,+}$ and $\bk_{3,+}$, and the other quadruplet selects its mirror image +$\bk_{2,-}$, $\bk_{2,-}$. Because there are 3 different components interacting in the two DIA-selected quadruplets, any discrete spectral component $(f_r,\theta)$ is actually involved in 6 quadruplets and directly exchanges energy with 12 other components $(f_r',\theta')$. Because the values of $f'_r$ and $\theta'$ do not fall exacly on other discrete components, the spectral density is interpolated using a bilinear interpolation, so that each source term value +$S_{nl}(\bk)$ contains the direct exchange of energy with 48 other discrete components. +we compute the three contributions that correspond to the situation in which $\bk$ takes the role of $\bk$,$\bk_{2,+}$, $\bk_{2,-}$, $\bk_{3,+}$ and $\bk_{3,-}$ in the quadruplet, namely the full source term is, without making explicit that bilinear interpolation, \begin{eqnarray} -S_{\mathrm{nl}}(\bk) &=& -2 \left[\delta S_{\mathrm{nl}}(\bk,\bk_2,\bk_3,+)+\delta S_{\mathrm{nl}}(\bk,\bk_2,\bk_3,-)\right] \nonumber \\ - & & + \delta S_{\mathrm{nl}}(\bk_4,\bk,\bk_5,+) + \delta S_{\mathrm{nl}}(\bk_6,\bk,\bk_7,-) \\ - & & + \delta S_{\mathrm{nl}}(\bk_8,\bk_9,\bk, +) + \delta S_{\mathrm{nl}}(\bk_{10},\bk_{11},\bk, -) . \label{eq:diasum} +S_{\mathrm{nl}}(\bk) &=& -2 \left[\delta S_{\mathrm{nl}}(\bk,\bk_{2,+},\bk_{3,+})+\delta S_{\mathrm{nl}}(\bk,\bk_{2,-},\bk_{3,-})\right] \nonumber \\ + & & + \delta S_{\mathrm{nl}}(\bk_4,\bk,\bk_5) + \delta S_{\mathrm{nl}}(\bk_6,\bk,\bk_7) \\ + & & + \delta S_{\mathrm{nl}}(\bk_8,\bk_9,\bk) + \delta S_{\mathrm{nl}}(\bk_{10},\bk_{11},\bk) . \label{eq:diasum} \end{eqnarray} -with elementary contributions given by +where the geometry of the quadruplet $(\bk_4,\bk_4,\bk,\bk_5)$ is obtained from that of $(\bk,\bk,\bk_{2,+},\bk_{3,+})$ by a dilation by a factor $(1+\lambda)^2$ and rotation by the angle $\delta_{\theta,2}$; $(\bk_6,\bk_6,\bk,\bk_7)$ has the same dilation but the opposite rotation; $(\bk_8,\bk_8,\bk_9,\bk)$ is dilated by a factor $(1-\lambda)^2$ and rotated by the angle $-\delta_{\theta,3}$: and $(\bk_{10},\bk_{10},\bk_{11},\bk)$ is dilated by the same factor and rotated by the opposite angle. + + +The elementary contributions $\delta S_{\mathrm{nl}}(\bk_l,\bk_m,\bk_n)$ are given by %----------------------------% % Nonlinear interactions DIA % %----------------------------% % eq:snl_dia \begin{equation} -\delta S_{\mathrm{nl}}(\bk,\bk_2,\bk_3,s) = \frac{C}{g^4} f_{r,1}^{11} \left [ F^2 \left ( \frac{F_{2,s}}{(1+\lambda_{nl})^4} + - \frac{F_{3,s}}{(1-\lambda_{nl})^4} \right ) - \frac{2 F F_{2,s} F_{3,s}}{(1-\lambda_{nl}^2)^4} \right] , +\delta S_{\mathrm{nl}}(\bk_l,\bk_m,\bk_n) = \frac{C}{g^4} f_{r,l}^{11} \left [ F_l^2 \left ( \frac{F_m}{(1+\lambda)^4} + + \frac{F_n}{(1-\lambda)^4} \right ) - \frac{2 F_l F_m F_n}{(1-\lambda^2)^4} \right] , \label{eq:snl_dia} \end{equation} -where $s=+$ or $s=-$ is a sign index, and the spectral densities are $F = F(f_{r} ,\theta)$, $F_{2,+} = F(f_{r,2} ,\theta + \delta_{\theta,2})$, $F_{2,-} = F(f_{r,2} ,\theta - \delta_{\theta,2})$, etc. +where the spectral densities are $F_l = F(f_{r,l} ,\theta_l)$, etc. $C$ is a proportionality constant that was tuned to reproduce the inverse energy cascade. Default values for different source term packages are presented in Table~\ref{tab:snl_par}. -As a result, when accounting for the two quadruplet configurations, the source term at $\bk$ includes the interactions with -10 other spectral components. Besides, because $f_{r,2}$ and $f_{r,3}$ nor $\theta_{2,\pm} $ and $\theta_{3,\pm} $ fall on discretized frequencies and directions, the spectral densities are bilinearly interpolated, which involves 4 discrete spectral components for each of these 10 components. - % tab:snl_par \begin{table} \begin{center} \begin{tabular}{|l|c|c|} \hline - & $\lambda_{nl}$ & $C$ \\ \hline + & $\lambda$ & $C$ \\ \hline ST6 & 0.25 & $3.00 \; 10^7$ \\ \hline \wam-3 & 0.25 & $2.78 \; 10^7$ \\ \hline ST4 (Ardhuin et al.)& 0.25 & $2.50 \; 10^7$ \\ \hline diff --git a/model/src/w3gridmd.F90 b/model/src/w3gridmd.F90 index f4df6d143..5af38fc12 100644 --- a/model/src/w3gridmd.F90 +++ b/model/src/w3gridmd.F90 @@ -840,7 +840,7 @@ MODULE W3GRIDMD ! #ifdef W3_ST4 INTEGER :: SWELLFPAR, SDSISO, SDSBRFDF, SINTABLE,& - TAUWBUG, VISCSTRESS + TAUWBUG REAL :: SDSBCHOICE REAL :: ZWND, ALPHA0, Z0MAX, BETAMAX, SINTHP,& ZALP, Z0RAT, TAUWSHELTER, SWELLF, & @@ -857,7 +857,7 @@ MODULE W3GRIDMD SDSBM0, SDSBM1, SDSBM2, SDSBM3, & SDSBM4, SDSFACMTF, SDSCUMP, SDSNUW, & SDSL, SDSMWD, SDSMWPOW, SPMSS, SDSNMTF, SINTAIL1, SINTAIL2, & - CUMSIGP + CUMSIGP, VISCSTRESS #endif ! #ifdef W3_ST6 diff --git a/model/src/w3iogrmd.F90 b/model/src/w3iogrmd.F90 index fadd237f6..f8723d812 100644 --- a/model/src/w3iogrmd.F90 +++ b/model/src/w3iogrmd.F90 @@ -1536,12 +1536,14 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & #endif IF (SINTAILPAR(1).GT.0.5) THEN WRITE (NDSM) DELUST, DELTAIL, DELTAUW, DELU, DELALP, & - TAUT, TAUHFT, TAUHFT2 + TAUT, TAUHFT + IF (TTAUWSHELTER.GT.0) WRITE (NDSM) TAUHFT2 #ifdef W3_ASCII - WRITE (NDSA,*) 'DELUST, DELTAIL, DELTAUW, DELU, DELALP,& - TAUT, TAUHFT, TAUHFT2:', & + WRITE (NDSA,*) 'DELUST, DELTAIL, DELTAUW, DELU, DELALP,& + TAUT, TAUHFT:', & DELUST, DELTAIL, DELTAUW, DELU, DELALP, & - TAUT, TAUHFT, TAUHFT2 + TAUT, TAUHFT + IF (TTAUWSHELTER.GT.0) WRITE (NDSA,*) 'TAUHFT2:', TAUHFT2 #endif END IF ELSE @@ -1561,7 +1563,8 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT & CALL INSIN4(.FALSE.) READ (NDSM,END=801,ERR=802,IOSTAT=IERR) & DELUST, DELTAIL, DELTAUW, DELU, DELALP, & - TAUT, TAUHFT, TAUHFT2 + TAUT, TAUHFT + IF (TTAUWSHELTER.GT.0) READ(NDSM,END=801,ERR=802,IOSTAT=IERR) TAUHFT2 END IF END IF #endif diff --git a/model/src/w3snl1md.F90 b/model/src/w3snl1md.F90 index 598b627ea..09c096d2b 100644 --- a/model/src/w3snl1md.F90 +++ b/model/src/w3snl1md.F90 @@ -825,6 +825,8 @@ SUBROUTINE W3SNLGQM(A,CG,WN,DEPTH,TSTOTn,TSDERn) USE CONSTANTS, ONLY: TPI USE W3GDATMD, ONLY: SIG, NK , NTH , DTH, XFR, FR1, GQTHRSAT, GQAMP + IMPLICIT NONE + REAL, intent(in) :: A(NTH,NK), CG(NK), WN(NK) REAL, intent(in) :: DEPTH REAL, intent(out) :: TSTOTn(NTH,NK), TSDERn(NTH,NK) @@ -883,8 +885,8 @@ SUBROUTINE W3SNLGQM(A,CG,WN,DEPTH,TSTOTn,TSDERn) ! Gamma_max=1.3 (JFMAX>NF) TO OBTAIN IMPROVED RESULTS ! Note by Fabrice Ardhuin: this appears to give the difference in tail benaviour with Gerbrant's WRT !======================================================================= - JFMIN= 1-INT(LOG(1.0D0)/LOG(RAISF)) - JFMAX=NF+INT(LOG(1.3D0)/LOG(RAISF)) + JFMIN=MAX(1-INT(LOG(1.0D0)/LOG(RAISF)),1) + JFMAX=MIN(NF+INT(LOG(1.3D0)/LOG(RAISF)),NK) ! !======================================================================= ! COMPUTES THE SPECTRUM THRESHOLD VALUES (BELOW WHICH QNL4 IS NOT @@ -1065,7 +1067,7 @@ SUBROUTINE W3SNLGQM(A,CG,WN,DEPTH,TSTOTn,TSDERn) TEMP=(TB_TPM(IQ_OM2,JT1,JF1)*(( F(JT1P2P,JFM2)*CF2 *F(JT1P3M,JFM3)*CF3)* & (F(JT,JFM0 )*CF0*TB_V14(JF1)+F(JT1P ,JFM1)*CF1) & -SP0*SP1P*(SP1P2P*V3_4+SP1P3M*V2_4))+T_2M3P*(AUX05*AUX01-AUX02*AUX06)) *CP0 - WRITE(995,'(5I3,3E12.3)') ICONF,JF,JT, F(JT,JFM0) + WRITE(995,'(3I3,3E12.3)') ICONF,JF,JT, F(JT,JFM0) TEMP=(Q_2P3M+Q_2M3P) *CP1 WRITE(995,'(5I3,3E12.3)') ICONF,JF,JT,JT1P, JFM1,AUX00 *CP1, F(JT1P,JFM1),TSTOT(JT1P,JFM1) WRITE(995,'(5I3,3E12.3)') ICONF,JF,JT,JT1P2P,JFM2,-Q_2P3M*CP2,F(JT1P2P,JFM2),TSTOT(JT1P2P,JFM2) @@ -1219,6 +1221,8 @@ FUNCTION COUPLE(XK1 ,YK1 ,XK2 ,YK2 ,XK3 ,YK3 ,XK4 ,YK4) !/ ------------------------------------------------------------------- / USE CONSTANTS, ONLY: GRAV ! + IMPLICIT NONE + DOUBLE PRECISION, INTENT(IN) :: XK1 , YK1 , XK2 , YK2 DOUBLE PRECISION, INTENT(IN) :: XK3 , YK3 DOUBLE PRECISION, INTENT(IN) :: XK4 , YK4 @@ -1305,6 +1309,7 @@ SUBROUTINE GAULEG (W_LEG ,X_LEG ,NPOIN) !/ ------------------------------------------------------------------- / !.....VARIABLES IN ARGUMENT ! """""""""""""""""""" + IMPLICIT NONE INTEGER , INTENT(IN) :: NPOIN DOUBLE PRECISION ,INTENT(INOUT) :: W_LEG(NPOIN) , X_LEG(NPOIN) ! @@ -1552,6 +1557,7 @@ SUBROUTINE INSNLGQM #ifdef W3_S CALL STRACE (IENT, 'INSNLGQM') #endif + IMPLICIT NONE !.....LOCAL VARIABLES INTEGER JF , JT , JF1 , JT1 , NF1P1 , IAUX , NT , NF , IK INTEGER IQ_TE1 , IQ_OM2 , LBUF , DIMBUF , IQ_OM1 , NQ_TE1 , NCONFM @@ -2084,10 +2090,7 @@ SUBROUTINE INSNLGQM AUX=0.0D0 DO JT1=1,GQNT1 DO IQ_OM2=1,GQNQ_OM2 - AAA=TB_FAC(IQ_OM2,JT1,JF1)*TB_TPM(IQ_OM2,JT1,JF1) - IF (AAA.GT.AUX) AUX=AAA - CCC=TB_FAC(IQ_OM2,JT1,JF1)*TB_TMP(IQ_OM2,JT1,JF1) - IF (CCC.GT.AUX) AUX=CCC + AUX=MAX(AUX,TB_FAC(IQ_OM2,JT1,JF1)*TB_TPM(IQ_OM2,JT1,JF1),TB_FAC(IQ_OM2,JT1,JF1)*TB_TMP(IQ_OM2,JT1,JF1)) ENDDO ENDDO MAXCLA(JF1)=AUX @@ -2099,6 +2102,7 @@ SUBROUTINE INSNLGQM DO JF1=1,GQNF1 IF (MAXCLA(JF1).GT.AUX) AUX=MAXCLA(JF1) ENDDO + TEST1=SEUIL1*AUX ! !.....Set to zero the coupling coefficients not used @@ -2128,7 +2132,9 @@ SUBROUTINE INSNLGQM ! !..... counts the fraction of the eliminated configurations ELIM=(1.D0-DBLE(NCONF)/DBLE(NCONFM))*100.D0 - ! WRITE(994,*) 'NCONF:',NCONF,ELIM +#ifdef W3_TGQM + WRITE(994,*) 'NCONF, ELIM FRACTION:',NCONF,ELIM +#endif END SUBROUTINE INSNLGQM !/ !/ End of module W3SNL1MD -------------------------------------------- / From 68377723e29f2afde808abffc6d9e3073ae86a49 Mon Sep 17 00:00:00 2001 From: Fabrice Ardhuin Date: Fri, 1 Dec 2023 15:54:59 +0100 Subject: [PATCH 13/15] Corrected type of VISCSTRESS from INTEGER to REAL --- model/src/w3gridmd.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/model/src/w3gridmd.F90 b/model/src/w3gridmd.F90 index f4df6d143..efac485a8 100644 --- a/model/src/w3gridmd.F90 +++ b/model/src/w3gridmd.F90 @@ -840,7 +840,7 @@ MODULE W3GRIDMD ! #ifdef W3_ST4 INTEGER :: SWELLFPAR, SDSISO, SDSBRFDF, SINTABLE,& - TAUWBUG, VISCSTRESS + TAUWBUG REAL :: SDSBCHOICE REAL :: ZWND, ALPHA0, Z0MAX, BETAMAX, SINTHP,& ZALP, Z0RAT, TAUWSHELTER, SWELLF, & @@ -853,7 +853,7 @@ MODULE W3GRIDMD SDSBR, SDSP, SDSBT, SDS4A, SDKOF, & SDSCOS, SDSDTH, SDSBCK, SDSABK, & SDSPBK, SDSBINT, SDSHCK, & - SDSBRF1, & + SDSBRF1, VISCSTRESS, & SDSBM0, SDSBM1, SDSBM2, SDSBM3, & SDSBM4, SDSFACMTF, SDSCUMP, SDSNUW, & SDSL, SDSMWD, SDSMWPOW, SPMSS, SDSNMTF, SINTAIL1, SINTAIL2, & From b00cd4eba94b35c58b9ae3c911d5559473bb8e15 Mon Sep 17 00:00:00 2001 From: ma05cde Date: Fri, 1 Dec 2023 15:27:15 +0000 Subject: [PATCH 14/15] correct typo --- model/src/w3gridmd.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/model/src/w3gridmd.F90 b/model/src/w3gridmd.F90 index 224144151..5af38fc12 100644 --- a/model/src/w3gridmd.F90 +++ b/model/src/w3gridmd.F90 @@ -853,7 +853,7 @@ MODULE W3GRIDMD SDSBR, SDSP, SDSBT, SDS4A, SDKOF, & SDSCOS, SDSDTH, SDSBCK, SDSABK, & SDSPBK, SDSBINT, SDSHCK, & - SDSBRF1, VISCSTRESS, & + SDSBRF1, & SDSBM0, SDSBM1, SDSBM2, SDSBM3, & SDSBM4, SDSFACMTF, SDSCUMP, SDSNUW, & SDSL, SDSMWD, SDSMWPOW, SPMSS, SDSNMTF, SINTAIL1, SINTAIL2, & From dbe8e6c548a29d39c19c034c7d6f27a376a14775 Mon Sep 17 00:00:00 2001 From: Fabrice Ardhuin Date: Thu, 4 Jan 2024 10:29:37 +0100 Subject: [PATCH 15/15] Fixed bug with undefined LEVTAIL used instead of LEVTAIL0 when table is used in W3SIN4 --- model/src/w3src4md.F90 | 13 ++++++++----- 1 file changed, 8 insertions(+), 5 deletions(-) diff --git a/model/src/w3src4md.F90 b/model/src/w3src4md.F90 index 07142f7ac..32eadaf82 100644 --- a/model/src/w3src4md.F90 +++ b/model/src/w3src4md.F90 @@ -850,9 +850,12 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & TEMP=TEMP+A(IS)*(MAX(COSWIND,0.))**3 END DO ! - LEVTAIL0= CONST0*TEMP ! LEVTAIL is sum of A(k,theta)*cos^3(theta-wind)*DTH*SIG^5/(g^2*2pi)*2*pi*SIG/CG - ! ! which is the same as sum of E(f,theta)*cos^3(theta-wind)*DTH*SIG^5/(g^2*2pi) - ! reminder: sum of E(f,theta)*DTH*SIG^5/(g^2*2pi) is 2*k^3*E(k) + LEVTAIL0= CONST0*TEMP ! LEVTAIL is sum over theta of A(k,theta)*cos^3(theta-wind)*DTH*SIG^5/(g^2*2pi)*2*pi*SIG/CG + ! which is the same as sum of E(f,theta)*cos^3(theta-wind)*DTH*SIG^5/(g^2*2pi) + ! reminder: sum of E(f,theta)*DTH*SIG^5/(g^2*2pi) is 2*k^3*E(k) +! +! Computation of stress supported by tail: uses table if SINTAILPAR(1)=1 , correspoding to SINTABLE = 1 +! IF (SINTAILPAR(1).LT.0.5) THEN ALLOCATE(W(JTOT)) W(2:JTOT-1)=1. @@ -896,7 +899,7 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & ZLOG = MIN(ALOG(ZMU),0.) ZBETA = CONST1*ZMU*ZLOG**4 ! - ! Optional addition of capillary wave peak + ! Optional addition of capillary wave peak if SINTAIL2=1 ! IF (SINTAILPAR(3).GT.0) THEN IF (USTR.LT.CM) THEN @@ -939,7 +942,7 @@ SUBROUTINE W3SIN4 (A, CG, K, U, USTAR, DRAT, AS, USDIR, Z0, CD, & DELJ1= MAX(0.,MIN (1. , XJ-FLOAT(J))) DELJ2=1. - DELJ1 IF (TTAUWSHELTER.GT.0) THEN - XK = LEVTAIL / DELTAIL + XK = LEVTAIL0/ DELTAIL I = MIN (ILEVTAIL-1, INT(XK)) DELK1= MIN (1. ,XK-FLOAT(I)) DELK2=1. - DELK1