WAVEWATCH III  beta 0.0.1
Functions/Subroutines
w3sbt9md Module Reference

Contains routines for computing dissipation by viscous fluid mud using Ng (2000). More...

Functions/Subroutines

subroutine w3sbt9 (AC, H_WDEPTH, S, D, IX, IY)
 Compute dissipation by viscous fluid mud using Ng (2000) (adapted from Erick Rogers code by Mark Orzech, NRL). More...
 
subroutine ng (SIGMA, H_WDEPTH, DTILDE, ZETA, SBLTM, GAMMA, WK, WKDR, DISS)
 Compute dissipation by viscous fluid mud using Ng (2000). More...
 
subroutine calc_nd (KWAVE, H_WDEPTH, SND2, ND)
 NA. More...
 

Detailed Description

Contains routines for computing dissipation by viscous fluid mud using Ng (2000).

Author
M. Orzech
W. E. Rogers
Date
21-Nov-2013
Copyright
Copyright 2009-2022 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.

Function/Subroutine Documentation

◆ calc_nd()

subroutine w3sbt9md::calc_nd ( real, intent(in)  KWAVE,
real, intent(in)  H_WDEPTH,
real, intent(in)  SND2,
real, intent(out)  ND 
)

NA.

Parameters
[in]KWAVE
[in]H_WDEPTH
[in]SND2
[out]ND
Author
NA
Date
NA

Definition at line 572 of file w3sbt9md.F90.

572  !/ ------------------------------------------------------------------- /
573 
574  IMPLICIT NONE
575  REAL, INTENT(IN) :: KWAVE
576  REAL, INTENT(IN) :: H_WDEPTH
577  REAL, INTENT(IN) :: SND2
578  REAL, INTENT(OUT) :: ND
579  REAL :: FAC1 ! LOCAL
580  REAL :: FAC2 ! LOCAL
581  REAL :: FAC3 ! LOCAL
582  REAL :: KND ! LOCAL
583 
584  knd = kwave*h_wdepth
585  fac1 = 2.*knd/sinh(2.*knd)
586  fac2 = snd2/knd
587  fac3 = 2.*fac2/(1.+fac2*fac2)
588  nd= fac1*(0.5/h_wdepth - kwave/fac3)
589 

Referenced by w3sbt9().

◆ ng()

subroutine w3sbt9md::ng ( real, intent(in)  SIGMA,
real, intent(in)  H_WDEPTH,
real, intent(in)  DTILDE,
real, intent(in)  ZETA,
real, intent(in)  SBLTM,
real, intent(in)  GAMMA,
real, intent(in)  WK,
real, intent(out)  WKDR,
real, intent(out)  DISS 
)

Compute dissipation by viscous fluid mud using Ng (2000).

Adapted from Erick Rogers code by Mark Orzech, NRL.

Parameters
[in]SIGMARadian frequency (rad).
[in]H_WDEPTHWater depth, denoted "h" in Ng (m).
[in]DTILDENormalized mud depth.
[in]ZETAThe ratio of stokes' boundary layer.
[in]SBLTMSbltm is what you get if you calculate sblt using the viscosity of the mud
[in]GAMMAThe gamma used in Ng page 238, density(water)/density(mud).
[in]WKUnmuddy wavenumber.
[out]WKDRMuddy wavenumber.
[out]DISSDissipation rate.
Author
E. Rogers
M. Orzech
Date
21-Nov-2013

Definition at line 430 of file w3sbt9md.F90.

