Ww3 diffraction

model/inp/ww3_grid.inp

@@ -317,45 +317,30 @@ $   limitation and the GSE alleviation.
 $
 $ Unstructured grids ------------------------------------------------ $
 $   UNST parameters       : Namelist UNST
-$                           UGBCCFL    : Turns on/off (TRUE/FALSE) the computation of the CFL number on the physical domain boundary. 
-$                                        If FALSE the explicit scheme can be much faster though stability is not guaranteed
-$                                        (default TRUE)
-$                           UGOBCAUTO  : TRUE: OBC points are taken from  type 15 elements (default)
-$                                        FALSE: OBC points must be listed in ww3_grid.inp
-$                           UGOBCDEPTH : Threshold ( < 0) depth (default -10) for OBC points if UGOBCAUTO is TRUE
-$                           UGOBCFILE  : File name of file for reading boudary (default 'unset') 
-$
-$                           The following are TRUE/FALSE variables and only one can be TRUE. 
-$                           By default, EXPFSN is TRUE and must be set to false to activate another option  
-$                           EXPFSN    : Activation of N scheme (default option)
+$                           UGOBCAUTO : TRUE: OBC points are taken from  type 15 elements
+$                                       FALSE: OBC points must be listed in ww3_grid.inp
+$                           UGOBCDEPTH: Threshold ( < 0) depth for OBC points if UGOBCAUTO is TRUE
+$                           EXPFSN    : Activation of N scheme
 $                           EXPFSPSI  : Activation of PSI scheme
 $                           EXPFSFCT  : Activation of FCT scheme
 $                           IMPFSN    : Activation of N implicit scheme
-$                           EXPTOTAL  : Activation of Block explicit N scheme solver
-$                           IMPTOTAL  : Activation of fully implicit scheme, non splitting
-$                             
-$                           The following TRUE/FALSE variables are only for IMPTOTAL=TRUE
-$                           IMPREFRACTION : Turn on implicit freq. shift (default FALSE) 
-$                           IMPFREQSHIFT  : Turn on implicit freq. shift terms (default FALSE)
-$                           IMPSOURCE     : Turn on implicit source terms (default FALSE)
-$
-$                           JGS_TERMINATE_MAXITER  : Terminate based on max number of iterations (TRUE/FALSE, default TRUE)
-$                           JGS_TERMINATE_DIFFERENCE  : Terminate based on the total change of the unweightet sum of wave action (TRUE/FALSE, default TRUE)
-$                           JGS_TERMINATE_NORM  : Terminate based on the norm of the solution (TRUE/FALSE, default FALSE)
-$                           JGS_USE_JACOBI  : Use Jacobi solver family (TRUE/FALSE, default TRUE)
-$                           JGS_BLOCK_GAUSS_SEIDEL  : Use Block Gauss Seidel method for imptotal instead of the conservative jacobi iterator. (TRUE/FALSE, default TRUE)
-$                           JGS_MAXITER  : Max. Number of solver iterations for JGS_TERMINATE_MAXITER (integer, default 100) 
-$                           JGS_PMIN  : % of grid points that do not need to converge during solver iteration (real, default 1) 
-$                           JGS_DIFF_THR  : Implicit solver threshold for JGS_TERMINATE_DIFFERENCE (real, default 1.0e-10)
-$                           JGS_NORM_THR  : Norm of the solution for JGS_TERMINATE_NORM (real, default 1.0e-20)
-$                           JGS_LIMITER  : TRUE: Use total (quasi-steady: limits whole equation) instead of local limiter (un-steady: limits only source terms)
-$                                          FALSE: default 
-$                           JGS_LIMITER_FUNC  : 1 - old limiter (default)
-$                                               2 - alternatnive limiter
-$                           SETUP_APPLY_WLV  : Compute wave setup (TRUE/FALSE, default TRUE) 
-$                           SOLVERTHR_SETUP  : Solver threshold for setup computations (default 1E-6)
-$                           CRIT_DEP_SETUP  : Critical depth for setup computations (default 0.1) 
-
+$                           IMPTOTAL  : Activation of fully implicit scheme | Non splitting
+$                           EXPTOTAL  : Turn on implicit refraction (only with imptotal)
+$                           IMPREFRACTION  : Turn on implicit freq. shift (only with imptotal)
+$                           IMPFREQSHIFT  : Turn on implicit freq. shift terms (only with imptotal)
+$                           IMPSOURCE  : Turn on implicit source terms (only with imptotal)
+$                           JGS_TERMINATE_MAXITER  : max. Number of iterations
+$                           JGS_TERMINATE_DIFFERENCE  : terminate based on the total change of wave action
+$                           JGS_TERMINATE_NORM  : terminate based on the norm of the solution
+$                           JGS_USE_JACOBI  : Use Jacobi solver for imptotal
+$                           JGS_BLOCK_GAUSS_SEIDEL  : Use Block Gauss Seidel method for imptotal
+$                           JGS_MAXITER  : max. Number of solver iterations
+$                           JGS_PMIN  : % of grid points that do not need to converge during solver iteration.
+$                           JGS_DIFF_THR  : implicit solver threshold for JGS_TERMINATE_DIFFERENCE
+$                           JGS_NORM_THR  : terminate based on the norm of the solution
+$                           SETUP_APPLY_WLV  : Compute wave setup (experimental)
+$                           SOLVERTHR_SETUP  : Solver threshold for setup computations
+$                           CRIT_DEP_SETUP  : Critical depths for setup computations
 $
 $ SMC grid propagation    : Namelist PSMC and default values
 $                           CFLSM  : Maximum CFL no. for propagation, 0.7

