Environmental Modeling Center Environmental Modeling Center Environmental Modeling Center United States Department of Commerce
 

WAVEWATCH III v. 3.14 errata, problems and fixes

The content provided on this page supports model development. These are not official NWS products and should not be relied upon for operational purposes. This web site is not subject to 24/7 support, and thus may be unavailable during system outages.

Operationally generated graphics of the wave fields (no spectra or source terms) are available from Model Analyses and Guidance.

Bulletin files are available from the Production FTP/HTTPS server. Look for gfs.YYYYMMDD/CC/wave/station/bulls.tCCz/gfswave.stationID.bull


This page presents errata to the WAVEWATCH III v. 3.14 manuals, and a list of identified problems and fixes. It will be updated as needed. Please let us know if you find any errors or problems so that we can add them to this page.

Errata

  • The example ww3_prep.inp file on page 86-88 does not include various new options. Use the new ww3_prep.inp file from the errata list below instead.
  • ....

Index of Problems and Fixes


USTAR initialization in flux computations :


Problem: In the flux computations in the subroutines W3FLX2 and W3FLX3 the friction velocity UST is declared with INTENT(OUT), but not initialized. In the main wave model code this is not an issue for most compilers, as the friction velocity is defined before the subroutine call. Some strict compilers, however, will initialize this as 0, resulting in run time model errors. Furthermore, these variables were not properly initialized in the output post-processors ww3_outp and gx_outp.
Fix: In subroutines W3FLX2 in file w3flx2md.ftn and W3FLX3 in file w3flx3md.ftn, replace the declarations

      REAL, INTENT(OUT)       :: UST, USTD, Z0, CD
with
      REAL, INTENT(INOUT)     :: UST
      REAL, INTENT(OUT)       :: USTD, Z0, CD

Fix: In the post-processor gx_outp.ftn, on lines 825 and following, add
!/ST3            TAUWY  = 0. 
!/ST3            LLWS(:,:)  = .TRUE. 
            USTAR  = 1. 
Fix: In the post-processor ww3_outp.ftn, on lines 1195 and following, add
!/ST2                ZWND   = ZWIND 
!/ST3                ZWND   = ZZWND 
!/ST0                USTAR  = 1. 
!/ST1                USTAR  = 1. 
!/ST2                USTAR  = 1. 
Identified by Jesus Portilla, Posted April 20, 2010
Back to index.

Obstruction in mosaic grids :


Problem: When using a mosaic of two-way nested grids, artificial obstructions may appear at the boundaries of inner grids, when obstructions are defined between grid points (FLAGTR = 1 or FLAGTR = 3 in ww3_grid.inp).
Fix: We have not yet pursued a software fix for this. For now, use no obstructions (FLAGTR = 0) of obstructions defined at the grid points instead (FLAGTR = 2 or FLAGTR = 4, the latter is the NCEP operational default approach).

Identified by Chris Bunney, Posted May 31, 2010
Back to index.

Errors in ww3_prep.inp :


Problem: The version of ww3_prep;inp distributed with the code and included in the manual does not include several new options.
Fix: Use the replacement ww3_prep.inp file available here.
May 31, 2010
Back to index.

Saving data in W3IOBC :


Problem: The subroutine W3IOBC processes boundary data files. When reading data for which the discrete spectral grid is different than the spectral grid for the model, this routine converts the spectra to the spectral model grid. This works well on the first reading of the file, but corresponding data for the input spectra are not save, and hence subsequent reading from the file fails.
Quick fix: *** THIS FIX IS TO BE USED WITH ww3_shel OR WITH ww3_multi ONLY WHEN ALL BOUNDARY DATA FILES HAVE THE SAME SPECTRAL RESOLUTION *** In subroutine W3IOBC in file w3iobcmd.ftn replace lines 196 - 204

!/ ------------------------------------------------------------------- /
!/ Local parameters
!/
      INTEGER                 :: IFILE, IERR, I, J, NKI, NTHI, IX, IY,&
                                 ISEA, IP, ISP, NPTS, ISOUT, IS, IGRD
!/T1      INTEGER               :: IK, ITH
!/S      INTEGER, SAVE           :: IENT = 0
      REAL                    :: XFRI, FR1I, TH1I
