ROMS
Loading...
Searching...
No Matches
mb_bbl.h
Go to the documentation of this file.
1 MODULE bbl_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Meinte Blaas !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This subroutine computes bottom stresses for combined waves and !
11! currents using the parametric approximation by Soulsby 1997: !
12! !
13! tauCW=tauC*[1+1.2*(tauW/(tauC+tauW))^3.2] !
14! !
15! and !
16! !
17! tauCWmax=SQRT([tauCW+tauW*COS(phiCW)]^2+[tauW*SIN(phiCW)]^2) !
18! !
19! where !
20! !
21! tauCW Combined wave-averaged stress (in current direction). !
22! tauC Stress due to currents, if waves would be absent. !
23! tauW Amplitude of stress due to waves without currents. !
24! tauCWmax Maximum combined wave-averaged stress. !
25! phiCW Angle between current and waves. !
26! !
27! References: !
28! !
29! Dyer 1986, Coastal and Estuarine Sediment Dynamics, Wiley, !
30! 342 pp. !
31! !
32! Harris and Wiberg 2001, Comp. and Geosci. 27, 675-690. !
33! !
34! Li and Amos 2001, Comp. and Geosci. 27, 619-645. !
35! !
36! Soulsby 1997, Dynamics of Marine Sands, Telford Publ., 249 pp. !
37! !
38! Soulsby 1995, Bed shear-stresses due to combined waves and !
39! currents, in: Stive et al: Advances in Coastal Morphodynamics, !
40! Wiley, 4.20-4.23 !
41! !
42! Wiberg and Harris 1994, J. Geophys. Res. 99(C4), 775-789. !
43! !
44!=======================================================================
45!
46 implicit none
47
48 PRIVATE
49 PUBLIC :: bblm
50
51 CONTAINS
52!
53!***********************************************************************
54 SUBROUTINE bblm (ng,tile)
55!***********************************************************************
56!
57 USE mod_param
58 USE mod_bbl
59 USE mod_forces
60 USE mod_grid
61 USE mod_ocean
62 USE mod_sedbed
63 USE mod_stepping
64!
65! Imported variable declarations.
66!
67 integer, intent(in) :: ng, tile
68!
69! Local variable declarations.
70!
71 character (len=*), parameter :: myfile = &
72 & __FILE__
73!
74# include "tile.h"
75!
76# ifdef PROFILE
77 CALL wclock_on (ng, inlm, 37, __line__, myfile)
78# endif
79 CALL mb_bbl_tile (ng, tile, &
80 & lbi, ubi, lbj, ubj, &
81 & imins, imaxs, jmins, jmaxs, &
82 & nrhs(ng), &
83 & grid(ng) % h, &
84 & grid(ng) % z_r, &
85 & grid(ng) % z_w, &
86 & grid(ng) % angler, &
87 & grid(ng) % ZoBot, &
88# ifdef MB_CALC_UB
89 & forces(ng) % Hwave, &
90# else
91 & forces(ng) % Uwave_rms, &
92# endif
93 & forces(ng) % Dwave, &
94 & forces(ng) % Pwave_bot, &
95 & ocean(ng) % rho, &
96 & ocean(ng) % u, &
97 & ocean(ng) % v, &
98 & sedbed(ng) % bottom, &
99 & bbl(ng) % Ubot, &
100 & bbl(ng) % Vbot, &
101 & bbl(ng) % Ur, &
102 & bbl(ng) % Vr, &
103 & bbl(ng) % bustrc, &
104 & bbl(ng) % bvstrc, &
105 & bbl(ng) % bustrw, &
106 & bbl(ng) % bvstrw, &
107 & bbl(ng) % bustrcwmax, &
108 & bbl(ng) % bvstrcwmax, &
109 & forces(ng) % bustr, &
110 & forces(ng) % bvstr)
111# ifdef PROFILE
112 CALL wclock_off (ng, inlm, 37, __line__, myfile)
113# endif
114!
115 RETURN
116 END SUBROUTINE bblm
117!
118!***********************************************************************
119 SUBROUTINE mb_bbl_tile (ng, tile, &
120 & LBi, UBi, LBj, UBj, &
121 & IminS, ImaxS, JminS, JmaxS, &
122 & nrhs, &
123 & h, z_r, z_w, angler, ZoBot, &
124# ifdef MB_CALC_UB
125 & Hwave, &
126# else
127 & Uwave_rms, &
128# endif
129 & Dwave, Pwave_bot, &
130 & rho, u, v, bottom, &
131 & Ubot, Vbot, Ur, Vr, &
132 & bustrc, bvstrc, &
133 & bustrw, bvstrw, &
134 & bustrcwmax, bvstrcwmax, &
135 & bustr, bvstr)
136!***********************************************************************
137!
138 USE mod_param
139 USE mod_scalars
140 USE mod_sediment
141!
142 USE bc_2d_mod
143# ifdef DISTRIBUTE
144 USE mp_exchange_mod, ONLY : mp_exchange2d
145# endif
146!
147! Imported variable declarations.
148!
149 integer, intent(in) :: ng, tile
150 integer, intent(in) :: LBi, UBi, LBj, UBj
151 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
152 integer, intent(in) :: nrhs
153!
154# ifdef ASSUMED_SHAPE
155 real(r8), intent(in) :: h(LBi:,LBj:)
156 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
157 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
158 real(r8), intent(in) :: angler(LBi:,LBj:)
159 real(r8), intent(in) :: ZoBot(LBi:,LBj:)
160# ifdef MB_CALC_UB
161 real(r8), intent(in) :: Hwave(LBi:,LBj:)
162# else
163 real(r8), intent(in) :: Uwave_rms(LBi:,LBj:)
164# endif
165 real(r8), intent(in) :: Dwave(LBi:,LBj:)
166 real(r8), intent(in) :: Pwave_bot(LBi:,LBj:)
167 real(r8), intent(in) :: rho(LBi:,LBj:,:)
168 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
169 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
170
171 real(r8), intent(inout) :: bottom(LBi:,LBj:,:)
172
173 real(r8), intent(out) :: Ubot(LBi:,LBj:)
174 real(r8), intent(out) :: Vbot(LBi:,LBj:)
175 real(r8), intent(out) :: Ur(LBi:,LBj:)
176 real(r8), intent(out) :: Vr(LBi:,LBj:)
177 real(r8), intent(out) :: bustrc(LBi:,LBj:)
178 real(r8), intent(out) :: bvstrc(LBi:,LBj:)
179 real(r8), intent(out) :: bustrw(LBi:,LBj:)
180 real(r8), intent(out) :: bvstrw(LBi:,LBj:)
181 real(r8), intent(out) :: bustrcwmax(LBi:,LBj:)
182 real(r8), intent(out) :: bvstrcwmax(LBi:,LBj:)
183 real(r8), intent(out) :: bustr(LBi:,LBj:)
184 real(r8), intent(out) :: bvstr(LBi:,LBj:)
185# else
186 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
187 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
188 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
189 real(r8), intent(in) :: angler(LBi:UBi,LBj:UBj)
190 real(r8), intent(in) :: ZoBot(LBi:UBi,LBj:UBj)
191# ifdef MB_CALC_UB
192 real(r8), intent(in) :: Hwave(LBi:UBi,LBj:UBj)
193# else
194 real(r8), intent(in) :: Uwave_rms(LBi:UBi,LBj:UBj)
195# endif
196 real(r8), intent(in) :: Dwave(LBi:UBi,LBj:UBj)
197 real(r8), intent(in) :: Pwave_bot(LBi:UBi,LBj:UBj)
198 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
199 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
200 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
201
202 real(r8), intent(inout) :: bottom(LBi:UBi,LBj:UBj,MBOTP)
203
204 real(r8), intent(out) :: Ubot(LBi:UBi,LBj:UBj)
205 real(r8), intent(out) :: Vbot(LBi:UBi,LBj:UBj)
206 real(r8), intent(out) :: Ur(LBi:UBi,LBj:UBj)
207 real(r8), intent(out) :: Vr(LBi:UBi,LBj:UBj)
208 real(r8), intent(out) :: bustrc(LBi:UBi,LBj:UBj)
209 real(r8), intent(out) :: bvstrc(LBi:UBi,LBj:UBj)
210 real(r8), intent(out) :: bustrw(LBi:UBi,LBj:UBj)
211 real(r8), intent(out) :: bvstrw(LBi:UBi,LBj:UBj)
212 real(r8), intent(out) :: bustrcwmax(LBi:UBi,LBj:UBj)
213 real(r8), intent(out) :: bvstrcwmax(LBi:UBi,LBj:UBj)
214 real(r8), intent(out) :: bustr(LBi:UBi,LBj:UBj)
215 real(r8), intent(out) :: bvstr(LBi:UBi,LBj:UBj)
216# endif
217!
218! Local variable declarations.
219!
220 integer :: i, ised, j
221
222 real(r8), parameter :: K1 = 0.6666666666_r8
223 real(r8), parameter :: K2 = 0.3555555555_r8
224 real(r8), parameter :: K3 = 0.1608465608_r8
225 real(r8), parameter :: K4 = 0.0632098765_r8
226 real(r8), parameter :: K5 = 0.0217540484_r8
227 real(r8), parameter :: K6 = 0.0065407983_r8
228
229 real(r8), parameter :: eps = 1.0e-10_r8
230
231 real(r8), parameter :: scf1 = 0.5_r8 * 1.39_r8
232 real(r8), parameter :: scf2 = 0.52_r8
233 real(r8), parameter :: scf3 = 2.0_r8 - scf2
234 real(r8), parameter :: scf4 = 1.2_r8
235 real(r8), parameter :: scf5 = 3.2_r8
236
237 real(r8) :: twopi = 2.0_r8 * pi
238
239 real(r8) :: RHbiomax = 0.006_r8 ! maximum biogenic ripple height
240 real(r8) :: RHmin = 0.001_r8 ! arbitrary minimum ripple height
241 real(r8) :: RLmin = 0.01_r8 ! arbitrary minimum ripple length
242
243 real(r8) :: Ab, Fwave_bot, Kbh, Kbh2, Kdh
244 real(r8) :: angleC, angleW, phiC, phiCW
245 real(r8) :: cff, cff1, cff2, d50, viscosity, wset
246 real(r8) :: RHbio, RHbiofac, RHmax, RLbio
247 real(r8) :: rhoWater, rhoSed
248 real(r8) :: rhgt, rlen, thetw
249 real(r8) :: Znot, ZnotC, Znot_bl, Znot_rip
250 real(r8) :: tau_bf, tau_ex, tau_en, tau_up, tau_w, tau_wb
251 real(r8) :: tau_c, tau_cb, tau_cs, tau_cw, tau_cwb, tau_cws
252
253 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Ub
254 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Ucur
255 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Vcur
256 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Zr
257 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Ur_mb
258 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Vr_mb
259 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Umag
260 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tauC
261 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tauW
262 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tauCW
263 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tauCWmax
264
265#include "set_bounds.h"
266!
267!-----------------------------------------------------------------------
268! Initalize stresses due to currents and waves.
269!-----------------------------------------------------------------------
270!
271 rhbiofac=1.0_r8/exp(4.11_r8)
272 DO j=jstrv-1,jend
273 DO i=istru-1,iend
274 tauc(i,j)=0.0_r8
275 taucw(i,j)=0.0_r8
276 taucwmax(i,j)=0.0_r8
277 END DO
278 END DO
279!
280!-----------------------------------------------------------------------
281! Set currents above bed.
282!-----------------------------------------------------------------------
283!
284 DO j=jstrv-1,jend+1
285 DO i=istru-1,iend+1
286 zr(i,j)=z_r(i,j,1)-z_w(i,j,0)
287 ur_mb(i,j)=u(i,j,1,nrhs)
288 vr_mb(i,j)=v(i,j,1,nrhs)
289 END DO
290 END DO
291!
292!=======================================================================
293! Compute bottom stresses.
294!=======================================================================
295!
296 DO j=jstrv-1,jend
297 DO i=istru-1,iend
298 rhbio=0.0_r8
299 rlbio=0.0_r8
300 rlen=bottom(i,j,irlen)
301 rhgt=bottom(i,j,irhgt)
302 rhowater=rho(i,j,1)+1000.0_r8
303 viscosity=0.0013_r8/rhowater ! kinematic viscosity
304!
305!-----------------------------------------------------------------------
306! Compute bed wave orbital velocity (m/s) and excursion amplitude (m)
307! from wind-induced waves. Use Dean and Dalrymple (1991) 6th-degree
308! polynomial to approximate wave number on shoaling water.
309!-----------------------------------------------------------------------
310!
311 fwave_bot=twopi/max(pwave_bot(i,j),0.05_r8)
312# ifdef MB_BBL_CALC_UB
313 kdh=h(i,j)*fwave_bot*fwave_bot/g
314 kbh2=kdh*kdh+ &
315 & kdh/(1.0_r8+kdh*(k1+kdh*(k2+kdh*(k3+kdh*(k4+ &
316 & kdh*(k5+k6*kdh))))))
317 kbh=sqrt(kbh2)
318!
319! Compute bed wave orbital velocity and excursion amplitude.
320!
321 ab=0.5_r8*hwave(i,j)/sinh(kbh)+eps
322 ub(i,j)=fwave_bot*ab
323# else
324 ub(i,j)=uwave_rms(i,j)
325 ab=ub(i,j)/fwave_bot+eps
326# endif
327!
328! Compute bottom current magnitude at RHO-points.
329!
330 ucur(i,j)=0.5_r8*(ur_mb(i,j)+ur_mb(i+1,j))
331 vcur(i,j)=0.5_r8*(vr_mb(i,j)+vr_mb(i,j+1))
332 umag(i,j)=sqrt(ucur(i,j)*ucur(i,j)+vcur(i,j)*vcur(i,j))+eps
333!
334! Compute angle between currents and waves (radians)
335!
336 IF (ucur(i,j).ne.0.0_r8) THEN
337 phic=atan2(vcur(i,j),ucur(i,j))
338 ELSE
339 phic=0.5_r8*pi*sign(1.0_r8,vcur(i,j))
340 END IF
341 phicw=1.5_r8*pi-dwave(i,j)-phic-angler(i,j)
342!
343!-----------------------------------------------------------------------
344! Determine skin roughness from sediment size. Set default logarithmic
345! profile consistent with current-only case.
346! Establish local median grain size for all calculations here and
347! determine local values of critical stresses.
348! Since most parameterizations have been derived ignoring multiple
349! grain sizes, we apply this single d50 also in the case of mixed beds.
350!-----------------------------------------------------------------------
351!
352 d50=bottom(i,j,isd50) ! (m)
353 tau_cb=bottom(i,j,itauc) ! (m2/s2)
354 wset=bottom(i,j,iwsed) ! (m/s)
355 rhosed=bottom(i,j,idens)/rhowater ! (nondimensional)
356!
357! Critical stress (m2/s2) for transition to sheet flow
358! (Li and Amos, 2001, eq. 11).
359!
360 tau_up=0.172_r8*(rhosed-1.0_r8)*g*(d50**0.624_r8)
361!
362! Threshold stress (m2/s2) for break off (Grant and Madsen,1982).
363!
364 tau_bf=0.79_r8*(viscosity**(-0.6_r8))* &
365 & (((rhosed-1.0_r8)*g)**0.3_r8)*(d50**0.9_r8)*tau_cb
366!
367! Set Znot for currents as maximum of user value or grain roughness.
368!
369 znotc=d50/12.0_r8
370 znot=max(zobot(i,j),znotc)
371!
372!-----------------------------------------------------------------------
373! Set tauC (m2/s2) acc. to current-only case (skin friction) [m/s]
374!-----------------------------------------------------------------------
375!
376 cff1=vonkar/log(zr(i,j)/znot)
377 cff2=min(cdb_max,max(cdb_min,cff1*cff1))
378 tauc(i,j)=cff2*umag(i,j)*umag(i,j)
379
380 cff1=vonkar/log(zr(i,j)/znotc)
381 tau_cs=cff1*cff1*umag(i,j)*umag(i,j)
382!
383!-----------------------------------------------------------------------
384! If significant waves (Ub > 0.01 m/s), wave-current interaction case
385! according to Soulsby (1995). Otherwise, tauCWmax = tauC for sediment
386! purposes.
387!-----------------------------------------------------------------------
388!
389 IF (ub(i,j).gt.0.01_r8) THEN
390!
391! Determine skin stresses (m2/s2) for pure waves and combined flow
392! using Soulsby (1995) approximation of the wave friction factor,
393! fw=2*scf1*(Znot/Ab)**scf2 so tauW=0.5fw*Ub**2.
394!
395 tau_w=scf1*((znotc*fwave_bot)**scf2)*(ub(i,j)**scf3)
396!
397! Wave-averaged, combined wave-current stress.(Eqn. 69, Soulsby, 1997).
398!
399 tau_cw=tau_cs* &
400 & (1.0_r8+scf4*((tau_w/(tau_w+tau_cs))**scf5))
401!
402! Maximum of combined wave-current skin stress (m2/s2) component for
403! sediment.(Eqn. 70, Soulsby, 1997).
404!
405 tau_cws=sqrt((tau_cw+tau_w*cos(phicw))**2+ &
406 & (tau_w*sin(phicw))**2)
407!! tauw(i,j)=tau_cws
408 taucwmax(i,j)=tau_cws
409 tauw(i,j)=tau_w
410!
411! Set combined stress for Znot.
412!
413 tau_w=scf1*((znot*fwave_bot)**scf2)*(ub(i,j)**scf3)
414 tau_cw=tauc(i,j)* &
415 & (1.0_r8+scf4*((tau_w/(tau_w+tauc(i,j)))**scf5))
416
417# ifdef MB_Z0BL
418!
419!-----------------------------------------------------------------------
420! Compute bedload roughness for ripple predictor and sediment purposes.
421! At high transport stages, friction depends on thickness of bedload
422! layer. Shear stress due to combined grain and bedload roughness is
423! used to predict ripples and onset of suspension (Li and Amos, 2001).
424!-----------------------------------------------------------------------
425!
426# ifdef MB_CALC_ZNOT
427 tau_ex=max((tau_cws-tau_cb),0.0_r8)
428 cff=(1.0_r8/((rhosed-1.0_r8)*g*d50))
429 znot_bl=17.4_r8*d50*(cff*tau_ex)**0.75_r8
430 znotc=znotc+znot_bl
431# endif
432!
433!-----------------------------------------------------------------------
434! Compute stresses (m2/s2)for sediment purposes, using grain and
435! bedload roughness.
436!-----------------------------------------------------------------------
437!
438 cff1=vonkar/log(zr(i,j)/znotc)
439 tau_c=cff1*cff1*umag(i,j)*umag(i,j)
440 tau_wb=scf1*((znotc*fwave_bot)**scf2)*(ub(i,j)**scf3)
441 tau_cw=tau_c*(1.0_r8+scf4*((tau_wb/(tau_wb+tau_c))**scf5))
442!
443! Maximum of combined wave-current stress (m2/s2) component for
444! sediment purposes.
445!
446 tau_cwb=sqrt((tau_cw+tau_wb*cos(phicw))**2+ &
447 & (tau_wb*sin(phicw))**2)
448 taucwmax(i,j)=tau_cwb
449 tauw(i,j)=tau_wb
450# endif
451
452# ifdef MB_Z0RIP
453!
454!-----------------------------------------------------------------------
455! Determine bedform roughness ripple height (m) and ripple length (m)
456! for sandy beds. Use structure according to Li and Amos (2001).
457!-----------------------------------------------------------------------
458!
459! Check median grain diameter
460!
461 IF (d50.ge.0.000063_r8) THEN
462!
463! Enhanced skin stress if pre-exisiting ripples (Nielsen, 1986).
464!
465 rhmax=0.8_r8*rlen/pi
466 rhgt=max(min(rhmax,rhgt),rhmin)
467 tau_en=max(tau_cws,tau_cws*(rlen/(rlen-pi*rhgt))**2)
468!
469 IF ((tau_cws.lt.tau_cb).and.(tau_en.ge.tau_cb)) THEN
470 rhgt=(19.6_r8*(sqrt(tau_cws/tau_cb))+20.9_r8)*d50
471 rlen=rhgt/0.12_r8 ! local transport
472 ELSE
473 IF ((tau_cws.ge.tau_cb).and.(tau_cwb.lt.tau_bf)) THEN
474 rhgt=(22.15_r8*(sqrt(tau_cwb/tau_cb))+6.38_r8)*d50
475 rlen=rhgt/0.12_r8 ! bed load in eq. range
476 ELSE
477 IF ((tau_cwb.ge.tau_bf).and.(tau_cwb.lt.tau_up)) THEN
478 rlen=535.0_r8*d50 ! break off regime
479 rhgt=0.15_r8*rlen*(sqrt(tau_up)-sqrt(tau_cwb))/ &
480 & (sqrt(tau_up)-sqrt(tau_bf ))
481 ELSE IF (tau_cwb.ge.tau_up) THEN
482 rlen=0.0_r8 ! sheet flow, plane bed
483 rhgt=0.0_r8
484 ELSE
485 rlen=bottom(i,j,irlen)
486 rhgt=bottom(i,j,irhgt)
487 END IF
488 END IF
489 END IF
490 END IF
491# endif
492
493# ifdef MB_Z0BIO
494!
495!-----------------------------------------------------------------------
496! Determine (biogenic) bedform roughness ripple height (m) and ripple
497! length (m) for silty beds, using Harris and Wiberg (2001).
498!-----------------------------------------------------------------------
499!
500! Use 10 cm default biogenic ripple length, RLbio (Wheatcroft 1994).
501!
502! NOTE: For mixed beds take average of silty and sandy bedform
503! roughnesses weighted according to silt and sand fractions.
504!
505 IF (d50.lt.0.000063_r8) THEN
506 rlbio=0.1_r8
507 thetw=tau_cws*(1.0_r8/((rhosed-1.0_r8)*g*d50))
508 rhbio=(thetw**(-1.67_r8))*rlbio*rhbiofac
509 rhgt=min(rhbio,rhbiomax)
510 rlen=rlbio
511 END IF
512# endif
513
514# if defined MB_Z0RIP || defined MB_Z0BIO
515!
516! Ripple roughness using Grant and Madsen (1982) roughness length.
517!
518# ifdef MB_CALC_ZNOT
519 znot_rip=0.92_r8*rhgt*rhgt/(max(rlen,rlmin))
520 znotc=znotc+znot_rip
521# endif
522!
523!-----------------------------------------------------------------------
524! Compute bottom stress (m2/s2) components based on total roughnes.
525!-----------------------------------------------------------------------
526!
527 cff1=vonkar/log(zr(i,j)/znotc)
528 tau_c=cff1*cff1*umag(i,j)*umag(i,j)
529 tau_w=scf1*((znotc*fwave_bot)**scf2)*(ub(i,j)**scf3)
530 tau_cw=tau_c* &
531 & (1.0_r8+scf4*((tau_w/(tau_w+tau_c))**scf5))
532!! tau_cwb=SQRT((tau_cw+tau_w*COS(phiCW))**2+ &
533!! & (tau_w*SIN(phiCW))**2)
534!! tauCWmax(i,j)=tau_cwb
535# endif
536!
537!-----------------------------------------------------------------------
538! Compute effective bottom shear velocity (m/s) relevant for flow and
539! eddy-diffusivities/viscosity.
540!-----------------------------------------------------------------------
541!
542!! tauC(i,j)=tau_cw
543 taucw(i,j)=tau_cw
544 tauw(i,j)=tau_w
545 ELSE IF (ub(i,j).le.0.01_r8) THEN
546!
547! If current-only, tauCWmax=tauC(skin) for use in sediment model; tauC
548! is determined using roughness due to current ripples.
549! In this limit, tauCWmax=tauCW=tauC
550!
551 taucwmax(i,j)=tauc(i,j)
552 tauw(i,j)=0.0_r8
553!! tauW(i,j)=tau_cs
554!! tauCW(i,j)=tau_cs
555
556# ifdef MB_Z0RIP
557!! IF (tauC(i,j).gt.tau_up) THEN
558 IF (tau_cs(i,j).gt.tau_up) THEN
559 rhgt=0.0_r8
560 rlen=0.0_r8
561!! ELSE IF (tauC(i,j).lt.tau_cb) THEN
562 ELSE IF (tau_cs(i,j).lt.tau_cb) THEN
563 rlen=bottom(i,j,irlen)
564 rhgt=bottom(i,j,irhgt)
565 ELSE
566 rlen=1000.0_r8*d50 ! Yalin (1964)
567 rhgt=0.0308_r8*(rlen**1.19_r8)
568!! rhgt=100.0_r8*0.074_r8* &
569!! & (0.01_r8*rlen)**1.19_r8 ! Allen (1970)
570 END IF
571# ifdef MB_CALC_ZNOT
572 znotc=znotc+0.92_r8*rhgt*rhgt/(max(rlen,rlmin))
573# endif
574# endif
575 cff1=vonkar/log(zr(i,j)/znotc)
576 cff2=min(cdb_max,max(cdb_min,cff1*cff1))
577!! tauc(i,j)=cff2*Umag(i,j)*Umag(i,j)
578 taucw(i,j)=cff2*umag(i,j)*umag(i,j)
579 END IF
580!
581!-----------------------------------------------------------------------
582! Load variables for output purposes.
583!-----------------------------------------------------------------------
584!
585 bottom(i,j,ibwav)=ab
586 bottom(i,j,irlen)=rlen
587 bottom(i,j,irhgt)=rhgt
588 bottom(i,j,izdef)=znot ! Zob(ng)
589 bottom(i,j,izapp)=znotc
590 END DO
591 END DO
592!
593!-----------------------------------------------------------------------
594! Compute kinematic bottom stress (m2/s2) components for flow due to
595! combined current and wind-induced waves.
596!-----------------------------------------------------------------------
597!
598 DO j=jstr,jend
599 DO i=istru,iend
600 anglec=ur_mb(i,j)/(0.5*(umag(i-1,j)+umag(i,j)))
601 bustr(i,j)=0.5_r8*(taucw(i-1,j)+taucw(i,j))*anglec
602# ifdef WET_DRY
603 cff2=0.75_r8*0.5_r8*(z_w(i-1,j,1)+z_w(i,j,1)- &
604 & z_w(i-1,j,0)-z_w(i,j,0))
605 bustr(i,j)=sign(1.0_r8,bustr(i,j))*min(abs(bustr(i,j)), &
606 & abs(u(i,j,1,nrhs))*cff2/dt(ng))
607# endif
608 END DO
609 END DO
610 DO j=jstrv,jend
611 DO i=istr,iend
612 anglec=vr_mb(i,j)/(0.5_r8*(umag(i,j-1)+umag(i,j)))
613 bvstr(i,j)=0.5_r8*(taucw(i,j-1)+taucw(i,j))*anglec
614# ifdef WET_DRY
615 cff2=0.75_r8*0.5_r8*(z_w(i,j-1,1)+z_w(i,j,1)- &
616 & z_w(i,j-1,0)-z_w(i,j,0))
617 bvstr(i,j)=sign(1.0_r8,bvstr(i,j))*min(abs(bvstr(i,j)), &
618 & abs(v(i,j,1,nrhs))*cff2/dt(ng))
619# endif
620 END DO
621 END DO
622 DO j=jstr,jend
623 DO i=istr,iend
624 anglec=ucur(i,j)/umag(i,j)
625 anglew=cos(1.5_r8*pi-dwave(i,j)-angler(i,j))
626 bustrc(i,j)=taucw(i,j)*anglec
627 bustrw(i,j)=tauw(i,j)*anglew
628 bustrcwmax(i,j)=taucwmax(i,j)*anglew
629 ubot(i,j)=ub(i,j)*anglew
630 ur(i,j)=ucur(i,j)
631!
632 anglec=vcur(i,j)/umag(i,j)
633 anglew=sin(1.5_r8*pi-dwave(i,j)-angler(i,j))
634 bvstrc(i,j)=taucw(i,j)*anglec
635 bvstrw(i,j)=tauw(i,j)*anglew
636 bvstrcwmax(i,j)=taucwmax(i,j)*anglew
637 vbot(i,j)=ub(i,j)*anglew
638 vr(i,j)=vcur(i,j)
639 END DO
640 END DO
641!
642! Apply periodic or gradient boundary conditions for output purposes.
643!
644 CALL bc_u2d_tile (ng, tile, &
645 & lbi, ubi, lbj, ubj, &
646 & bustr)
647 CALL bc_v2d_tile (ng, tile, &
648 & lbi, ubi, lbj, ubj, &
649 & bvstr)
650 CALL bc_r2d_tile (ng, tile, &
651 & lbi, ubi, lbj, ubj, &
652 & bustrc)
653 CALL bc_r2d_tile (ng, tile, &
654 & lbi, ubi, lbj, ubj, &
655 & bvstrc)
656 CALL bc_r2d_tile (ng, tile, &
657 & lbi, ubi, lbj, ubj, &
658 & bustrw)
659 CALL bc_r2d_tile (ng, tile, &
660 & lbi, ubi, lbj, ubj, &
661 & bvstrw)
662 CALL bc_u2d_tile (ng, tile, &
663 & lbi, ubi, lbj, ubj, &
664 & bustrcwmax)
665 CALL bc_r2d_tile (ng, tile, &
666 & lbi, ubi, lbj, ubj, &
667 & bvstrcwmax)
668 CALL bc_r2d_tile (ng, tile, &
669 & lbi, ubi, lbj, ubj, &
670 & ubot)
671 CALL bc_r2d_tile (ng, tile, &
672 & lbi, ubi, lbj, ubj, &
673 & vbot)
674 CALL bc_r2d_tile (ng, tile, &
675 & lbi, ubi, lbj, ubj, &
676 & ur)
677 CALL bc_r2d_tile (ng, tile, &
678 & lbi, ubi, lbj, ubj, &
679 & vr)
680 CALL bc_r2d_tile (ng, tile, &
681 & lbi, ubi, lbj, ubj, &
682 & bottom(:,:,ibwav))
683 CALL bc_r2d_tile (ng, tile, &
684 & lbi, ubi, lbj, ubj, &
685 & bottom(:,:,irhgt))
686 CALL bc_r2d_tile (ng, tile, &
687 & lbi, ubi, lbj, ubj, &
688 & bottom(:,:,irlen))
689 CALL bc_r2d_tile (ng, tile, &
690 & lbi, ubi, lbj, ubj, &
691 & bottom(:,:,izdef))
692 CALL bc_r2d_tile (ng, tile, &
693 & lbi, ubi, lbj, ubj, &
694 & bottom(:,:,izapp))
695
696#ifdef DISTRIBUTE
697 CALL mp_exchange2d (ng, tile, inlm, 4, &
698 & lbi, ubi, lbj, ubj, &
699 & nghostpoints, &
700 & ewperiodic(ng), nsperiodic(ng), &
701 & bustr, bvstr, bustrc, bvstrc)
702 CALL mp_exchange2d (ng, tile, inlm, 4, &
703 & lbi, ubi, lbj, ubj, &
704 & nghostpoints, &
705 & ewperiodic(ng), nsperiodic(ng), &
706 & bustrw, bvstrw, bustrcwmax, bvstrcwmax)
707 CALL mp_exchange2d (ng, tile, inlm, 4, &
708 & lbi, ubi, lbj, ubj, &
709 & nghostpoints, &
710 & ewperiodic(ng), nsperiodic(ng), &
711 & ubot, vbot, ur, vr)
712 CALL mp_exchange2d (ng, tile, inlm, 3, &
713 & lbi, ubi, lbj, ubj, &
714 & nghostpoints, &
715 & ewperiodic(ng), nsperiodic(ng), &
716 & bottom(:,:,ibwav), &
717 & bottom(:,:,irhgt), &
718 & bottom(:,:,irlen))
719 CALL mp_exchange2d (ng, tile, inlm, 2, &
720 & lbi, ubi, lbj, ubj, &
721 & nghostpoints, &
722 & ewperiodic(ng), nsperiodic(ng), &
723 & bottom(:,:,izdef), &
724 & bottom(:,:,izapp))
725#endif
726!
727 RETURN
728 END SUBROUTINE mb_bbl_tile
729
730 END MODULE bbl_mod
bbl_mod
Definition mb_bbl.h:1
bbl_mod::bblm
subroutine, public bblm(ng, tile)
Definition mb_bbl.h:55
bbl_mod::mb_bbl_tile
subroutine mb_bbl_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, h, z_r, z_w, angler, zobot, hwave, uwave_rms, dwave, pwave_bot, rho, u, v, bottom, ubot, vbot, ur, vr, bustrc, bvstrc, bustrw, bvstrw, bustrcwmax, bvstrcwmax, bustr, bvstr)
Definition mb_bbl.h:136
bc_2d_mod
Definition bc_2d.F:2
bc_2d_mod::bc_r2d_tile
subroutine bc_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:44
bc_2d_mod::bc_v2d_tile
subroutine bc_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:345
bc_2d_mod::bc_u2d_tile
subroutine bc_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:167
mod_bbl
Definition mod_bbl.F:2
mod_bbl::bbl
type(t_bbl), dimension(:), allocatable bbl
Definition mod_bbl.F:62
mod_forces
Definition mod_forces.F:2
mod_forces::forces
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
mod_grid
Definition mod_grid.F:2
mod_grid::grid
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
mod_ocean
Definition mod_ocean.F:2
mod_ocean::ocean
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
mod_param
Definition mod_param.F:2
mod_param::inlm
integer, parameter inlm
Definition mod_param.F:662
mod_param::nghostpoints
integer nghostpoints
Definition mod_param.F:710
mod_scalars
Definition mod_scalars.F:2
mod_scalars::cdb_min
real(dp) cdb_min
Definition mod_scalars.F:772
mod_scalars::vonkar
real(dp) vonkar
Definition mod_scalars.F:469
mod_scalars::dt
real(dp), dimension(:), allocatable dt
Definition mod_scalars.F:269
mod_scalars::ewperiodic
logical, dimension(:), allocatable ewperiodic
Definition mod_scalars.F:510
mod_scalars::nsperiodic
logical, dimension(:), allocatable nsperiodic
Definition mod_scalars.F:511
mod_scalars::g
real(dp) g
Definition mod_scalars.F:464
mod_scalars::cdb_max
real(dp) cdb_max
Definition mod_scalars.F:773
mod_scalars::pi
real(dp), parameter pi
Definition mod_scalars.F:813
mod_sedbed
Definition sedbed_mod.h:1
mod_sedbed::sedbed
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
mod_sediment
Definition sediment_mod.h:1
mod_sediment::iwsed
integer, parameter iwsed
Definition sediment_mod.h:126
mod_sediment::irlen
integer, parameter irlen
Definition sediment_mod.h:128
mod_sediment::irhgt
integer, parameter irhgt
Definition sediment_mod.h:129
mod_sediment::ibwav
integer, parameter ibwav
Definition sediment_mod.h:130
mod_sediment::idens
integer, parameter idens
Definition sediment_mod.h:125
mod_sediment::isd50
integer, parameter isd50
Definition sediment_mod.h:124
mod_sediment::izapp
integer, parameter izapp
Definition sediment_mod.h:132
mod_sediment::izdef
integer, parameter izdef
Definition sediment_mod.h:131
mod_sediment::itauc
integer, parameter itauc
Definition sediment_mod.h:127
mod_stepping
Definition mod_stepping.F:2
mod_stepping::nrhs
integer, dimension(:), allocatable nrhs
Definition mod_stepping.F:70
mp_exchange_mod
Definition mp_exchange.F:2
mp_exchange_mod::mp_exchange2d
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
Definition mp_exchange.F:294
wclock_off
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
wclock_on
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3