430  !/
431  !/ +-----------------------------------+
432  !/ | WAVEWATCH III NOAA/NCEP |
433  !/ | E. Rogers and M. Orzech |
434  !/ | FORTRAN 90 |
435  !/ | Last update : 21-Nov-2013 |
436  !/ +-----------------------------------+
437  !
438  !/ 28-Jul-2011 : Origination. ( version 4.01 )
439  !/ 21-Nov-2013 : Preparing distribution version. ( version 4.11 )
440  !/
441  ! 1. Purpose :
442  !
443  ! Compute dissipation by viscous fluid mud using Ng (2000)
444  ! (adapted from Erick Rogers code by Mark Orzech, NRL).
445  !
446  ! 2. Method :
447  !
448  ! 3. Parameters :
449  !
450  ! Parameter list
451  ! ----------------------------------------------------------------
452  ! SIGMA Real I radian frequency (rad)
453  ! H_WDEPTH Real I water depth
454  ! DTILDE Real I normalized mud depth
455  ! ZETA Real I zeta as used in Ng
456  ! SBLTM Real I mud Stokes boundary layer thickness
457  ! GAMMA Real I gamma as used in Ng
458  ! WK Real I wavenumber w/out mud
459  ! WKDR Real O wavenumber w/mud
460  ! DISS Real O dissipation rate
461  ! ----------------------------------------------------------------
462  !
463  ! 4. Subroutines used :
464  !
465  ! None.
466  !
467  ! 5. Called by :
468  !
469  ! Name Type Module Description
470  ! ----------------------------------------------------------------
471  ! W3SBT9 Subr. W3SBT9MD Main routine (all freqs)
472  ! ----------------------------------------------------------------
473  !
474  ! 6. Error messages :
475  !
476  ! None.
477  !
478  ! 7. Remarks :
479  ! Calculations for the "B coefficients" came from a code by Jim Kaihatu
480  !
481  ! 8. Structure :
482  !
483  ! See source code.
484  !
485  ! 9. Switches :
486  !
487  ! None.
488  !
489  ! 10. Source code :
490  !
491  !/ ------------------------------------------------------------------- /
492  !/ ------------------------------------------------------------------- /
493  !/
494 
495  !
496  IMPLICIT NONE
497 
498  ! INPUT VARIABLES :
499  REAL, INTENT(IN) :: SIGMA ! radian frequency (rad)
500  REAL, INTENT(IN) :: H_WDEPTH! water depth, denoted "h" in Ng (m)
501  REAL, INTENT(IN) :: DTILDE ! normalized mud depth = mud depth / sbltm,
502  ! delta is the sblt= sqrt(2*visc/sigma)
503  REAL, INTENT(IN) :: ZETA ! this is zeta as used in Ng pg. 238. it is
504  ! the ratio of stokes' boundary layer
505  ! thicknesses, or sbltm/delta_w
506  REAL, INTENT(IN) :: GAMMA ! this is the gamma used in Ng pg. 238.
507  ! this is density(water)/density(mud)
508  REAL, INTENT(IN) :: SBLTM ! sbltm is what you get if you calculate
509  ! sblt using the viscosity of the mud,
510  ! sbltm=sqrt(2*visc_m/sigma)
511  ! .....also delta_m
512  REAL, INTENT(IN) :: WK ! unmuddy wavenumber
513 
514  ! OUTPUT VARIABLES :
515  REAL, INTENT(OUT) :: WKDR ! muddy wavenumber
516  REAL, INTENT(OUT) :: DISS ! dissipation rate
517 
518  ! LOCAL VARIABLES :
519  REAL :: B1 ! an Ng coefficient
520  REAL :: B2 ! an Ng coefficient
521  REAL :: B3 ! an Ng coefficient
522  REAL :: BR ! an Ng coefficient
523  REAL :: BI ! an Ng coefficient
524  REAL :: BRP ! an Ng coefficient
525  REAL :: BIP ! an Ng coefficient
526  REAL :: DM ! MUD DEPTH, ADDED JUNE 2 2006
527 
528 
529  dm=dtilde*sbltm ! DTILDE=DM/SBLTM
530  ! NOW CALCULATE Ng's B coefficients : see Ng pg 238
531  b1=gamma*(-2.0*gamma**2+2.0*gamma-1.-zeta**2)*sinh(dtilde)* &
532  cosh(dtilde)-gamma**2*zeta*((cosh(dtilde))**2+ &
533  (sinh(dtilde))**2)-(gamma-1.)**2*zeta*((cosh(dtilde))**2 &
534  *(cos(dtilde))**2+(sinh(dtilde))**2*(sin(dtilde))**2)-2.0 &
535  *gamma*(1.-gamma)*(zeta*cosh(dtilde)+gamma*sinh(dtilde)) &
536  *cos(dtilde)
537 
538  b2=gamma*(-2.0*gamma**2+2.0*gamma-1.+zeta**2)*sin(dtilde)* &
539  cos(dtilde) -2.0*gamma*(1.-gamma)*(zeta*sinh(dtilde)+gamma &
540  *cosh(dtilde))*sin(dtilde)
541 
542  b3=(zeta*cosh(dtilde)+gamma*sinh(dtilde))**2*(cos(dtilde))**2 &
543  +(zeta*sinh(dtilde)+gamma*cosh(dtilde))**2*(sin(dtilde))**2
544 
545  br=wk*sbltm*(b1-b2)/(2.0*b3)+gamma*wk*dm
546 
547  bi=wk*sbltm*(b1+b2)/(2.0*b3)
548  brp=b1/b3 ! "B_R PRIME"
549  bip=b2/b3 ! "B_I PRIME"
550 
551  ! now calculate dissipation rate and wavenumber
552  diss=-sbltm*(brp+bip)*wk**2/(sinh(2.0*wk*h_wdepth)+2.0*wk*h_wdepth)
553  wkdr=wk-br*wk/(sinh(wk*h_wdepth)*cosh(wk*h_wdepth)+wk*h_wdepth)
554 
555  RETURN
556 