!/T1      REAL                  :: HS, HS0
with
!/ ------------------------------------------------------------------- /
!/ Local parameters
!/
      INTEGER, SAVE           :: IFILE, IERR, I, J, NKI, NTHI, IX, IY,&
                                 ISEA, IP, ISP, NPTS, ISOUT, IS, IGRD
!/T1      INTEGER               :: IK, ITH
!/S      INTEGER, SAVE           :: IENT = 0
      REAL, SAVE              :: XFRI, FR1I, TH1I
!/T1      REAL                  :: HS, HS0
Full fix: Add relevant parameters to be saved to the data structure in w3odatmd.ftn to be saved for each separate grid individually. In file w3odatmd.ftn, make the following additions (green text) around lines 335,
 !/
      TYPE OTYPE5
        INTEGER               :: NBI, NBI2, NFBPO, NBO(0:9),          &
                                 NBO2(0:9), NDSL(9), NKI, NTHI
!/MPI        INTEGER               :: NRQBP, NRQBP2
        INTEGER, POINTER      :: IPBPI(:,:), ISBPI(:),                &
                                 IPBPO(:,:), ISBPO(:)
!/MPI        INTEGER, POINTER      :: IRQBP1(:), IRQBP2(:)
        REAL                  :: XFRI, FR1I, TH1I
        REAL, POINTER         :: XBPI(:), YBPI(:), RDBPI(:,:),        &
                                 XBPO(:), YBPO(:), RDBPO(:,:),        &
                                 ABPI0(:,:), ABPIN(:,:), ABPOS(:,:),  &
around line 427,
 !/
      TYPE OTYPE5
        INTEGER               :: NBI, NBI2, NFBPO, NBO(0:9),          &
                                 NBO2(0:9), NDSL(9), NKI, NTHI
!/MPI        INTEGER               :: NRQBP, NRQBP2
        INTEGER, POINTER      :: IPBPI(:,:), ISBPI(:),                &
                                 IPBPO(:,:), ISBPO(:)
!/MPI        INTEGER, POINTER      :: IRQBP1(:), IRQBP2(:)
        REAL                  :: XFRI, FR1I, TH1I
        REAL, POINTER         :: XBPI(:), YBPI(:), RDBPI(:,:),        &
                                 XBPO(:), YBPO(:), RDBPO(:,:),        &
                                 ABPI0(:,:), ABPIN(:,:), ABPOS(:,:),  &
and around line 1330
      NBO    => OUTPTS(IMOD)%OUT5%NBO
      NBO2   => OUTPTS(IMOD)%OUT5%NBO2
      NDSL   => OUTPTS(IMOD)%OUT5%NDSL
      NKI    => OUTPTS(IMOD)%OUT5%NKI
      NTHI   => OUTPTS(IMOD)%OUT5%NTHI
      XFRI   => OUTPTS(IMOD)%OUT5%XFRI
      FR1I   => OUTPTS(IMOD)%OUT5%FR1I
      TH1I   => OUTPTS(IMOD)%OUT5%TH1I
      FLBPI  => OUTPTS(IMOD)%OUT5%FLBPI
      FLBPO  => OUTPTS(IMOD)%OUT5%FLBPO
      FILER  => OUTPTS(IMOD)%OUT5%FILER
Furthermore, in the file w3iobcmd,ftn, around line 170 add
      USE W3ADATMD, ONLY: CG
      USE W3ODATMD, ONLY: NDSE, NDST, IAPROC, NAPROC, NAPERR, NAPBPT, &
                          NBI, NBI2, NFBPO, NBO, NBO2, NDSL,          &
                          NKI, NTHI, XFRI, FR1I, TH1I,                &
                          IPBPI, ISBPI, XBPI, YBPI, RDBPI,            &
                          IPBPO, ISBPO, XBPO, YBPO, RDBPO,            &
                          ABPI0, ABPIN, ABPOS, FLBPI, FILER, FILEW,   &
and around line 190 remove (red text)
!/ ------------------------------------------------------------------- /
!/ Local parameters
!/
      INTEGER                 :: IFILE, IERR, I, J, NKI, NTHI, IX, IY,&
                                 ISEA, IP, ISP, NPTS, ISOUT, IS, IGRD
