#include "w3macros.h" !/ ------------------------------------------------------------------- / MODULE W3SRC2MD !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | FORTRAN 90 | !/ | Last update : 29-May-2009 | !/ +-----------------------------------+ !/ !/ 04-Feb-2000 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ 21-Feb-2004 : Multiple model version. ( version 3.06 ) !/ 03-Jul-2006 : Extract stress computation. ( version 3.09 ) !/ 13-Apr-2007 : EMEAN in W3SPR2 par list. ( version 3.11 ) !/ 29-May-2009 : Preparing distribution version. ( version 3.14 ) !/ !/ Copyright 2009 National Weather Service (NWS), !/ National Oceanic and Atmospheric Administration. All rights !/ reserved. WAVEWATCH III is a trademark of the NWS. !/ No unauthorized use without permission. !/ ! 1. Purpose : ! ! Tolman and Chalikov (1996) input and dissipation source terms. ! Bundled with interpolation tables. ! ! 2. Variables and types : ! ! Interpolation tables : ! ! Name Type Scope Description ! ---------------------------------------------------------------- ! NRSIGA I.P. Public Array dimension (SIGA). ! NRDRAG I.P. Public Array dimension (drag coefficient). ! SIGAMX R.P. Public Maximum nondiensional frequency. ! DRAGMX R.P. Public Maximum drag coefficient. ! DSIGA Real Public Table increment. ! DDRAG Real Public Id. ! BETATB R.A. Public Interpolation table. ! ---------------------------------------------------------------- ! ! 3. Subroutines and functions : ! ! Name Type Scope Description ! ---------------------------------------------------------------- ! W3SPR2 Subr. Public Mean parameters from spectrum. ! W3SIN2 Subr. Public Input source term. ! W3SDS2 Subr. Public Dissipation source term. ! INPTAB Subr. Public Interpolation table for wind-wave ! interaction parameter. ! W3BETA R.F. INPTAB Id. function. ! ---------------------------------------------------------------- ! ! 4. Subroutines and functions used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ( !/S ) ! PRT2DS Subr. W3ARRYMD Print plot of spectra. ( !/T0 ) ! OUTMAT Subr. W3WRRYMD Print out 2D matrix. ( !/T1 ) ! ... Data W3DISPMD Interpolation tables to solve ! dispersion relation. ! ---------------------------------------------------------------- ! ! 5. Remarks : ! ! 6. Switches : ! ! !/S Enable subroutine tracing. ! !/T(n) Test output, see subroutines. ! ! 7. Source code : ! !/ ------------------------------------------------------------------- / PUBLIC !/ !/ Interpolation table !/ INTEGER, PARAMETER, PRIVATE :: NRSIGA = 400 INTEGER, PARAMETER, PRIVATE :: NRDRAG = 20 REAL, PARAMETER, PRIVATE :: SIGAMX = 40. REAL, PARAMETER, PRIVATE :: DRAGMX = 1.E-2 ! REAL, PRIVATE :: DSIGA, DDRAG, & BETATB(-NRSIGA:NRSIGA+1,NRDRAG+1) !/ CONTAINS !/ ------------------------------------------------------------------- / SUBROUTINE W3SPR2 (A, CG, WN, DEPTH, FPI, U, USTAR, & EMEAN, FMEAN, WNMEAN, AMAX, ALFA, FP ) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | D.Chalikov | !/ | FORTRAN 90 | !/ | Last update : 13-Apr-2007 | !/ +-----------------------------------+ !/ !/ 06-Dec-1996 : Final version 1.18 / FORTRAN 77 version. !/ 16-Nov-1999 : Add itteration to section 5. for removal of W3APR2. !/ 04-Feb-2000 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ 21-Dec-2004 : Multiple model version. ( version 3.06 ) !/ 03-Jul-2006 : Extract stress computation. ( version 3.09 ) !/ 13-Apr-2007 : EMEAN in parameter list. ( version 3.11 ) ! ! 1. Purpose : ! ! Calculate mean wave parameters for the use in the source term ! routines. (Tolman and Chalikov) ! ! 2. Method : ! ! See source term routines. ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! A R.A. I Action density spectrum. ! CG R.A. I Group velocities. ! WN R.A. I Wavenumbers. ! DEPTH Real I Water depth. ! FPI Real I Peak input frequency. ! U Real I Wind speed. ! USTAR Real I Friction velocity. ! EMEAN Real O Total energy (variance). ! FMEAN Real O Mean frequency. ! WNMEAN Real O Mean wavenumber. ! AMAX Real O Maximum of action spectrum. ! ALFA R.A. O Phillips' constant. ! FP Real O Peak frequency. ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3SRCE Subr. W3SRCEMD Source term integration. ! W3EXPO Subr. N/A Point output post-processor. ! GXEXPO Subr. N/A GrADS point output post-processor. ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/S Enable subroutine tracing. ! !/T Enable test output. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE CONSTANTS USE W3GDATMD, ONLY: NK, NTH, DTH, SIG, DDEN, FTE, FTF, FTWN, & NITTIN, ZWIND, CINXSI !/T USE W3ODATMD, ONLY: NDST !/S USE W3SERVMD, ONLY: STRACE USE W3DISPMD, ONLY: NAR1D, DFAC, N1MAX, ECG1, EWN1, DSIE !/ IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ REAL, INTENT(IN) :: A(NTH,NK), CG(NK), WN(NK), DEPTH, & FPI, U, USTAR REAL, INTENT(OUT) :: EMEAN, FMEAN, WNMEAN, AMAX, & ALFA(NK), FP !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ INTEGER :: IK, ITH, I1, ITT !/S INTEGER, SAVE :: IENT = 0 REAL :: EBAND, FPISTR, EB(NK), UST !/ !/ ------------------------------------------------------------------- / !/ !/S CALL STRACE (IENT, 'W3SPR2') ! UST = MAX ( 0.0001 , USTAR ) ! EMEAN = 0. FMEAN = 0. WNMEAN = 0. AMAX = 0. ! ! 1. Integral over directions and maximum --------------------------- * ! DO IK=1, NK EB(IK) = 0. DO ITH=1, NTH EB(IK) = EB(IK) + A(ITH,IK) AMAX = MAX ( AMAX , A(ITH,IK) ) END DO END DO ! ! 2. Integrate over directions -------------------------------------- * ! DO IK=1, NK ALFA(IK) = 2. * DTH * SIG(IK) * EB(IK) * WN(IK)**3 EB(IK) = EB(IK) * DDEN(IK) / CG(IK) EMEAN = EMEAN + EB(IK) FMEAN = FMEAN + EB(IK) / SIG(IK) WNMEAN = WNMEAN + EB(IK) / SQRT(WN(IK)) END DO ! ! 3. Add tail beyond discrete spectrum and get mean pars ------------ * ! ( DTH * SIG absorbed in FTxx ) ! EBAND = EB(NK) / DDEN(NK) EMEAN = EMEAN + EBAND * FTE FMEAN = FMEAN + EBAND * FTF WNMEAN = WNMEAN + EBAND * FTWN ! FMEAN = TPIINV * EMEAN / MAX ( 1.E-7 , FMEAN ) WNMEAN = ( EMEAN / MAX ( 1.E-7 , WNMEAN ) )**2 ! ! 4. Estimate peak frequency from FPI ------------------------------- * ! FPISTR = MAX ( 0.008 , FPI * UST / GRAV ) FP = ( 3.6E-4 + 0.92*FPISTR - 6.3E-10/FPISTR**3 )/UST*GRAV FP = FP * TPIINV ! RETURN !/ !/ End of W3SPR2 ----------------------------------------------------- / !/ END SUBROUTINE W3SPR2 !/ ------------------------------------------------------------------- / SUBROUTINE W3SIN2 ( A, CG, K, U, UDIR, CD, Z0, FPI, S, D ) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | D.Chalikov | !/ | FORTRAN 90 | !/ | Last update : 21-Feb-2004 | !/ +-----------------------------------+ !/ !/ 14-Jan-1997 : Final FORTRAN 77 ( version 1.18 ) !/ 04-Feb-2000 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ 21-Feb-2004 : Multiple model version. ( version 3.06 ) !/ ! 1. Purpose : ! ! Calculate input source term. ! ! 2. Method : ! ! Tolman and Chalikov (1996), see manual. ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! A R.A. I Action density spectrum (1-D). ! CG R.A. I Group velocities for k-axis of spectrum. ! K R.A. I Wavenumber for entire spectrum (1-D). ! U Real I Wind speed at reference height. ! UDIR Real I Direction of U. ! CD Real I Drag coefficient at wind level ZWIND. ! Z0 Real I Corresponding z0. ! FPI R.A. O Input 'peak' frequency. ! S R.A. O Source term (1-D version). ! D R.A. O Diagonal term of derivative (1-D version). ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! PRT2DS Subr. W3ARRYMD Print plot of spectra. ! OUTMAT Subr. W3WRRYMD Print out 2D matrix. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3SRCE Subr. W3SRCEMD Source term integration. ! W3EXPO Subr. N/A Point output post-processor. ! GXEXPO Subr. N/A GrADS point output post-processor. ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! 7. Remarks : ! ! - Actual height of wind speed does not need to be 10 m, but is ! given by ZWIND. ! - Abs(cos) > 0.0087 to asure continuity in beta. Corresponds ! to shift of up to half a degree. ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/S Enable subroutine tracing. ! !/T Enable general test output. ! !/T0 Print arrays. ! !/T1 Calculation of diagonal within spectrum ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE CONSTANTS USE W3GDATMD, ONLY: NK, NTH, NSPEC, XFR, DDEN, SIG, SIG2, & ESIN, ECOS, FTE, FTTR, FPIMIN, ZWIND, & FACTI1, FACTI2, FSWELL !/T USE W3ODATMD, ONLY: NDST !/S USE W3SERVMD, ONLY: STRACE !/T0 USE W3ARRYMD, ONLY: PRT2DS !/T1 USE W3ARRYMD, ONLY: OUTMAT !/ IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ REAL, INTENT(IN) :: A(NSPEC), CG(NK), K(NSPEC), U, UDIR, & CD, Z0 REAL, INTENT(OUT) :: S(NSPEC), D(NSPEC), FPI !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ INTEGER :: IS, IK, IOMA, ICL, NKFILT, NKFIL2 !/S INTEGER, SAVE :: IENT = 0 !/T0 INTEGER ITH REAL :: COSU, SINU, COSFAC, LAMBDA, ULAM, & CLAM, OMA, M0, M1, RD1, RD2, BETA, & FACLN1, FACLN2, USTAR, TRANS, FPISTR,& FP1STR, FP1, SIN1A(NK) REAL, PARAMETER :: TRANSF = 0.75 REAL, PARAMETER :: PEAKFC = 0.8 !/T0 REAL :: DOUT(NK,NTH) !/ !/ ------------------------------------------------------------------- / !/ !/S CALL STRACE (IENT, 'W3SIN2') ! !/T WRITE (NDST,9000) DSIGA, DDRAG, U, UDIR*RADE, CD, Z0 ! ! 1. Preparations ! COSU = COS(UDIR) SINU = SIN(UDIR) ! ! 2. Loop over spectrum ! !/T2 WRITE (NDST,9020) ! FACLN1 = U / LOG(ZWIND/Z0) FACLN2 = LOG(Z0) ! DO IS=1, NSPEC COSFAC = ECOS(IS)*COSU + ESIN(IS)*SINU COSFAC = SIGN ( MAX ( 0.0087 , ABS(COSFAC) ) , COSFAC ) LAMBDA = TPI / ( K(IS) * ABS(COSFAC) ) ULAM = FACLN1 * ( LOG(LAMBDA) - FACLN2 ) CLAM = CD * ( U / ULAM )**2 OMA = K(IS) * ULAM * COSFAC / SIG2(IS) IOMA = INT ( OMA/DSIGA ) + & MIN ( 0 , INT ( SIGN ( -1.1 , OMA ) ) ) ICL = INT ( CLAM/DDRAG ) RD1 = OMA/DSIGA - REAL(IOMA) RD2 = CLAM/DDRAG - REAL(ICL) IOMA = MAX ( -NRSIGA , MIN ( NRSIGA , IOMA ) ) ICL = MAX ( 1 , MIN ( NRDRAG , ICL ) ) BETA = (1.-RD1) * (1.-RD2) * BETATB( IOMA , ICL ) & + RD1 * (1.-RD2) * BETATB(IOMA+1, ICL ) & + (1.-RD1) * RD2 * BETATB( IOMA ,ICL+1) & + RD1 * RD2 * BETATB(IOMA+1,ICL+1) D(IS) = BETA * SIG2(IS) S(IS) = A(IS) * D(IS) !/T2 WRITE (NDST,9021) IS, COSFAC, LAMBDA, ULAM, CLAM*1.E3, & !/T2 OMA, BETA*1.E4 END DO ! ! 3. Calculate FPI ! DO IK=1, NK SIN1A(IK) = 0. DO IS=(IK-1)*NTH+1, IK*NTH SIN1A(IK) = SIN1A(IK) + MAX ( 0. , S(IS) ) END DO END DO ! M0 = 0. M1 = 0. DO IK=1, NK SIN1A(IK) = SIN1A(IK) * DDEN(IK) / ( CG(IK) * SIG(IK)**3 ) M0 = M0 + SIN1A(IK) M1 = M1 + SIN1A(IK)/SIG(IK) END DO ! SIN1A(NK) = SIN1A(NK) / DDEN(NK) M0 = M0 + SIN1A(NK) * FTE M1 = M1 + SIN1A(NK) * FTTR IF ( M1 .LT. 1E-20 ) THEN FPI = XFR * SIG(NK) ELSE FPI = M0 / M1 END IF ! ! 4. Filter for swell ! USTAR = U * SQRT(CD) FPISTR = MAX ( FPIMIN , FPI * USTAR / GRAV ) FP1STR = 3.6E-4 + 0.92*FPISTR - 6.3E-10/FPISTR**3 FP1 = PEAKFC * FP1STR * GRAV / USTAR ! NKFILT = MIN ( NK , INT(FACTI2+FACTI1*LOG(FP1)) ) NKFIL2 = MIN ( NK , INT(FACTI2+FACTI1*LOG(TRANSF*FP1)) ) NKFIL2 = MAX ( 0 , NKFIL2 ) ! DO IS=1, NKFIL2*NTH D(IS) = MAX ( D(IS) , FSWELL*D(IS) ) S(IS) = A(IS) * D(IS) END DO ! DO IK=NKFIL2+1, NKFILT TRANS = ( SIG(IK)/FP1 - TRANSF ) / (1.-TRANSF) DO IS=(IK-1)*NTH+1, IK*NTH D(IS) = (1.-TRANS)*MAX(D(IS),FSWELL*D(IS)) + TRANS*D(IS) S(IS) = A(IS) * D(IS) END DO END DO ! ! ... Test output of arrays ! !/T0 DO IK=1, NK !/T0 DO ITH=1, NTH !/T0 DOUT(IK,ITH) = D(ITH+(IK-1)*NTH) !/T0 END DO !/T0 END DO ! !/T0 CALL PRT2DS (NDST, NK, NK, NTH, DOUT, SIG(1:), ' ', 1., & !/T0 0.0, 0.001, 'Diag Sin', ' ', 'NONAME') ! !/T1 CALL OUTMAT (NDST, D, NTH, NTH, NK, 'diag Sin') ! RETURN ! ! Formats ! !/T 9000 FORMAT (' TEST W3SIN2 : DSIGA,DDRAG,U,UDIR,CD,Z0(IN) : '/ & !/T ' ',F8.4,F9.6,F7.2,F6.1,F8.5,F8.5) ! !/T2 9020 FORMAT (' TEST W3SIN2 : IS, COS, LAMBDA, ULAM, CLAM*1E3, ', & !/T2 'OMA, BETA*1E4') !/T2 9021 FORMAT (6X,I6,F7.2,1X,F6.1,2(1X,F5.2),2(1X,F6.2)) !/ !/ End of W3SIN2 ----------------------------------------------------- / !/ END SUBROUTINE W3SIN2 !/ ------------------------------------------------------------------- / SUBROUTINE W3SDS2 (A, CG, K, FPI, USTAR, ALFA, S, D) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | FORTRAN 90 | !/ | Last update : 21-Feb-2004 | !/ +-----------------------------------+ !/ !/ 12-Jun-1996 : Final FORTRAN 77 ( version 1.18 ) !/ 04-Feb-2000 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ 23-Apr-2002 : Erick Rogers' fix ( version 2.19 ) !/ 21-Feb-2004 : Multiple model version. ( version 3.06 ) !/ ! 1. Purpose : ! ! Calculate whitecapping source term and diagonal term of der. ! ! 2. Method : ! ! Tolman and Chalikov (1995). ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! A R.A. I Input action density spectrum. ! CG R.A. I Group velocity array. ! K R.A. I Wavenumber array. ! FPI Real I 'Peak frequency' of input (rad/s). ! USTAR Real I Friction velocity (m/s). ! ALFA R.A. I Phillips' constant. ! S R.A. O Source term (1-D version). ! D R.A. O Diagonal term of derivative (1-D version). ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! PRT2DS Subr. W3ARRYMD Print plot of spectra. ! OUTMAT Subr. W3WRRYMD Print out 2D matrix. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3SRCE Subr. W3SRCEMD Source term integration. ! W3EXPO Subr. N/A Point output post-processor. ! GXEXPO Subr. N/A GrADS point output post-processor. ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/S Enable subroutine tracing. ! !/T Enable general test output. ! !/T0 Print arrays. ! !/T1 Print filter and constituents. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE CONSTANTS USE W3GDATMD, ONLY: NK, NTH, SIG, DDEN, DTH, FTE, FPIMIN, & FACTI1, FACTI2, XF1, XF2, XFH, SDSALN, & CDSA0, CDSA1, CDSA2, CDSB0, CDSB1, CDSB2, & CDSB3 !/T USE W3ODATMD, ONLY: NDST !/S USE W3SERVMD, ONLY: STRACE !/T0 USE W3ARRYMD, ONLY: PRT2DS !/T1 USE W3ARRYMD, ONLY: OUTMAT ! IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ REAL, INTENT(IN) :: A(NTH,NK), CG(NK), K(NK), FPI, & USTAR, ALFA(NK) REAL, INTENT(OUT) :: S(NTH,NK), D(NTH,NK) !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ INTEGER :: IK, ITH, IKHW !/S INTEGER, SAVE :: IENT = 0 REAL :: FHW, XHW, FPIT, PHI, AF1, AF2, & AFILT, BFILT, CDIST, FILT, POW, & CDISH, CDISP, HW, EHIGH, EBD(NK) !/T REAL POWMAX !/T0 REAL DOUT(NK,NTH) !/ !/ ------------------------------------------------------------------- / !/ !/S CALL STRACE (IENT, 'W3SDS2') ! !/T WRITE (NDST,9000) FPI, USTAR ! ! 1. Preparations ! 1.a HW ! FHW = XFH*FPI XHW = FACTI2 + FACTI1*LOG(FHW) IKHW = MIN ( NK , INT ( XHW + 0.5 ) ) DO IK=IKHW, NK EBD(IK) = 0. DO ITH=1, NTH EBD(IK) = EBD(IK) + A(ITH,IK) END DO END DO ! IF ( FHW .LT. SIG(NK+1) ) THEN XHW = 1. - MOD ( XHW + 0.5 , 1. ) IF ( IKHW .EQ. NK ) XHW = MAX ( 0. , XHW - 0.5 ) HW = XHW * EBD(IKHW)*DDEN(IKHW)/CG(IKHW) DO IK=IKHW+1, NK HW = HW + EBD(IK)*DDEN(IK)/CG(IK) END DO HW = 4. * SQRT ( HW + EBD(NK)/CG(NK)*FTE ) ELSE EHIGH = EBD(NK)/CG(NK) * SIG(NK)*DTH * (SIG(NK)/FHW)**5 HW = 4. * SQRT ( 0.25 * FHW * EHIGH ) END IF ! ! 1.b PHI ! FPIT = MAX ( FPIMIN , FPI*TPIINV*USTAR/GRAV ) PHI = CDSB0 + CDSB1*FPIT + CDSB2/FPIT**CDSB3 ! ! 1.c Set-up filter ! AF2 = XF2*FPI AF1 = XF1*FPI BFILT = 1. / ( AF2 - AF1 ) AFILT = - BFILT * AF1 ! ! 1.d Constants ! CDIST = - 2. * USTAR * HW * PHI CDISH = G2PI3I * USTAR**2 CDISP = G1PI1I * USTAR ! ! 2. Combined diagonal factor ! !/T2 WRITE (NDST,9020) !/T POWMAX = 0. DO IK=1, NK FILT = MIN ( 1., MAX ( 0. , AFILT + BFILT*SIG(IK) )) POW = MIN ( 25. , CDSA1 / ( CDISP*SIG(IK) )**CDSA2 ) IF ( FILT .GT. 0. ) THEN D(1,IK) = (1.-FILT) * CDIST * K(IK)**2 & - FILT * CDSA0 * CDISH * SIG(IK)**3 & * (ALFA(IK)/SDSALN)**POW ELSE D(1,IK) = (1.-FILT) * CDIST * K(IK)**2 END IF !/T POWMAX = MAX(POW*FILT,POWMAX) !/T2 WRITE (NDST,9021) IK, FILT, ALFA(IK)/SDSALN, & !/T2 CDIST*PHI*K(IK)**2, CDSA0*CDISH*SIG(IK)**3 & !/T2 * (ALFA(IK)/SDSALN)**POW, D(1,IK) END DO ! !/T WRITE (NDST,9010) AF1, AF2, AFILT, BFILT, POWMAX ! ! 3. 2-D diagonal array ! DO IK=1, NK DO ITH=2, NTH D(ITH,IK) = D(1,IK) END DO END DO ! S = D * A ! ! ... Test output of arrays ! !/T0 DO IK=1, NK !/T0 DO ITH=1, NTH !/T0 DOUT(IK,ITH) = D(ITH,IK) !/T0 END DO !/T0 END DO ! !/T0 CALL PRT2DS (NDST, NK, NK, NTH, DOUT, SIG(1:), ' ', 1., & !/T0 0.0, 0.001, 'Diag Sds', ' ', 'NONAME') ! !/T1 CALL OUTMAT (NDST, D, NTH, NTH, NK, 'diag Sds') ! RETURN ! ! Formats ! !/T 9000 FORMAT (' TEST W3SDS2 : FPI, USTAR : ',2F8.3) !/T 9010 FORMAT (' TEST W3SDS2 : AF1-2, A-BFILT, PMAX : ',4F7.3,E10.3) !/T2 9020 FORMAT (' TEST W3SDS2 : IK, FILT, ALFA, DDST, DDSH, DDS') !/T2 9021 FORMAT (' ',I6,2F7.3,3E11.3) !/ !/ End of W3SDS2 ----------------------------------------------------- / !/ END SUBROUTINE W3SDS2 !/ ------------------------------------------------------------------- / SUBROUTINE INPTAB !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | FORTRAN 90 | !/ | Last update : 21-Feb-2004 | !/ +-----------------------------------+ !/ !/ 03-Jun-1996 : Final version 1.18 / FORTRAN 77 version. !/ 06-Dec-1999 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ 21-Feb-2004 : Multiple model version. ( version 3.06 ) !/ ! 1. Purpose : ! ! Generate an interpolation table for the air-sea interaction ! parameter of Chalikov and Belevich (1993). ! ! 2. Method : ! ! The size of the table is set in parameter statements, the range ! is set by the input parameters of this routine. The first counter ! of the table corresponds to the nondimensional frequency ! ! SIGMA Ul ! SIGA = ---------- COS ( THETA - THETA ) (1) ! g wind ! ! The second counter of the table represents the drag coefficient. ! The maximum values of both parameters are passed to the routine ! through the parameter list. ! ! 3. Parameters : ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! STRACE Subr. W3SERVMD Subroutine tracing. ! W3BETA Func. Internal Function to calculate the ! interaction parameter. ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3IOGR Subr. W3IOGRMD Model definition IO routine. ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/S Enable subroutine tracing. ! !/T Enable test output. ! !/T0 Print table. ! !/T1 Estimate maximum errors. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE CONSTANTS USE W3ODATMD, ONLY: NDST !/S USE W3SERVMD, ONLY: STRACE !/ IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ INTEGER :: ISIGA, IDRAG !/S INTEGER, SAVE :: IENT = 0 !/T0 INTEGER :: I1 !/T1 INTEGER :: IE1 REAL :: SIGA, DRAG !/T0 REAL :: BMIN, BMAX !/T1 REAL :: ENORM, ERR(NRDRAG) !/ !/ ------------------------------------------------------------------- / !/ !/S CALL STRACE (IENT, 'INPTAB') ! ! 1. Determine range and increments of table ------------------------ * ! DSIGA = SIGAMX / REAL(NRSIGA) DDRAG = DRAGMX / REAL(NRDRAG) ! !/T WRITE (NDST,9000) SIGAMX, DSIGA, DRAGMX, DDRAG ! ! 2. Fill table ----------------------------------------------------- * ! DO ISIGA=-NRSIGA,NRSIGA+1 SIGA = REAL(ISIGA) * DSIGA DO IDRAG=1, NRDRAG+1 DRAG = REAL(IDRAG) * DDRAG BETATB(ISIGA,IDRAG) = W3BETA ( SIGA, DRAG , NDST ) END DO END DO ! ! 3. Test output ---------------------------------------------------- * ! !/T0 WRITE (NDST,9010) !/T0 I1 = MIN (35,NRDRAG) !/T0 DO ISIGA=-NRSIGA,NRSIGA !/T0 SIGA = REAL(ISIGA) * DSIGA !/T0 BMIN = 0. !/T0 BMAX = 0. !/T0 DO IDRAG=1, NRDRAG !/T0 BMIN = MIN ( BMIN , BETATB(ISIGA,IDRAG) ) !/T0 BMAX = MAX ( BMAX , BETATB(ISIGA,IDRAG) ) !/T0 END DO !/T0 BMAX = MAX ( BMAX , -BMIN ) !/T0 WRITE (NDST,9011) ISIGA, SIGA, BMAX, & !/T0 (NINT(BETATB(ISIGA,IDRAG)/BMAX*100.),IDRAG=1,I1) !/T0 IF (I1.LT.NRDRAG) WRITE (NDST,9012) & !/T0 (NINT(BETATB(ISIGA,IDRAG)/BMAX*100.),IDRAG=I1+1,NRDRAG) !/T0 END DO ! !/T1 WRITE (NDST,9020) !/T1 IE1 = MIN (30,NRDRAG-1) !/T1 ENORM = 1000. / ABS(BETATB(0,NRDRAG)) !/T1 DO ISIGA=-NRSIGA,NRSIGA !/T1 SIGA = REAL(ISIGA) * DSIGA !/T1 IF ( ABS(SIGA) .LT. 5.01 ) THEN !/T1 DO IDRAG=1, NRDRAG-1 !/T1 DRAG = ( REAL(IDRAG) + 0.5 ) * DDRAG !/T1 ERR(IDRAG) = - W3BETA (SIGA,DRAG,NDST) + 0.5 * & !/T1 ( BETATB(ISIGA,IDRAG) + BETATB(ISIGA,IDRAG+1) ) !/T1 END DO !/T1 WRITE (NDST,9021) ISIGA, SIGA, & !/T1 (NINT(ENORM*ERR(IDRAG)),IDRAG=1,IE1) !/T1 IF (IE1.LT.NRDRAG-1) WRITE (NDST,9022) & !/T1 (NINT(ENORM*ERR(IDRAG)),IDRAG=IE1+1,NRDRAG-1) !/T1 ENDIF !/T1 END DO ! !/T1 WRITE (NDST,9030) !/T1 IE1 = MIN (30,NRDRAG) !/T1 ENORM = 1000. / ABS(BETATB(0,NRDRAG)) !/T1 DO ISIGA=-NRSIGA,NRSIGA-1 !/T1 SIGA = ( REAL(ISIGA) + 0.5 ) * DSIGA !/T1 IF ( ABS(SIGA) .LT. 5.01 ) THEN !/T1 DO IDRAG=1, NRDRAG !/T1 DRAG = REAL(IDRAG) * DDRAG !/T1 ERR(IDRAG) = - W3BETA (SIGA,DRAG,NDST) + 0.5 * & !/T1 ( BETATB(ISIGA,IDRAG) + BETATB(ISIGA+1,IDRAG) ) !/T1 END DO !/T1 WRITE (NDST,9031) ISIGA, SIGA, & !/T1 (NINT(ENORM*ERR(IDRAG)),IDRAG=1,IE1) !/T1 IF (IE1.LT.NRDRAG) WRITE (NDST,9032) & !/T1 (NINT(ENORM*ERR(IDRAG)),IDRAG=IE1+1,NRDRAG) !/T1 ENDIF !/T1 END DO ! RETURN ! ! Formats ! !/T 9000 FORMAT ( ' TEST INPTAB : SIGAMX, DSIGA : ',F6.2,F8.2/ & !/T ' DRAGMX, DDRAG : ',F8.4,F9.5) ! !/T0 9010 FORMAT (/' TEST INPTAB : TABLE, NORMALIZED WITH ', & !/T0 'BETATB(ISIGA,NRDRAG)'/ & !/T0 ' ISIGA, SIGA, BETA_MAX, TABLE (x100)') !/T0 9011 FORMAT (1X,I4,F7.2,F6.4,1X,35I3) !/T0 9012 FORMAT (19X,35I3) ! !/T1 9020 FORMAT (/' TEST INPTAB : ERROR DUE TO DRAG, NORMALIZED ', & !/T1 'WITH BETATB(ISIGA,NRDRAG)'/ & !/T1 ' ISIGA, SIGA, TABLE (x1000)') !/T1 9021 FORMAT (1X,I4,F7.2,35I3) !/T1 9022 FORMAT (12X,35I3) ! !/T1 9030 FORMAT (/' TEST INPTAB : ERROR DUE TO SIGA, NORMALIZED WITH ', & !/T1 'BETATB(ISIGA,NRDRAG)'/ & !/T1 ' ISIGA, SIGA, TABLE (x1000)') !/T1 9031 FORMAT (1X,I4,F7.2,35I3) !/T1 9032 FORMAT (12X,35I3) !/ !/ Internal function W3BETA !/ CONTAINS !/ ------------------------------------------------------------------- / REAL FUNCTION W3BETA ( OMA , CL , NDST ) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | H. L. Tolman | !/ | D.Chalikov | !/ | FORTRAN 90 | !/ | Last update : 21-Feb-2004 | !/ +-----------------------------------+ !/ !/ 06-Dec-1996 : Final version 1.18 / FORTRAN 77 version. !/ 06-Dec-1999 : Upgrade to FORTRAN 90 ( version 2.00 ) !/ 21-Feb-2004 : Multiple model version. ( version 3.06 ) !/ ! 1. Purpose : ! ! Calculate wind-wave interaction parameter beta. ! ! 2. Method : ! ! Chalikov and Belevich (1992), see also manual. ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! W3BETA Real O Wind-wave interaction parameter multiplied ! by density ratio. ! OMA Real I Non-dimensional apparent frequency. ! ! OMA = OMEGA | U | cos(theta-theta ) / g ! l w ! ! CL Real I Drag coefficient at height l ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! 5. Called by : ! ! 6. Error messages : ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/S Enable subroutine tracing. ! !/T0 Enable test output. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list !/ INTEGER, INTENT(IN) :: NDST REAL, INTENT(IN) :: OMA, CL !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ !/S INTEGER, SAVE :: IENT = 0 REAL :: OM1, OM2, A0, A1, A2, A3, A4, A5, & A6, A7, A8, A9, A10 !/ !/ ------------------------------------------------------------------- / !/ !/S CALL STRACE (IENT, 'W3BETA') ! !/T0 WRITE (NDST,9000) OMA, CL ! ! calculate Omegas ! OM1 = 1.075 + 75.*CL OM2 = 1.2 + 300.*CL ! ! calculate factors a ! A1 = 0.25 + 395.*CL A2 = 0.35 + 150.*CL A4 = 0.3 + 300.*CL A9 = 0.35 + 240.*CL A10 = -0.06 + 470.*CL ! A5 = A4 * OM1 A0 = 0.25 * A5**2 / A4 A3 = (A0-A2-A1) / (A0+A4+A5) A6 = A0 * (1.-A3) A7 = (A9*(OM2-1)**2+A10) / (OM2-OM1) A8 = A7 * OM1 ! !/T0 WRITE (NDST,9001) OM1, OM2 !/T0 WRITE (NDST,9002) A0, A1, A2, A3, A4, A5, A6, A7, A8, A9, A10 ! ! calculate beta * 1.e4 ! IF ( OMA .LT. -1. ) THEN W3BETA = -A1 * OMA**2 - A2 ELSE IF (OMA .LT. OM1/2.) THEN W3BETA = A3 * OMA * ( A4 * OMA - A5 ) - A6 ELSE IF (OMA .LT. OM1) THEN W3BETA = OMA * ( A4 * OMA - A5 ) ELSE IF (OMA .LT. OM2) THEN W3BETA = A7 * OMA - A8 ELSE W3BETA = A9 * (OMA-1.)**2 + A10 END IF ! ! beta * dwat / dair ! W3BETA = W3BETA * 1.E-4 !/T0 WRITE (NDST,9003) W3BETA ! RETURN ! ! Formats ! !/T0 9000 FORMAT ( ' TEST W3BETA : INPUT : ',2E10.3) !/T0 9001 FORMAT ( ' TEST W3BETA : OM1-2 : ',2E10.3) !/T0 9002 FORMAT ( ' TEST W3BETA : A0-10 : ',5E10.3/ & !/T0 ' ',6E10.3) !/T0 9003 FORMAT ( ' TEST W3BETA : BETA : ',E10.3) !/ !/ End of W3BETA ----------------------------------------------------- / !/ END FUNCTION W3BETA !/ !/ End of INPTAB ----------------------------------------------------- / !/ END SUBROUTINE INPTAB !/ !/ End of module W3SRC2MD -------------------------------------------- / !/ END MODULE W3SRC2MD