Referenced by w3sbt9().

◆ w3sbt9()

subroutine w3sbt9md::w3sbt9 ( real, dimension(nspec), intent(in)  AC,
real, intent(in)  H_WDEPTH,
real, dimension(nspec), intent(out)  S,
real, dimension(nspec), intent(out)  D,
integer, intent(in)  IX,
integer, intent(in)  IY 
)

Compute dissipation by viscous fluid mud using Ng (2000) (adapted from Erick Rogers code by Mark Orzech, NRL).

Parameters
[in]ACAction density.
[in]H_WDEPTHMean water depth.
[out]SSource term (1-D version).
[out]DDiagonal term of derivative (1-D version).
[in]IX
[in]IY
Author
M. Orzech
W. E. Rogers
Date
21-Nov-2013

Definition at line 119 of file w3sbt9md.F90.

119  !/
120  !/ +-----------------------------------+
121  !/ | WAVEWATCH III NOAA |
122  !/ | M. Orzech NRL |
123  !/ | W. E. Rogers NRL |
124  !/ | FORTRAN 90 |
125  !/ | Last update : 21-Nov-2013 |
126  !/ +-----------------------------------+
127 
128  !/ 28-Jul-2011 : Origination. ( version 4.01 )
129  !/ 21-Nov-2013 : Preparing distribution version. ( version 4.11 )
130  !/
131  ! 1. Purpose :
132  !
133  ! Compute dissipation by viscous fluid mud using Ng (2000)
134  ! (adapted from Erick Rogers code by Mark Orzech, NRL).
135  !
136  ! 2. Method :
137  !
138  ! 3. Parameters :
139  !
140  ! Parameter list
141  ! ----------------------------------------------------------------
142  ! H_WDEPTH Real I Mean water depth.
143  ! S R.A. O Source term (1-D version).
144  ! D R.A. O Diagonal term of derivative (1-D version).
145  ! ----------------------------------------------------------------
146  !
147  ! 4. Subroutines used :
148  !
149  ! Name Type Module Description
150  ! ----------------------------------------------------------------
151  ! STRACE Subr. W3SERVMD Subroutine tracing.
152  ! CALC_ND
153  ! NG
154  ! ----------------------------------------------------------------
155  !
156  ! 5. Called by :
157  !
158  ! Name Type Module Description
159  ! ----------------------------------------------------------------
160  ! W3SRCE Subr. W3SRCEMD Source term integration.
161  ! W3EXPO Subr. N/A Point output post-processor.
162  ! GXEXPO Subr. N/A GrADS point output post-processor.
163  ! ----------------------------------------------------------------
164  !
165  ! 6. Error messages :
166  !
167  ! None.
168  !
169  ! 7. Remarks :
170  !
171  ! Cg_mud calculation could be improved by using dsigma/dk instead
172  ! of n*C. The latter is a "naive" method and its accuracy has
173  ! not been confirmed.
174  !
175  ! 8. Structure :
176  !
177  ! See source code.
178  !
179  ! 9. Switches :
180  !
181  ! !/S Enable subroutine tracing.
182  !
183  ! 10. Source code :
184  !
185  !/ ------------------------------------------------------------------- /
186  USE w3gdatmd, ONLY: nk,sig,nspec,mapwn
187  USE w3idatmd, ONLY: mudt, mudv, mudd, inflags1
188  USE constants, ONLY: pi,grav,dwat,nu_water
189  USE w3servmd, ONLY: extcde
190  USE w3odatmd, ONLY: ndse
191 #ifdef W3_S
192  USE w3servmd, ONLY: strace
193 #endif
194  !/
195  IMPLICIT NONE
196  !/
197  !/ ------------------------------------------------------------------- /
198  !/ Parameter list
199  !/
200  REAL, INTENT(IN) :: H_WDEPTH ! WATER DEPTH, DENOTED "H" IN NG (M)
201  REAL, INTENT(IN) :: AC(NSPEC) ! ACTION DENSITY
202  INTEGER, INTENT(IN) :: IX, IY
203  REAL, INTENT(OUT) :: S(NSPEC), D(NSPEC)
204  !/
205  !/ ------------------------------------------------------------------- /
206  !/ Local parameters
207  !/
208 #ifdef W3_S
209  INTEGER, SAVE :: IENT = 0
210 #endif
211 
212  ! LOCAL VARIABLES
213  REAL :: DMW(NK)
214  REAL :: ROOTDG
215  REAL :: SND
216  REAL :: SND2
217  REAL :: WGD
218  REAL :: CWAVE
219  REAL :: KD_ROCK
220  REAL :: CG_MUD
221  REAL :: K_MUD
222  REAL :: NWAVE_MUD
223  REAL :: ND_MUD
224  REAL :: SMUDWD(NK) ! DISSIPATION DUE TO MUD
225  REAL :: CG_ROCK
226  REAL :: K_ROCK
227  REAL :: NWAVE_ROCK
228  REAL :: ND_ROCK
229  REAL :: KINVISM ! := THE KINEMATIC VISCOSITY OF THE MUD
230  REAL :: KINVISW ! := KINEMATIC VISCOSITY OF WATER
231  REAL :: RHOW ! := DENSITY OF WATER
232  REAL :: RHOM ! := DENSITY OF MUD
233  REAL :: DM ! := DEPTH OF MUD LAYER
234  REAL :: ZETA ! := THIS IS ZETA AS USED IN NG PG. 238. IT IS THE
235  ! RATIO OF STOKES' BOUNDARY LAYER THICKNESSES,
236  ! OR DELTA_M/DELTA_W
237  REAL :: GAMMA ! := THIS IS THE GAMMA USED IN NG PG. 238. THIS IS
238  ! DENSITY(WATER)/DENSITY(MUD)
239  REAL :: SBLTW ! := A FUNCTION OF VISCOSITY AND FREQ
240  REAL :: SBLTM ! := A FUNCTION OF VISCOSITY AND FREQ
241  REAL :: DTILDE ! := NORMALIZED MUD DEPTH = MUD DEPTH / DELTA_M,
242  ! DELTA IS THE SBLT= SQRT(2*VISC/SIGMA)
243  REAL :: ZTMP
244  REAL :: KDCUTOFF
245  REAL :: KD
246 
247  INTEGER :: IS
248 
249  LOGICAL :: INAN
250 
251  !/ ------------------------------------------------------------------- /
252  !/
253 #ifdef W3_S
254  CALL strace (ient, 'W3SBT9')
255 #endif
256  !
257  ! 0. Initializations ------------------------------------------------ *
258  !
259  ! Ng (2000), Waves over soft muds.
260  ! Based on code for SWAN created by Erick Rogers.
261  ! Adapted for WW3 by Mark Orzech, Nov 2012.
262 
263  ! Initialize properties from mud fields if available
264  IF (inflags1(-2))THEN
265  rhom = mudd(ix,iy)
266  ELSE
267  WRITE(ndse,*)'RHOM NOT SET'
268  CALL extcde ( 1 )
269  ENDIF
270  IF (inflags1(-1)) THEN
271  dm = mudt(ix,iy)
272  ELSE
273  WRITE(ndse,*)'DM NOT SET'
274  CALL extcde ( 2 )
275  ENDIF
276  IF (inflags1(0)) THEN
277  kinvism = mudv(ix,iy)
278  ELSE
279  WRITE(ndse,*)'KINVISM NOT SET'
280  CALL extcde ( 3 )
281  ENDIF
282 
283  rootdg = sqrt(h_wdepth/grav)
284  wgd = rootdg*grav
285  DO is = 1, nk
286  ! SND is dimensionless frequency
287  snd = sig(is) * rootdg
288  IF (snd .GE. 2.5) THEN
289  ! ******* DEEP WATER *******
290  k_rock = sig(is) * sig(is) / grav
291  cg_rock = 0.5 * grav / sig(is)
292  nwave_rock = 0.5
293  nd_rock = 0.
294  ELSE IF (snd.LT.1.e-6) THEN
295  ! *** VERY SHALLOW WATER ***
296  k_rock = snd/h_wdepth
297  cg_rock = wgd
298  nwave_rock = 1.
299  nd_rock = 0.
300  ELSE
301 
302  snd2 = snd*snd
303  cwave = sqrt(grav*h_wdepth/(snd2+1./(1.+0.666*snd2 &
304  +0.445*snd2**2 -0.105*snd2**3+0.272*snd2**4)))
305  k_rock = sig(is)/cwave
306 
307  CALL calc_nd(k_rock,h_wdepth,snd2,nd_rock)
308 
309  nwave_rock = 0.5*(1.0+2.0*k_rock*h_wdepth/sinh(2.0*k_rock*h_wdepth))
310  cg_rock= nwave_rock*cwave
311 
312  snd2=0
313  cwave=0
314 
315  ENDIF
316 
317  kdcutoff = 10.0 ! hardwired (same as w3sbt8md)
318 
319  ! now that kh is known, we can use a definition of "deep" that is
320  ! consistent with the definition used in sbot
321  k_mud=0.0
322  dmw(is)=0.0
323  kd_rock = k_rock * h_wdepth
324  ! KD_ROCK is used to determine whether we make the mud calculation
325  IF((kd_rock.LT.kdcutoff).AND.(dm.GT.1.0e-5))THEN
326  kinvisw=nu_water
327  rhow=dwat
328  zeta=sqrt(kinvism/kinvisw)
329  gamma=rhow/rhom
330  sbltw=sqrt(2.0*kinvisw/sig(is))
331  sbltm=sqrt(2.0*kinvism/sig(is))
332 
333  dtilde=dm/sbltm
334  CALL ng(sig(is),h_wdepth,dtilde,zeta,sbltm,gamma,k_rock,k_mud, &
335  dmw(is))
336 
337  ELSE ! IF ( KD_ROCK .LT. KDCUTOFF ) THEN
338  k_mud=k_rock
339  END IF ! IF ( KD_ROCK .LT. KDCUTOFF ) THEN
340 
341  ! calculate cg_mud, nwave_mud here
342  cwave=sig(is)/k_mud
343 
344  ztmp=2.0*k_mud*h_wdepth
345  IF(ztmp.LT.70)THEN
346  ztmp=sinh(ztmp)
347  ELSE
348  ztmp=1.0e+30
349  ENDIF
350  nwave_mud=0.5*(1.0+2.0*k_mud*h_wdepth/ztmp)
351 
352  cg_mud=nwave_mud*cwave
353  snd2 = snd*snd
354 
355  CALL calc_nd(k_mud,h_wdepth,snd2,nd_mud)
356 
357  snd2=0
358  cwave=0
359 
360  ! If we wanted to include the effects of mud on the real part of the
361  ! wavnumber (as we do in SWAN), this is where we would do it.
362  ! Set output variables k_out, cg_out, nwave_out, nd_out, dmw.
363  !kinematics IF(MUD)THEN !
364  !kinematics K_OUT(IS) =K_MUD
365  !kinematics CG_OUT(IS) =CG_MUD
366  !kinematics NWAVE_OUT(IS)=NWAVE_MUD
367  !kinematics ND_OUT(IS) =ND_MUD
368  !kinematics ELSE ! USE ROCKY WAVENUMBER,ETC.
369  !kinematics K_OUT(IS) =K_ROCK
370  !kinematics CG_OUT(IS) =CG_ROCK
371  !kinematics NWAVE_OUT(IS)=NWAVE_ROCK
372  !kinematics ND_OUT(IS) =ND_ROCK
373  !kinematics DMW(IS)=0.0
374  !kinematics ENDIF
375 
376  kd = k_mud * h_wdepth
377  IF ( kd .LT. kdcutoff ) THEN
378  ! note that "IS" here is for the 1d spectrum
379  smudwd(is)=2.0*dmw(is)*cg_mud
380  END IF
381 
382  ! NaN check:
383  inan = .NOT. ( dmw(is) .GE. -huge(dmw(is)) .AND. dmw(is) &
384  .LE. huge(dmw(is)) )
385  IF (inan) THEN
386  WRITE(*,'(/1A/)') 'W3SBT9 ERROR -- DMW(IS) IS NAN'
387  WRITE(*,*)'W3SBT9: RHOM, DM, KINVISM = ',rhom, dm, kinvism
388  WRITE(*,*)'W3SBT9: IS,NK = ',is,nk
389  WRITE(*,*)'W3SBT9: H_WDEPTH,KD,KDCUTOFF = ',h_wdepth,kd, kdcutoff
390  WRITE(*,*)'W3SBT9: K_MUD,CG_MUD,NWAVE_MUD = ',k_mud,cg_mud,nwave_mud
391  CALL extcde (1)
392  END IF
393 
394  END DO ! DO IS = 1, NK
395 
396  ! *** store the results in the DIAGONAL arrays D and S ***
397  DO is = 1,nspec
398  ! note that "IS" here is for the directional spectrum (2d)
399  d(is) = -smudwd(mapwn(is))
400  END DO
401 
402  s = d * ac
403 
404  RETURN
405 