!/T1      INTEGER               :: IK, ITH
!/S      INTEGER, SAVE           :: IENT = 0
      REAL                    :: XFRI, FR1I, TH1I
!/T1      REAL                  :: HS, HS0
      REAL, ALLOCATABLE       :: TMPSPC(:,:)
      LOGICAL                 :: FLOK

Identified by Sander Hulst.
Quick fix posted May 31, 2010
Full fix posted July 29. 2010
Back to index.

Declaration error for restart file under MPI :


Problem: A declaration error was found in the auxiliary array VAAUX used to gather all spectra into a single processor when writing the restart file under MPI. On some systems, particularly with different word lengths for integers and reals, this will break the restart file.
Fix: In file w3odatmd.ftn replace lines around 329 and line 422. This needs to be done twice, although at the second location the original line is formatted a little different.


!/MPI        INTEGER, POINTER      :: IRQRS(:), IRQRSS(:), VAAUX(:,:,:)

with

!/MPI        INTEGER, POINTER      :: IRQRS(:), IRQRSS(:)
!/MPI        REAL, POINTER         :: VAAUX(:,:,:)

Identified by Yalin Fan, Posted June 14, 2010
Back to index.

Freezing of wind fields and wave height noise :


Problem: The operational web sites for NOAA WAVEWATCH III have been featuring frozen wind fields and associated noise in wave height fields in the global multi-grid wave model starting around March 2010 and severely aggravated as of July 27, 2010. This operational model uses a modified version of ww3_multi.ftn that allows for computations to progress while wind files are incrementally updated, i.e., this model can thus run side-by-side with the weather model forcing it in a precursor to full dynamic coupling of the weather and wave model. The spurious features in this wave model were traced back to issues with the file system while dynamically updating the input wind files for the wave model. This problem was solved for the Aug. 17 12z model cycle by adding various tests to the scripting and model driver ww3_multi.ftn for this model.
Fix: No fix needed for public distribution version of WAVEWATCH III, because the corresponding version of the wave model driver ww3_multi is not publicly distributed.

Identified by many, Posted Aug. 17, 2010
Back to index.

Initialization of spectra for processing of boundary data :


Problem: For some applications of the ww3_multi.ftn wave model driver and in some hardware/compiler combinations the model crashes in the flux computation routines W3FLX2 or W3FLX3. This problem was traced back to NaNs being introduced in the spectrum at boundary points of nested grids. The NaNs were introduced due to lack of initialization of arrays, but on most systems do no materialize as the non-initialized arrays are multiplied by 0 on first usage.
Fix: In file w3updtmd.ftn add the following lines in subroutine W3UBPT around line 802

!
          BBPI0(:,0) = 0.
          BBPIN(:,0) = 0.
          ABPI0(:,0) = 0.
          ABPIN(:,0) = 0.
!
          DO IBI=1, NBI
            ISEA   = ISBPI(IBI)
Identified by several, Posted Aug. 18, 2010
Back to index.

Managing MAPSTA in inundation grids :


Problem: When using grids with inundation capabilities in the multi-grid model driver ww3_multi, the status map MAPSTA is not properly processed in the wave model routine WMWAVE in wmwavemd,ftn, resulting in dry grid point accidentally being switched back on with zero water depth. For some compilers, this will lead to NaNs in the refraction routines and subsequent model failure, typically somewhere in the source term or flux computations.
Fix: In file wmwavemd.ftn add the following line around line 990 in section 9.b of subroutine WMWAVE. Beware that this line is already present several lines earlier before the call to subroutine W3WAVE.

                    IF ( FLGHG1 .AND. .NOT.FLGHG2 .AND.               &
                         GRDHGH(I,0).GT.0 ) THEN
                        MAPST2 = MAPST2 + 8*MAPMSK
                        MAPSTA = ABS(MAPSTA)
                        DO IX=1, NX
                          DO IY=1, NY
                            IF ( MAPST2(IY,IX) .GT. 0 )               &
Posted Aug. 27, 2010
Back to index.

Possible issue with depth treatment in W3SRCE :