model/src/CMakeLists.txt

@@ -47,16 +47,15 @@ target_include_directories(ww3_lib
 
 # Set compiler flags.
 if(CMAKE_Fortran_COMPILER_ID MATCHES "^(Intel)$")
-  set(compile_flags -no-fma -ip -g -traceback -i4 -real-size 32 -fp-model precise
-    -assume byterecl -fno-alias -fno-fnalias)
+  set(compile_flags)
   # -sox only works on Linux
   if(LINUX)
     list(APPEND compile_flags -sox)
   endif()
-  set(compile_flags_release -O3)
+  set(compile_flags_release)
   # SHELL: prefix fixes CMake attempting to de-duplicate the repeated uses of 'all' in -warn, -debug, -check
   # See https://cmake.org/cmake/help/latest/command/target_compile_options.html#option-de-duplication
-  set(compile_flags_debug -O0 "SHELL:-debug all" "SHELL:-warn all" "SHELL:-check all" -check noarg_temp_created -fp-stack-check -heap-arrays -traceback -fpe0)
+  set(compile_flags_debug)
 
   if(APPLE)
     # The linker on macOS does not include `common symbols` (usually module variables without a default value) by default

model/src/w3dispmd.F90

@@ -343,7 +343,7 @@ SUBROUTINE WAVNU2 (W,H,K,CG,EPS,NMAX,ICON)
     !/
   END SUBROUTINE WAVNU2
   !/
-  PURE SUBROUTINE WAVNU3 (SI,H,K,CG)
+  PURE SUBROUTINE WAVNU3 (SI,H,K,CG,C)
     !/
     !/                  +-----------------------------------+
     !/                  | WAVEWATCH III           NOAA/NCEP |
@@ -420,7 +420,7 @@ PURE SUBROUTINE WAVNU3 (SI,H,K,CG)
     !/ Parameter list
     !/
     REAL, INTENT(IN)        :: SI, H
-    REAL, INTENT(OUT)       :: K, CG
+    REAL, INTENT(OUT)       :: K, CG, C
     !/
     !/ ------------------------------------------------------------------- /
     !/ Local parameters
@@ -445,6 +445,7 @@ PURE SUBROUTINE WAVNU3 (SI,H,K,CG)
     KH  = KH0 * (1 + KH0**1.09 * 1./EXP(MIN(KDMAX,TMP))) / SQRT(TANH(MIN(KDMAX,KH0)))
     K   = KH/H
     CG  = 0.5*(1+(2*KH/SINH(MIN(KDMAX,2*KH))))*SI/K
+    C   = CG/K
     !
     RETURN
     !/
@@ -533,15 +534,15 @@ PURE SUBROUTINE WAVNU_LOCAL (SIG,DW,WNL,CGL)
     !/
     !/ ------------------------------------------------------------------- /
     !/ Local parameters
-    REAL                    :: DEPTH
+    REAL                    :: DEPTH, CP
     !
     !/
     !/ End of WAVNU_LOCAL------------------------------------------------- /
     !/
 
     DEPTH  = MAX ( DMIN , DW)
     !
-    CALL WAVNU3(SIG,DEPTH,WNL,CGL)
+    CALL WAVNU3(SIG,DEPTH,WNL,CGL,CP)
 
   END SUBROUTINE WAVNU_LOCAL
   !/

model/src/w3gdatmd.F90

@@ -758,6 +758,9 @@ MODULE W3GDATMD
     INTEGER*2, POINTER   :: IOBP_LOC(:)
     INTEGER*1, POINTER   :: IOBDP_LOC(:)
     INTEGER*1, POINTER   :: IOBPA_LOC(:)
+    REAL, POINTER      :: DIFRM(:)
+    REAL, POINTER      :: DIFRX(:)
+    REAL, POINTER      :: DIFRY(:)
 #endif
 
     REAL(8), POINTER      :: LEN(:,:),SI(:), IEN(:,:)
@@ -1043,6 +1046,7 @@ MODULE W3GDATMD
     LOGICAL :: FSFREQSHIFT  = .FALSE.
     LOGICAL :: FSSOURCE     = .FALSE.
     LOGICAL :: FSBCCFL      = .FALSE.
+    LOGICAL :: L2DM         = .TRUE. 
     LOGICAL :: DO_CHANGE_WLV
     REAL(8) :: SOLVERTHR_STP
     REAL(8) :: CRIT_DEP_STP
@@ -1058,6 +1062,8 @@ MODULE W3GDATMD
     REAL*8  :: B_JGS_DIFF_THR
     REAL*8  :: B_JGS_NORM_THR
     INTEGER :: B_JGS_NLEVEL
+    LOGICAL :: B_JGS_LDIFR
+    INTEGER :: B_JGS_IDIFR 
     LOGICAL :: B_JGS_SOURCE_NONLINEAR
   END TYPE SCHM
   !
@@ -1116,6 +1122,9 @@ MODULE W3GDATMD
   INTEGER*2, POINTER   :: IOBP_LOC(:)
   INTEGER*1, POINTER   :: IOBDP_LOC(:)
   INTEGER*1, POINTER   :: IOBPA_LOC(:)
+  REAL, POINTER        :: DIFRM(:)
+  REAL, POINTER        :: DIFRX(:)
+  REAL, POINTER        :: DIFRY(:)
 #endif
 
   REAL(8), POINTER        :: IEN(:,:), LEN(:,:), SI(:)
@@ -1403,6 +1412,7 @@ MODULE W3GDATMD
   LOGICAL, POINTER  :: FSN,FSPSI,FSFCT,FSNIMP,FSTOTALIMP,FSTOTALEXP
   LOGICAL, POINTER  :: FSREFRACTION, FSFREQSHIFT, FSSOURCE, FSBCCFL
   LOGICAL, POINTER  :: DO_CHANGE_WLV
+  LOGICAL, POINTER  :: L2DM 
   REAL(8), POINTER  :: SOLVERTHR_STP
   REAL(8), POINTER  :: CRIT_DEP_STP
   LOGICAL, POINTER :: B_JGS_TERMINATE_MAXITER
@@ -1417,6 +1427,8 @@ MODULE W3GDATMD
   REAL(8), POINTER :: B_JGS_DIFF_THR
   REAL(8), POINTER :: B_JGS_NORM_THR
   INTEGER, POINTER :: B_JGS_NLEVEL
+  LOGICAL, POINTER :: B_JGS_LDIFR
+  INTEGER, POINTER :: B_JGS_IDIFR
   LOGICAL, POINTER :: B_JGS_SOURCE_NONLINEAR
   !/
   !/ Data aliasing for structure SICP(S)
@@ -2840,6 +2852,7 @@ SUBROUTINE W3SETG ( IMOD, NDSE, NDST )
     FSREFRACTION => MPARS(IMOD)%SCHMS%FSREFRACTION
     FSFREQSHIFT => MPARS(IMOD)%SCHMS%FSFREQSHIFT
     FSSOURCE => MPARS(IMOD)%SCHMS%FSSOURCE
+    L2DM => MPARS(IMOD)%SCHMS%L2DM
     DO_CHANGE_WLV => MPARS(IMOD)%SCHMS%DO_CHANGE_WLV
     SOLVERTHR_STP => MPARS(IMOD)%SCHMS%SOLVERTHR_STP
     CRIT_DEP_STP => MPARS(IMOD)%SCHMS%CRIT_DEP_STP
@@ -2851,6 +2864,8 @@ SUBROUTINE W3SETG ( IMOD, NDSE, NDST )
     B_JGS_BLOCK_GAUSS_SEIDEL => MPARS(IMOD)%SCHMS%B_JGS_BLOCK_GAUSS_SEIDEL
     B_JGS_MAXITER => MPARS(IMOD)%SCHMS%B_JGS_MAXITER
     B_JGS_LIMITER_FUNC => MPARS(IMOD)%SCHMS%B_JGS_LIMITER_FUNC
+    B_JGS_LDIFR => MPARS(IMOD)%SCHMS%B_JGS_LDIFR
+    B_JGS_IDIFR => MPARS(IMOD)%SCHMS%B_JGS_IDIFR
     B_JGS_PMIN => MPARS(IMOD)%SCHMS%B_JGS_PMIN
     B_JGS_DIFF_THR => MPARS(IMOD)%SCHMS%B_JGS_DIFF_THR
     B_JGS_NORM_THR => MPARS(IMOD)%SCHMS%B_JGS_NORM_THR

model/src/w3gridmd.F90

@@ -926,12 +926,15 @@ MODULE W3GRIDMD
   LOGICAL :: IMPREFRACTION
   LOGICAL :: IMPFREQSHIFT
   LOGICAL :: IMPSOURCE
+  LOGICAL :: LREAD2DM
   LOGICAL :: SETUP_APPLY_WLV
   INTEGER :: JGS_MAXITER
   INTEGER :: nbSel
   INTEGER :: UNSTSCHEMES(6)
   INTEGER :: UNSTSCHEME
   INTEGER :: JGS_NLEVEL
+  LOGICAL :: JGS_LDIFR
+  INTEGER :: JGS_IDIFR
   REAL*8  :: JGS_PMIN
   REAL*8  :: JGS_DIFF_THR
   REAL*8  :: JGS_NORM_THR
@@ -1087,7 +1090,7 @@ MODULE W3GRIDMD
        UGBCCFL, EXPFSN, EXPFSPSI, EXPFSFCT,       &
        IMPFSN, IMPTOTAL, EXPTOTAL,                &
        IMPREFRACTION, IMPFREQSHIFT,               &
-       IMPSOURCE,                                 &
+       IMPSOURCE, LREAD2DM,                       &
        JGS_TERMINATE_MAXITER,                     &
        JGS_TERMINATE_DIFFERENCE,                  &
        JGS_TERMINATE_NORM,                        &
@@ -1102,7 +1105,8 @@ MODULE W3GRIDMD
        JGS_NLEVEL,                                &
        JGS_SOURCE_NONLINEAR,                      &
        SETUP_APPLY_WLV, SOLVERTHR_SETUP,          &
-       CRIT_DEP_SETUP
+       CRIT_DEP_SETUP, JGS_LDIFR, JGS_IDIFR
+
   NAMELIST /MISC/ CICE0, CICEN, LICE, XSEED, FLAGTR, XP, XR, &
        XFILT, PMOVE, IHM, HSPM, WSM, WSC, FLC, FMICHE, &
        RWNDC, FACBERG, NOSW, GSHIFT, WCOR1, WCOR2,     &
@@ -2441,6 +2445,7 @@ SUBROUTINE W3GRID()
     IMPREFRACTION = .FALSE.
     IMPFREQSHIFT = .FALSE.
     IMPSOURCE = .FALSE.
+    LREAD2DM   = .FALSE. 
     SETUP_APPLY_WLV = .TRUE.
     SOLVERTHR_SETUP=1E-6
     CRIT_DEP_SETUP=0.1
@@ -2457,9 +2462,12 @@ SUBROUTINE W3GRID()
     JGS_NORM_THR = 1.E-20
     JGS_NLEVEL = 0
     JGS_SOURCE_NONLINEAR = .FALSE.
+    JGS_LDIFR = .FALSE.
+    JGS_IDIFR = 1
     ! read data from the unstructured devoted namelist
     CALL READNL ( NDSS, 'UNST', STATUS )
 
+    L2DM             = LREAD2DM
     B_JGS_USE_JACOBI = JGS_USE_JACOBI
     B_JGS_TERMINATE_MAXITER = JGS_TERMINATE_MAXITER
     B_JGS_TERMINATE_DIFFERENCE = JGS_TERMINATE_DIFFERENCE
@@ -2472,6 +2480,8 @@ SUBROUTINE W3GRID()
     B_JGS_DIFF_THR = JGS_DIFF_THR
     B_JGS_NORM_THR = JGS_NORM_THR
     B_JGS_NLEVEL = JGS_NLEVEL
+    B_JGS_LDIFR  = JGS_LDIFR
+    B_JGS_IDIFR  = JGS_IDIFR
     B_JGS_SOURCE_NONLINEAR = JGS_SOURCE_NONLINEAR
 
     nbSel=0
@@ -3318,10 +3328,10 @@ SUBROUTINE W3GRID()
       WRITE (NDSO,2953) CFLTM, WDTHCG, WDTHTH
 #endif
       !
-      WRITE (NDSO,2956) UGBCCFL, UGOBCAUTO, UGOBCDEPTH,TRIM(UGOBCFILE), &
+      WRITE (NDSO,*) UGBCCFL, UGOBCAUTO, UGOBCDEPTH,TRIM(UGOBCFILE), &
            EXPFSN, EXPFSPSI, EXPFSFCT, IMPFSN, EXPTOTAL,&
            IMPTOTAL, IMPREFRACTION, IMPFREQSHIFT,      &
-           IMPSOURCE, SETUP_APPLY_WLV,                 &
+           IMPSOURCE, LREAD2DM, SETUP_APPLY_WLV,        &
            JGS_TERMINATE_MAXITER,                      &
            JGS_TERMINATE_DIFFERENCE,                   &
            JGS_TERMINATE_NORM,                         &
@@ -4000,7 +4010,11 @@ SUBROUTINE W3GRID()
       !
       ! Reading depths on unstructured grid (this also sets number of mesh points, NX)
       !
-      CALL READMSH(NDSG,FNAME)
+      IF (L2DM) THEN
+        CALL READ2DM(NDSG,FNAME)
+      ELSE
+        CALL READMSH(NDSG,FNAME)
+      ENDIF
       ALLOCATE(ZBIN(NX, NY),OBSX(NX,NY),OBSY(NX,NY))
       ZBIN(:,1) = VSC * ZB(:)
       !
@@ -4440,6 +4454,9 @@ SUBROUTINE W3GRID()
          CALL READMSHOBC(NDSG,UGOBCFILE,TMPSTA,UGOBCOK)
     IF ((GTYPE.EQ.UNGTYPE).AND.UGOBCAUTO.AND.(.NOT.UGOBCOK))  &
          CALL UG_GETOPENBOUNDARY(TMPSTA,ZBIN,UGOBCDEPTH)
+    IF ((GTYPE.EQ.UNGTYPE).AND.L2DM.AND.(.NOT.UGOBCOK))  &
+         CALL READ2DM_TMPSTA(23956,trim("meshbnd.2dm"),TMPSTA)     
+
     !
     ! 8.b Determine where to get the data
     !

model/src/w3initmd.F90

@@ -421,7 +421,7 @@ SUBROUTINE W3INIT ( IMOD, IsMulti, FEXT, MDS, MTRACE, ODAT, FLGRD,  FLGR2, FLGD,
     use yowNodepool, only: npa
     use yowRankModule, only : rank
 #endif
-    USE W3GDATMD, ONLY: GTYPE, UNGTYPE
+    USE W3GDATMD, ONLY: GTYPE, UNGTYPE, IOBDP_LOC
 #ifdef W3_PDLIB
     USE PDLIB_W3PROFSMD, ONLY : PDLIB_MAPSTA_INIT, SET_IOBDP_PDLIB, PDLIB_IOBP_INIT, SET_IOBPA_PDLIB
     USE PDLIB_W3PROFSMD, ONLY : BLOCK_SOLVER_INIT, BLOCK_SOLVER_EXPLICIT_INIT, PDLIB_INIT, DEALLOCATE_PDLIB_GLOBAL
@@ -499,7 +499,7 @@ SUBROUTINE W3INIT ( IMOD, IsMulti, FEXT, MDS, MTRACE, ODAT, FLGRD,  FLGR2, FLGD,
 #endif
     REAL                    :: DTTST, DEPTH, FRACOS
     REAL                    :: FACTOR
-    REAL                    :: WLVeff
+    REAL                    :: WLVeff, CTMP
 #ifdef W3_T
     REAL, ALLOCATABLE      :: XOUT(:,:)
 #endif
@@ -1336,6 +1336,7 @@ SUBROUTINE W3INIT ( IMOD, IsMulti, FEXT, MDS, MTRACE, ODAT, FLGRD,  FLGR2, FLGD,
 #ifdef W3_PDLIB
     IF ( IAPROC .LE. NAPROC ) THEN
       CALL SET_IOBDP_PDLIB
+      WRITE(*,*) MAXVAL(IOBDP_LOC), MINVAL(IOBDP_LOC) 
     ENDIF
 #endif
 
@@ -1386,7 +1387,7 @@ SUBROUTINE W3INIT ( IMOD, IsMulti, FEXT, MDS, MTRACE, ODAT, FLGRD,  FLGR2, FLGD,
         !
         !         Calculate wavenumbers and group velocities.
 #ifdef W3_PDLIB
-        CALL WAVNU3(SIG(IK),DEPTH,WN(IK,IS),CG(IK,IS))
+        CALL WAVNU3(SIG(IK),DEPTH,WN(IK,IS),CG(IK,IS), CTMP)
 #else
         CALL WAVNU1(SIG(IK),DEPTH,WN(IK,IS),CG(IK,IS))
 #endif

model/src/w3iogrmd.F90

@@ -796,6 +796,8 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT &
              B_JGS_DIFF_THR,                                  &
              B_JGS_NORM_THR,                                  &
              B_JGS_NLEVEL,                                    &
+             B_JGS_LDIFR,                                     & 
+             B_JGS_IDIFR,                                     &
              B_JGS_SOURCE_NONLINEAR
 #ifdef W3_ASCII
         WRITE (NDSA,*)                                        &
@@ -995,6 +997,8 @@ SUBROUTINE W3IOGR ( INXOUT, NDSM, IMOD, FEXT &
              B_JGS_DIFF_THR,                                  &
              B_JGS_NORM_THR,                                  &
              B_JGS_NLEVEL,                                    &
+             B_JGS_LDIFR,                                     &
+             B_JGS_IDIFR,                                     &
              B_JGS_SOURCE_NONLINEAR
         IF (.NOT. GUGINIT) THEN
           CALL W3DIMUG ( IGRD, NTRI, NX, COUNTOT, NNZ, NDSE, NDST )