References calc_nd(), constants::dwat, w3servmd::extcde(), constants::grav, w3idatmd::inflags1, w3gdatmd::mapwn, w3odatmd::ndse, ng(), w3gdatmd::nk, w3gdatmd::nspec, constants::nu_water, constants::pi, w3gdatmd::sig, and w3servmd::strace().

Referenced by w3srcemd::w3srce().

constants::pi
real, parameter pi
PI Value of Pi.
Definition: constants.F90:71
w3idatmd::inflags1
logical, dimension(:), pointer inflags1
Definition: w3idatmd.F90:260
w3gdatmd::sig
real, dimension(:), pointer sig
Definition: w3gdatmd.F90:1234
w3sbt9md::calc_nd
subroutine calc_nd(KWAVE, H_WDEPTH, SND2, ND)
NA.
Definition: w3sbt9md.F90:572
w3odatmd::ndse
integer, pointer ndse
Definition: w3odatmd.F90:456
w3sbt9md::ng
subroutine ng(SIGMA, H_WDEPTH, DTILDE, ZETA, SBLTM, GAMMA, WK, WKDR, DISS)
Compute dissipation by viscous fluid mud using Ng (2000).
Definition: w3sbt9md.F90:430
w3servmd
Definition: w3servmd.F90:3
constants::nu_water
real, parameter nu_water
NU_WATER Kinematic viscosity of water (m2/s).
Definition: constants.F90:67
w3odatmd
Definition: w3odatmd.F90:3
constants::dwat
real, parameter dwat
DWAT Density of water (kg/m3).
Definition: constants.F90:62
w3servmd::strace
subroutine strace(IENT, SNAME)
Definition: w3servmd.F90:148
w3idatmd
Define data structures to set up wave model input data for several models simultaneously.
Definition: w3idatmd.F90:16
w3gdatmd::mapwn
integer, dimension(:), pointer mapwn
Definition: w3gdatmd.F90:1231
constants
Define some much-used constants for global use (all defined as PARAMETER).
Definition: constants.F90:20
w3gdatmd
Definition: w3gdatmd.F90:16
constants::grav
real, parameter grav
GRAV Acc.
Definition: constants.F90:61