Possible problem: In this routine, the depth is allowed to go to 0, whereas in the propagation routines and all other computations, the depth is not allowed to be smaller than DMIN. Although we have not had any reports of this giving problems with model runs, we suggest the following preventive correction to the code.
Fix: In the subroutine W3SRCE in the file w3srcemd.ftn apply the following changes. In the header of the routine, replace

      SUBROUTINE W3SRCE ( IX, IY, IMOD, SPEC, ALPHA, WN1, CG1, DEPTH, &
with
      SUBROUTINE W3SRCE ( IX, IY, IMOD, SPEC, ALPHA, WN1, CG1, D_INP, &
In the use statement for the grid data structure replace
      USE W3GDATMD, ONLY: NK, NTH, NSPEC, SIG, TH, DTMAX, DTMIN,      &
                          FACTI1, FACTI2, FACSD, FACHFA, FACP,        &
with
      USE W3GDATMD, ONLY: NK, NTH, NSPEC, SIG, TH, DMIN, DTMAX,       &
                          DTMIN, FACTI1, FACTI2, FACSD, FACHFA, FACP, &
In the declaration of the parameter list, replace
      REAL, INTENT(IN)        :: WN1(NK), CG1(NK), DEPTH, U10ABS,     &
with
      REAL, INTENT(IN)        :: WN1(NK), CG1(NK), D_INP, U10ABS,     &
In the declaration of variables internal to the subroutine add
      REAL                    :: DTTOT, FHIGH, DT, AFILT, DAMAX, AFAC,&
                                 HDT, ZWND, FP, DEPTH
Finally, around line 367, add the lines
!
!/T      FLTEST = .TRUE.
!
      DEPTH  = MAX ( DMIN , D_INP )
!
!/LN0      VSLN = 0.
!/SEED      VSLN = 0.
!/ST0      VSIN = 0.
Posted Aug. 27, 2010
Back to index.

Bug in combining wind seas in the partitioning routine :

Problem: In w3partmd.f90, the subroutine PTMEAN returns NP, the number of the NP_MAX partitions for which HS >= HSPMIN, but it doesn't return information about *which* partitions were selected. The code that combines wind-seas assumes that the NP partitions selected are the *first* NP of the NP_MAX partitions, which is not necessarily the case.

Fix:

  1. In subroutine W3PART, declare a new local integer array:
        INTEGER :: PMAP(DIMXP)
    
  2. Modify calls to PTMEAN from
         CALL PTMEAN(NP_MAX,IMO,ZP,DEPTH,UABS,UDIR,WN,NP,XP,DIMXP)
    
    to
         CALL PTMEAN(NP_MAX,IMO,ZP,DEPTH,UABS,UDIR,WN,NP,XP,DIMXP,PMAP)
    
  3. Replace lines 236-242:
         IPW    = INDEX(1)
         DO IP=2, NWS
           IPT    = INDEX(IP)
           DO ISP=1, NSPEC
             IF ( IMO(ISP) .EQ. IPT ) IMO(ISP) = IPW
             END DO
           END DO  
    
    with
         IPW    = PMAP(INDEX(1))
         DO IP=2, NWS
           IPT    = PMAP(INDEX(IP))
           DO ISP=1, NSPEC
             IF ( IMO(ISP) .EQ. IPT ) IMO(ISP) = IPW
             END DO
           END DO
    
  4. In subroutine PTMEAN, add an OUT parameter PMAP as above:
         INTEGER, INTENT(OUT) :: PMAP(DIMXP)
    
  5. and modify lines 1064-1065:
         NPO       = NPO + 1
         XP(1,NPO) = HS
    
    to
         NPO       = NPO + 1
         IF (IP.GT.0) PMAP(NPO) = IP
         XP(1,NPO) = HS
    
Identified by Mark Szyszka
Posted December 2, 2010
Back to index.

Bug: Out of bounds error message for point data :

Problem: When point output for the global grid (grid that wraps around) lies between NX and 1, we get an out of bounds error in the initialization routine.

Fix: On line 313 in wmiopomd.ftn add the following line

      IF ( GLOBAL .AND. IX1 .EQ. NX ) IXS = 1-NX 


Posted January 17, 2011
Back to index.