ROMS
Loading...
Searching...
No Matches
sed_bedload.F File Reference
#include "cppdefs.h"
#include "tile.h"
#include "set_bounds.h"
Include dependency graph for sed_bedload.F:

Go to the source code of this file.

Macros

#define SLOPE_NEMETH
 
#define BSTRESS_UPWIND
 

Functions/Subroutines

program __sed_bedload_f__
 
subroutine sed_bedload (ng, tile)
 
subroutine sed_bedload_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, pm, pn, rmask, umask, vmask, rmask_wet, z_w, bustrc, bvstrc, bustrw, bvstrw, bustrcwmax, bvstrcwmax, dwave, pwave_bot, bustr, bvstr, t, hwave, lwave, angler, bed_thick, h, om_r, om_u, on_r, on_v, bedldu, bedldv, bed, bed_frac, bed_mass, bottom)
 

Macro Definition Documentation

◆ BSTRESS_UPWIND

#define BSTRESS_UPWIND

◆ SLOPE_NEMETH

#define SLOPE_NEMETH

Function/Subroutine Documentation

◆ __sed_bedload_f__()

program __sed_bedload_f__

Definition at line 2 of file sed_bedload.F.

References sed_bedload().

Here is the call graph for this function:

◆ sed_bedload()

subroutine __sed_bedload_f__::sed_bedload ( integer, intent(in) ng,
integer, intent(in) tile )
private

Definition at line 46 of file sed_bedload.F.

47!***********************************************************************
48!
49 USE mod_param
50 USE mod_forces
51 USE mod_grid
52 USE mod_ocean
53 USE mod_sedbed
54 USE mod_stepping
55# ifdef BBL_MODEL
56 USE mod_bbl
57# endif
58!
59! Imported variable declarations.
60!
61 integer, intent(in) :: ng, tile
62!
63! Local variable declarations.
64!
65 character (len=*), parameter :: MyFile = &
66 & __FILE__
67!
68# include "tile.h"
69!
70# ifdef PROFILE
71 CALL wclock_on (ng, inlm, 16, __line__, myfile)
72# endif
73 CALL sed_bedload_tile (ng, tile, &
74 & lbi, ubi, lbj, ubj, &
75 & imins, imaxs, jmins, jmaxs, &
76 & nstp(ng), nnew(ng), &
77 & grid(ng) % pm, &
78 & grid(ng) % pn, &
79# ifdef MASKING
80 & grid(ng) % rmask, &
81 & grid(ng) % umask, &
82 & grid(ng) % vmask, &
83# endif
84# ifdef WET_DRY
85 & grid(ng) % rmask_wet, &
86# endif
87 & grid(ng) % z_w, &
88# ifdef BBL_MODEL
89 & bbl(ng) % bustrc, &
90 & bbl(ng) % bvstrc, &
91 & bbl(ng) % bustrw, &
92 & bbl(ng) % bvstrw, &
93 & bbl(ng) % bustrcwmax, &
94 & bbl(ng) % bvstrcwmax, &
95 & forces(ng) % Dwave, &
96 & forces(ng) % Pwave_bot, &
97# endif
98 & forces(ng) % bustr, &
99 & forces(ng) % bvstr, &
100 & ocean(ng) % t, &
101# if defined BEDLOAD_SOULSBY
102 & forces(ng) % Hwave, &
103 & forces(ng) % Lwave, &
104 & grid(ng) % angler, &
105# endif
106# if defined SED_MORPH
107 & sedbed(ng) % bed_thick, &
108# endif
109# if defined BEDLOAD_MPM || defined BEDLOAD_SOULSBY
110 & grid(ng) % h, &
111 & grid(ng) % om_r, &
112 & grid(ng) % om_u, &
113 & grid(ng) % on_r, &
114 & grid(ng) % on_v, &
115 & sedbed(ng) % bedldu, &
116 & sedbed(ng) % bedldv, &
117# endif
118 & sedbed(ng) % bed, &
119 & sedbed(ng) % bed_frac, &
120 & sedbed(ng) % bed_mass, &
121 & sedbed(ng) % bottom)
122# ifdef PROFILE
123 CALL wclock_off (ng, inlm, 16, __line__, myfile)
124# endif
125!
126 RETURN
type(t_bbl), dimension(:), allocatable bbl
Definition mod_bbl.F:62
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter inlm
Definition mod_param.F:662
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine sed_bedload_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, pm, pn, rmask, umask, vmask, rmask_wet, z_w, bustrc, bvstrc, bustrw, bvstrw, bustrcwmax, bvstrcwmax, dwave, pwave_bot, bustr, bvstr, t, hwave, lwave, angler, bed_thick, h, om_r, om_u, on_r, on_v, bedldu, bedldv, bed, bed_frac, bed_mass, bottom)
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3

References mod_bbl::bbl, mod_forces::forces, mod_grid::grid, mod_param::inlm, mod_stepping::nnew, mod_stepping::nstp, mod_ocean::ocean, sed_bedload_tile(), mod_sedbed::sedbed, wclock_off(), and wclock_on().

Referenced by __sed_bedload_f__(), and sediment_mod::sediment().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ sed_bedload_tile()

subroutine __sed_bedload_f__::sed_bedload_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nstp,
integer, intent(in) nnew,
real(r8), dimension(lbi:,lbj:), intent(in) pm,
real(r8), dimension(lbi:,lbj:), intent(in) pn,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:), intent(in) rmask_wet,
real(r8), dimension(lbi:,lbj:,0:), intent(in) z_w,
real(r8), dimension(lbi:,lbj:), intent(in) bustrc,
real(r8), dimension(lbi:,lbj:), intent(in) bvstrc,
real(r8), dimension(lbi:,lbj:), intent(in) bustrw,
real(r8), dimension(lbi:,lbj:), intent(in) bvstrw,
real(r8), dimension(lbi:,lbj:), intent(in) bustrcwmax,
real(r8), dimension(lbi:,lbj:), intent(in) bvstrcwmax,
real(r8), dimension(lbi:,lbj:), intent(in) dwave,
real(r8), dimension(lbi:,lbj:), intent(in) pwave_bot,
real(r8), dimension(lbi:,lbj:), intent(in) bustr,
real(r8), dimension(lbi:,lbj:), intent(in) bvstr,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) t,
real(r8), dimension(lbi:,lbj:), intent(in) hwave,
real(r8), dimension(lbi:,lbj:), intent(in) lwave,
real(r8), dimension(lbi:,lbj:), intent(in) angler,
real(r8), dimension(lbi:,lbj:,:), intent(inout) bed_thick,
real(r8), dimension(lbi:,lbj:), intent(in) h,
real(r8), dimension(lbi:,lbj:), intent(in) om_r,
real(r8), dimension(lbi:,lbj:), intent(in) om_u,
real(r8), dimension(lbi:,lbj:), intent(in) on_r,
real(r8), dimension(lbi:,lbj:), intent(in) on_v,
real(r8), dimension(lbi:,lbj:,:), intent(inout) bedldu,
real(r8), dimension(lbi:,lbj:,:), intent(inout) bedldv,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) bed,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) bed_frac,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) bed_mass,
real(r8), dimension(lbi:,lbj:,:), intent(inout) bottom )
private

Definition at line 130 of file sed_bedload.F.

163!***********************************************************************
164!
165 USE mod_param
166 USE mod_ncparam
167 USE mod_scalars
168 USE mod_sediment
169!
170 USE bc_3d_mod, ONLY : bc_r3d_tile
171# ifdef BEDLOAD
173# endif
174# ifdef DISTRIBUTE
176# endif
177!
178! Imported variable declarations.
179!
180 integer, intent(in) :: ng, tile
181 integer, intent(in) :: LBi, UBi, LBj, UBj
182 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
183 integer, intent(in) :: nstp, nnew
184!
185# ifdef ASSUMED_SHAPE
186 real(r8), intent(in) :: pm(LBi:,LBj:)
187 real(r8), intent(in) :: pn(LBi:,LBj:)
188# ifdef MASKING
189 real(r8), intent(in) :: rmask(LBi:,LBj:)
190 real(r8), intent(in) :: umask(LBi:,LBj:)
191 real(r8), intent(in) :: vmask(LBi:,LBj:)
192# endif
193# ifdef WET_DRY
194 real(r8), intent(in) :: rmask_wet(LBi:,LBj:)
195# endif
196 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
197# ifdef BBL_MODEL
198 real(r8), intent(in) :: bustrc(LBi:,LBj:)
199 real(r8), intent(in) :: bvstrc(LBi:,LBj:)
200 real(r8), intent(in) :: bustrw(LBi:,LBj:)
201 real(r8), intent(in) :: bvstrw(LBi:,LBj:)
202 real(r8), intent(in) :: bustrcwmax(LBi:,LBj:)
203 real(r8), intent(in) :: bvstrcwmax(LBi:,LBj:)
204 real(r8), intent(in) :: Dwave(LBi:,LBj:)
205 real(r8), intent(in) :: Pwave_bot(LBi:,LBj:)
206# endif
207 real(r8), intent(in) :: bustr(LBi:,LBj:)
208 real(r8), intent(in) :: bvstr(LBi:,LBj:)
209# if defined BEDLOAD_SOULSBY
210 real(r8), intent(in) :: Hwave(LBi:,LBj:)
211 real(r8), intent(in) :: Lwave(LBi:,LBj:)
212 real(r8), intent(in) :: angler(LBi:,LBj:)
213# endif
214# if defined SED_MORPH
215 real(r8), intent(inout):: bed_thick(LBi:,LBj:,:)
216# endif
217# if defined BEDLOAD_MPM || defined BEDLOAD_SOULSBY
218 real(r8), intent(in) :: h(LBi:,LBj:)
219 real(r8), intent(in) :: om_r(LBi:,LBj:)
220 real(r8), intent(in) :: om_u(LBi:,LBj:)
221 real(r8), intent(in) :: on_r(LBi:,LBj:)
222 real(r8), intent(in) :: on_v(LBi:,LBj:)
223 real(r8), intent(inout) :: bedldu(LBi:,LBj:,:)
224 real(r8), intent(inout) :: bedldv(LBi:,LBj:,:)
225# endif
226 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
227 real(r8), intent(inout) :: bed(LBi:,LBj:,:,:)
228 real(r8), intent(inout) :: bed_frac(LBi:,LBj:,:,:)
229 real(r8), intent(inout) :: bed_mass(LBi:,LBj:,:,:,:)
230 real(r8), intent(inout) :: bottom(LBi:,LBj:,:)
231# else
232 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
233 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
234# ifdef MASKING
235 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
236 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
237 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
238# endif
239# ifdef WET_DRY
240 real(r8), intent(in) :: rmask_wet(LBi:UBi,LBj:UBj)
241# endif
242 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
243# ifdef BBL_MODEL
244 real(r8), intent(in) :: bustrc(LBi:UBi,LBj:UBj)
245 real(r8), intent(in) :: bvstrc(LBi:UBi,LBj:UBj)
246 real(r8), intent(in) :: bustrw(LBi:UBi,LBj:UBj)
247 real(r8), intent(in) :: bvstrw(LBi:UBi,LBj:UBj)
248 real(r8), intent(in) :: bustrcwmax(LBi:UBi,LBj:UBj)
249 real(r8), intent(in) :: bvstrcwmax(LBi:UBi,LBj:UBj)
250 real(r8), intent(in) :: Dwave(LBi:UBi,LBj:UBj)
251 real(r8), intent(in) :: Pwave_bot(LBi:UBi,LBj:UBj)
252# endif
253 real(r8), intent(in) :: bustr(LBi:UBi,LBj:UBj)
254 real(r8), intent(in) :: bvstr(LBi:UBi,LBj:UBj)
255# if defined BEDLOAD_SOULSBY
256 real(r8), intent(in) :: Hwave(LBi:UBi,LBj:UBj)
257 real(r8), intent(in) :: Lwave(LBi:UBi,LBj:UBj)
258 real(r8), intent(in) :: angler(LBi:UBi,LBj:UBj)
259# endif
260# if defined SED_MORPH
261 real(r8), intent(inout):: bed_thick(LBi:UBi,LBj:UBj,3)
262# endif
263# if defined BEDLOAD_MPM || defined BEDLOAD_SOULSBY
264 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
265 real(r8), intent(in) :: om_r(LBi:UBi,LBj:UBj)
266 real(r8), intent(in) :: om_u(LBi:UBi,LBj:UBj)
267 real(r8), intent(in) :: on_r(LBi:UBi,LBj:UBj)
268 real(r8), intent(in) :: on_v(LBi:UBi,LBj:UBj)
269 real(r8), intent(inout) :: bedldu(LBi:UBi,LBj:UBj,NST)
270 real(r8), intent(inout) :: bedldv(LBi:UBi,LBj:UBj,NST)
271# endif
272 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
273 real(r8), intent(inout) :: bed(LBi:UBi,LBj:UBj,Nbed,MBEDP)
274 real(r8), intent(inout) :: bed_frac(LBi:UBi,LBj:UBj,Nbed,NST)
275 real(r8), intent(inout) :: bed_mass(LBi:UBi,LBj:UBj,Nbed,1:2,NST)
276 real(r8), intent(inout) :: bottom(LBi:UBi,LBj:UBj,MBOTP)
277# endif
278!
279! Local variable declarations.
280!
281 integer :: i, ised, j, k
282
283 real(r8), parameter :: eps = 1.0e-14_r8
284
285 real(r8) :: cff, cff1, cff2, cff3, cff4, cff5
286
287 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tau_w
288# ifdef BSTRESS_UPWIND
289 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tau_wX
290 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tau_wE
291# endif
292# ifdef BEDLOAD
293 real(r8) :: a_slopex, a_slopey, sed_angle
294 real(r8) :: bedld, bedld_mass, dzdx, dzdy
295 real(r8) :: smgd, smgdr, osmgd, Umag
296 real(r8) :: rhs_bed, Ua, Ra, phi, Clim
297
298 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
299 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
300 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX_r
301 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE_r
302# endif
303# if defined BEDLOAD_MPM
304 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: angleu
305 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: anglev
306# endif
307# if defined BEDLOAD_SOULSBY
308 real(r8) :: theta_mean, theta_wav, w_asym
309 real(r8) :: theta_max, theta_max1, theta_max2
310 real(r8) :: phi_x1, phi_x2, phi_x, phi_y, Dstp
311 real(r8) :: bedld_x, bedld_y, tau_cur, waven, wavec
312
313 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: phic
314 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: phicw
315 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tau_wav
316 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tau_mean
317 real(r8), parameter :: kdmax = 100.0_r8
318# endif
319
320# include "set_bounds.h"
321!
322!-----------------------------------------------------------------------
323! Compute maximum bottom stress for MPM bedload or suspended load.
324!-----------------------------------------------------------------------
325!
326# if defined BEDLOAD_MPM || defined SUSPLOAD
327# ifdef BBL_MODEL
328 DO j=jstr-1,jend+1
329 DO i=istr-1,iend+1
330 tau_w(i,j)=sqrt(bustrcwmax(i,j)*bustrcwmax(i,j)+ &
331 & bvstrcwmax(i,j)*bvstrcwmax(i,j))
332# ifdef WET_DRY
333 tau_w(i,j)=tau_w(i,j)*rmask_wet(i,j)
334# endif
335 END DO
336 END DO
337# else
338 DO j=jstrm1,jendp1
339 DO i=istrm1,iendp1
340# ifdef BSTRESS_UPWIND
341 cff1=0.5_r8*(1.0_r8+sign(1.0_r8,bustr(i+1,j)))
342 cff2=0.5_r8*(1.0_r8-sign(1.0_r8,bustr(i+1,j)))
343 cff3=0.5_r8*(1.0_r8+sign(1.0_r8,bustr(i ,j)))
344 cff4=0.5_r8*(1.0_r8-sign(1.0_r8,bustr(i ,j)))
345 tau_wx(i,j)=cff3*(cff1*bustr(i,j)+ &
346 & cff2*0.5_r8*(bustr(i,j)+bustr(i+1,j)))+ &
347 & cff4*(cff2*bustr(i+1,j)+ &
348 & cff1*0.5_r8*(bustr(i,j)+bustr(i+1,j)))
349 cff1=0.5_r8*(1.0_r8+sign(1.0_r8,bvstr(i,j+1)))
350 cff2=0.5_r8*(1.0_r8-sign(1.0_r8,bvstr(i,j+1)))
351 cff3=0.5_r8*(1.0_r8+sign(1.0_r8,bvstr(i,j)))
352 cff4=0.5_r8*(1.0_r8-sign(1.0_r8,bvstr(i,j)))
353 tau_we(i,j)=cff3*(cff1*bvstr(i,j)+ &
354 & cff2*0.5_r8*(bvstr(i,j)+bvstr(i,j+1)))+ &
355 & cff4*(cff2*bvstr(i,j+1)+ &
356 & cff1*0.5_r8*(bvstr(i,j)+bvstr(i,j+1)))
357# endif
358 tau_w(i,j)=0.5_r8*sqrt((bustr(i,j)+bustr(i+1,j))* &
359 & (bustr(i,j)+bustr(i+1,j))+ &
360 & (bvstr(i,j)+bvstr(i,j+1))* &
361 & (bvstr(i,j)+bvstr(i,j+1)))
362# ifdef WET_DRY
363 tau_w(i,j)=tau_w(i,j)*rmask_wet(i,j)
364# endif
365 END DO
366 END DO
367# endif
368# endif
369
370# ifdef BEDLOAD
371!
372!-----------------------------------------------------------------------
373! Compute bedload sediment transport.
374!-----------------------------------------------------------------------
375!
376! Compute some constant bed slope parameters.
377!
378 sed_angle=dtan(33.0_r8*pi/180.0_r8)
379!
380! Compute angle between currents and waves (radians).
381!
382 DO j=jstrm1,jendp1
383 DO i=istrm1,iendp1
384# if defined BEDLOAD_SOULSBY
385!
386! Compute angle between currents and waves, measure CCW from current
387! direction toward wave vector.
388!
389 IF (bustrc(i,j).eq.0.0_r8) THEN
390 phic(i,j)=0.5_r8*pi*sign(1.0_r8,bvstrc(i,j))
391 ELSE
392 phic(i,j)=atan2(bvstrc(i,j),bustrc(i,j))
393 ENDIF
394 phicw(i,j)=1.5_r8*pi-dwave(i,j)-phic(i,j)-angler(i,j)
395!
396! Compute stress components at rho points.
397!
398 tau_cur=sqrt(bustrc(i,j)*bustrc(i,j)+ &
399 & bvstrc(i,j)*bvstrc(i,j))
400 tau_wav(i,j)=sqrt(bustrw(i,j)*bustrw(i,j)+ &
401 & bvstrw(i,j)*bvstrw(i,j))
402 tau_mean(i,j)=tau_cur*(1.0_r8+1.2_r8*((tau_wav(i,j)/ &
403 & (tau_cur+tau_wav(i,j)+eps))**3.2_r8))
404!
405# elif defined BEDLOAD_MPM
406 cff1=0.5_r8*(bustr(i,j)+bustr(i+1,j))
407 cff2=0.5_r8*(bvstr(i,j)+bvstr(i,j+1))
408 umag=sqrt(cff1*cff1+cff2*cff2)+eps
409 angleu(i,j)=cff1/umag
410 anglev(i,j)=cff2/umag
411# endif
412 END DO
413 END DO
414!
415 DO ised=ncs+1,nst
416 smgd=(srho(ised,ng)/rho0-1.0_r8)*g*sd50(ised,ng)
417 osmgd=1.0_r8/smgd
418 smgdr=sqrt(smgd)*sd50(ised,ng)*srho(ised,ng)
419!
420 DO j=jstrm1,jendp1
421 DO i=istrm1,iendp1
422# ifdef BEDLOAD_SOULSBY
423!
424! Compute wave asymmetry factor, based on Fredosoe and Deigaard.
425!
426 dstp=z_w(i,j,n(ng))+h(i,j)
427 waven=2.0_r8*pi/(lwave(i,j)+eps)
428 wavec=sqrt(g/waven*tanh(waven*dstp))
429 cff4=min(waven*dstp,kdmax)
430 cff1=-0.1875_r8*wavec*(waven*dstp)**2/(sinh(cff4))**4
431 cff2=0.125_r8*g*hwave(i,j)**2/(wavec*dstp+eps)
432 cff3=pi*hwave(i,j)/(pwave_bot(i,j)*sinh(cff4)+eps)
433 w_asym=max(min((cff1-cff2)/cff3,0.2_r8),0.0_r8)
434 w_asym=0.0_r8
435!
436! Compute nondimensional stresses.
437!
438 theta_wav=tau_wav(i,j)*osmgd+eps
439 theta_mean=tau_mean(i,j)*osmgd
440!
441 cff1=theta_wav*(1.0_r8+w_asym)
442 cff2=theta_wav*(1.0_r8-w_asym)
443 theta_max1=sqrt((theta_mean+ &
444 & cff1*cos(phicw(i,j)))**2+ &
445 & (cff1*sin(phicw(i,j)))**2)
446 theta_max2=sqrt((theta_mean+ &
447 & cff2*cos(phicw(i,j)+pi))**2+ &
448 & (cff2*sin(phicw(i,j)+pi))**2)
449 theta_max=max(theta_max1,theta_max2)
450!
451! Motion initiation factor.
452!
453 cff3=0.5_r8*(1.0_r8+sign(1.0_r8, &
454 & theta_max/tau_ce(ised,ng)-1.0_r8))
455!
456! Calculate bed loads in direction of current and perpendicular
457! direction.
458!
459 phi_x1=12.0_r8*sqrt(theta_mean)* &
460 & max((theta_mean-tau_ce(ised,ng)),0.0_r8)
461 phi_x2=12.0_r8*(0.9534_r8+0.1907*cos(2.0_r8*phicw(i,j)))* &
462 & sqrt(theta_wav)*theta_mean+ &
463 & 12.0_r8*(0.229_r8*w_asym*theta_wav**1.5_r8* &
464 & cos(phicw(i,j)))
465! phi_x=MAX(phi_x1,phi_x2) ! <- original
466 IF (abs(phi_x2).gt.phi_x1) THEN
467 phi_x=phi_x2
468 ELSE
469 phi_x=phi_x1
470 END IF
471 bedld_x=phi_x*smgdr*cff3
472!
473 cff5=theta_wav**1.5_r8+1.5_r8*(theta_mean**1.5_r8)
474 phi_y=12.0_r8*0.1907_r8*theta_wav*theta_wav* &
475 & (theta_mean*sin(2.0_r8*phicw(i,j))+1.2_r8*w_asym* &
476 & theta_wav*sin(phicw(i,j)))/cff5*cff3
477 bedld_y=phi_y*smgdr
478!
479! Partition bedld into xi and eta directions, still at rho points.
480! (FX_r and FE_r have dimensions of kg).
481!
482 fx_r(i,j)=(bedld_x*cos(phic(i,j))-bedld_y*sin(phic(i,j)))* &
483 & on_r(i,j)*dt(ng)
484 fe_r(i,j)=(bedld_x*sin(phic(i,j))+bedld_y*cos(phic(i,j)))* &
485 & om_r(i,j)*dt(ng)
486# elif defined BEDLOAD_MPM
487# ifdef BSTRESS_UPWIND
488!
489! Magnitude of bed load at rho points. Meyer-Peter Muller formulation.
490! bedld has dimensions of kg m-1 s-1. Use partitions of stress
491! from upwind direction, still at rho points.
492! (FX_r and FE_r have dimensions of kg).
493!
494 bedld=8.0_r8*(max((abs(tau_wx(i,j))*osmgd-0.047_r8), &
495 & 0.0_r8)**1.5_r8)*smgdr* &
496 & sign(1.0_r8,tau_wx(i,j))
497 fx_r(i,j)=bedld*on_r(i,j)*dt(ng)
498 bedld=8.0_r8*(max((abs(tau_we(i,j))*osmgd-0.047_r8), &
499 & 0.0_r8)**1.5_r8)*smgdr* &
500 & sign(1.0_r8,tau_we(i,j))
501 fe_r(i,j)=bedld*om_r(i,j)*dt(ng)
502# else
503!
504! Magnitude of bed load at rho points. Meyer-Peter Muller formulation.
505! (BEDLD has dimensions of kg m-1 s-1).
506!
507 bedld=8.0_r8*(max((tau_w(i,j)*osmgd-0.047_r8), &
508 & 0.0_r8)**1.5_r8)*smgdr
509!
510! Partition bedld into xi and eta directions, still at rho points.
511! (FX_r and FE_r have dimensions of kg).
512!
513 fx_r(i,j)=angleu(i,j)*bedld*on_r(i,j)*dt(ng)
514 fe_r(i,j)=anglev(i,j)*bedld*om_r(i,j)*dt(ng)
515# endif
516# endif
517!
518! Correct for along-direction slope. Limit slope to 0.9*sed angle.
519!
520 cff1=0.5_r8*(1.0_r8+sign(1.0_r8,fx_r(i,j)))
521 cff2=0.5_r8*(1.0_r8-sign(1.0_r8,fx_r(i,j)))
522# if defined SLOPE_NEMETH
523 dzdx=(h(i+1,j)-h(i ,j))/om_u(i+1,j)*cff1+ &
524 & (h(i-1,j)-h(i ,j))/om_u(i ,j)*cff2
525 dzdy=(h(i,j+1)-h(i,j ))/on_v(i,j+1)*cff1+ &
526 & (h(i,j-1)-h(i,j ))/on_v(i ,j)*cff2
527# ifdef BEDLOAD_MPM
528 cff=abs(tau_w(i,j))
529# else
530 cff=abs(tau_mean(i,j))
531# endif
532 a_slopex=0.3_r8*cff**0.5_r8*0.002_r8*dzdx+ &
533 & 0.3_r8*cff**1.5_r8*3.330_r8*dzdx
534 a_slopey=0.3_r8*cff**0.5_r8*0.002_r8*dzdy+ &
535 & 0.3_r8*cff**1.5_r8*3.330_r8*dzdy
536!
537! Add contriubiton of bed slope to bed load transport fluxes.
538!
539 fx_r(i,j)=fx_r(i,j)+a_slopex
540 fe_r(i,j)=fe_r(i,j)+a_slopey
541# elif defined SLOPE_LESSER
542 dzdx=min(((h(i+1,j)-h(i ,j))/om_u(i+1,j)*cff1+ &
543 & (h(i ,j)-h(i-1,j))/om_u(i ,j)*cff2),0.52_r8)* &
544 & sign(1.0_r8,fx_r(i,j))
545 dzdy=min(((h(i,j+1)-h(i,j ))/on_v(i,j+1)*cff1+ &
546 & (h(i,j )-h(i,j-1))/on_v(i ,j)*cff2),0.52_r8)* &
547 & sign(1.0_r8,fe_r(i,j))
548 cff=datan(dzdx)
549 a_slopex=sed_angle/(cos(cff)*(sed_angle-dzdx))
550 cff=datan(dzdy)
551 a_slopey=sed_angle/(cos(cff)*(sed_angle-dzdy))
552!
553! Add contriubiton of bed slope to bed load transport fluxes.
554!
555 fx_r(i,j)=fx_r(i,j)*a_slopex
556 fe_r(i,j)=fe_r(i,j)*a_slopey
557# endif
558!
559!
560# ifdef SED_MORPH
561!
562! Apply morphology factor.
563!
564 fx_r(i,j)=fx_r(i,j)*morph_fac(ised,ng)
565 fe_r(i,j)=fe_r(i,j)*morph_fac(ised,ng)
566
567# endif
568!
569! Apply bedload transport rate coefficient. Also limit
570! bedload to the fraction of each sediment class.
571!
572 fx_r(i,j)=fx_r(i,j)*bedload_coeff(ng)*bed_frac(i,j,1,ised)
573 fe_r(i,j)=fe_r(i,j)*bedload_coeff(ng)*bed_frac(i,j,1,ised)
574!
575! Limit bed load to available bed mass.
576!
577 bedld_mass=abs(fx_r(i,j))+abs(fe_r(i,j))+eps
578 fx_r(i,j)=min(abs(fx_r(i,j)), &
579 & bed_mass(i,j,1,nstp,ised)* &
580 & om_r(i,j)*on_r(i,j)*abs(fx_r(i,j))/ &
581 & bedld_mass)* &
582 & sign(1.0_r8,fx_r(i,j))
583 fe_r(i,j)=min(abs(fe_r(i,j)), &
584 & bed_mass(i,j,1,nstp,ised)* &
585 & om_r(i,j)*on_r(i,j)*abs(fe_r(i,j))/ &
586 & bedld_mass)* &
587 & sign(1.0_r8,fe_r(i,j))
588 END DO
589 END DO
590!
591! Apply boundary conditions (gradient).
592!
593 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
594 IF (domain(ng)%Western_Edge(tile)) THEN
595 DO j=jstrm1,jendp1
596 fx_r(istr-1,j)=fx_r(istr,j)
597 fe_r(istr-1,j)=fe_r(istr,j)
598 END DO
599 END IF
600 END IF
601 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
602 IF (domain(ng)%Eastern_Edge(tile)) THEN
603 DO j=jstrm1,jendp1
604 fx_r(iend+1,j)=fx_r(iend,j)
605 fe_r(iend+1,j)=fe_r(iend,j)
606 END DO
607 END IF
608 END IF
609!
610 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
611 IF (domain(ng)%Southern_Edge(tile)) THEN
612 DO i=istrm1,iendp1
613 fx_r(i,jstr-1)=fx_r(i,jstr)
614 fe_r(i,jstr-1)=fe_r(i,jstr)
615 END DO
616 END IF
617 END IF
618 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
619 IF (domain(ng)%Northern_Edge(tile)) THEN
620 DO i=istrm1,iendp1
621 fx_r(i,jend+1)=fx_r(i,jend)
622 fe_r(i,jend+1)=fe_r(i,jend)
623 END DO
624 END IF
625 END IF
626!
627 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
628 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
629 IF (domain(ng)%SouthWest_Corner(tile)) THEN
630 fx_r(istr-1,jstr-1)=0.5_r8*(fx_r(istr ,jstr-1)+ &
631 & fx_r(istr-1,jstr ))
632 fe_r(istr-1,jstr-1)=0.5_r8*(fe_r(istr ,jstr-1)+ &
633 & fe_r(istr-1,jstr ))
634 END IF
635 END IF
636
637 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
638 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
639 IF (domain(ng)%SouthEast_Corner(tile)) THEN
640 fx_r(iend+1,jstr-1)=0.5_r8*(fx_r(iend ,jstr-1)+ &
641 & fx_r(iend+1,jstr ))
642 fe_r(iend+1,jstr-1)=0.5_r8*(fe_r(iend ,jstr-1)+ &
643 & fe_r(iend+1,jstr ))
644 END IF
645 END IF
646
647 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
648 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
649 IF (domain(ng)%NorthWest_Corner(tile)) THEN
650 fx_r(istr-1,jend+1)=0.5_r8*(fx_r(istr-1,jend )+ &
651 & fx_r(istr ,jend+1))
652 fe_r(istr-1,jend+1)=0.5_r8*(fe_r(istr-1,jend )+ &
653 & fe_r(istr ,jend+1))
654 END IF
655 END IF
656
657 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
658 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
659 IF (domain(ng)%NorthEast_Corner(tile)) THEN
660 fx_r(iend+1,jend+1)=0.5_r8*(fx_r(iend+1,jend )+ &
661 & fx_r(iend ,jend+1))
662 fe_r(iend+1,jend+1)=0.5_r8*(fe_r(iend+1,jend )+ &
663 & fe_r(iend ,jend+1))
664 END IF
665 END IF
666!
667! Upwind shift FX_r and FE_r to u and v points.
668!
669 DO j=jstr-1,jend+1
670 DO i=istr,iend+1
671 cff1=0.5_r8*(1.0_r8+sign(1.0_r8,fx_r(i,j)))
672 cff2=0.5_r8*(1.0_r8-sign(1.0_r8,fx_r(i,j)))
673 fx(i,j)=0.5_r8*(1.0_r8+sign(1.0_r8,fx_r(i-1,j)))* &
674 & (cff1*fx_r(i-1,j)+ &
675 & cff2*0.5_r8*(fx_r(i-1,j)+fx_r(i,j)))+ &
676 & 0.5_r8*(1.0_r8-sign(1.0_r8,fx_r(i-1,j)))* &
677 & (cff2*fx_r(i ,j)+ &
678 & cff1*0.5_r8*(fx_r(i-1,j)+fx_r(i,j)))
679# ifdef MASKING
680 fx(i,j)=fx(i,j)*umask(i,j)
681# endif
682 END DO
683 END DO
684 DO j=jstr,jend+1
685 DO i=istr-1,iend+1
686 cff1=0.5_r8*(1.0_r8+sign(1.0_r8,fe_r(i,j)))
687 cff2=0.5_r8*(1.0_r8-sign(1.0_r8,fe_r(i,j)))
688 fe(i,j)=0.5_r8*(1.0_r8+sign(1.0_r8,fe_r(i,j-1)))* &
689 & (cff1*fe_r(i,j-1)+ &
690 & cff2*0.5_r8*(fe_r(i,j-1)+fe_r(i,j)))+ &
691 & 0.5_r8*(1.0_r8-sign(1.0_r8,fe_r(i,j-1)))* &
692 & (cff2*fe_r(i ,j)+ &
693 & cff1*0.5_r8*(fe_r(i,j-1)+fe_r(i,j)))
694# ifdef MASKING
695 fe(i,j)=fe(i,j)*vmask(i,j)
696# endif
697 END DO
698 END DO
699!
700! Limit fluxes to prevent bottom from breaking thru water surface.
701!
702! DO j=Jstr,Jend
703! DO i=Istr,Iend
704! cff1=1.0_r8/(Srho(ised,ng)*(1.0_r8-bed(i,j,1,iporo)))
705! rhs_bed=(FX(i+1,j)-FX(i,j)+ &
706! & FE(i,j+1)-FE(i,j))*pm(i,j)*pn(i,j)
707! cff2=MAX(rhs_bed*cff1+h(i,j)-Dcrit(ng),0.0_r8)
708! cff=cff2/ABS(cff2+eps)
709! FX(i ,j )=MAX(FX(i ,j ),0.0_r8)*cff+ &
710! & MIN(FX(i ,j ),0.0_r8)
711! FX(i+1,j )=MAX(FX(i+1,j ),0.0_r8)+ &
712! & MIN(FX(i+1,j ),0.0_r8)*cff
713! FE(i ,j )=MAX(FE(i ,j ),0.0_r8)*cff+ &
714! & MIN(FE(i ,j ),0.0_r8)
715! FE(i ,j+1)=MAX(FE(i ,j+1),0.0_r8)+ &
716! & MIN(FE(i ,j+1),0.0_r8)*cff
717! END DO
718! END DO
719!
720! Apply boundary conditions (gradient).
721!
722 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
723 IF (domain(ng)%Western_Edge(tile)) THEN
724 IF (lbc(iwest,istvar(idsed(ised)),ng)%closed) THEN
725 DO j=jstr-1,jend+1
726 fx(istr,j)=0.0_r8
727 END DO
728 END IF
729 END IF
730 END IF
731 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
732 IF (domain(ng)%Eastern_Edge(tile)) THEN
733 IF (lbc(ieast,istvar(idsed(ised)),ng)%closed) THEN
734 DO j=jstr-1,jend+1
735 fx(iend+1,j)=0.0_r8
736 END DO
737 END IF
738 END IF
739 END IF
740!
741 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
742 IF (domain(ng)%Southern_Edge(tile)) THEN
743 IF (lbc(isouth,istvar(idsed(ised)),ng)%closed) THEN
744 DO i=istr-1,iend+1
745 fe(i,jstr)=0.0_r8
746 END DO
747 END IF
748 END IF
749 END IF
750 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
751 IF (domain(ng)%Northern_Edge(tile)) THEN
752 IF (lbc(inorth,istvar(idsed(ised)),ng)%closed) THEN
753 DO i=istr-1,iend+1
754 fe(i,jend+1)=0.0_r8
755 END DO
756 END IF
757 END IF
758 END IF
759!
760! Determine flux divergence and evaluate change in bed properties.
761!
762 DO j=jstr,jend
763 DO i=istr,iend
764 cff=(fx(i+1,j)-fx(i,j)+ &
765 & fe(i,j+1)-fe(i,j))*pm(i,j)*pn(i,j)
766 bed_mass(i,j,1,nnew,ised)=max(bed_mass(i,j,1,nstp,ised)- &
767 & cff,0.0_r8)
768# if !defined SUSPLOAD
769 DO k=2,nbed
770 bed_mass(i,j,k,nnew,ised)=bed_mass(i,j,k,nstp,ised)
771 END DO
772# endif
773 bed(i,j,1,ithck)=max(bed(i,j,1,ithck)- &
774 & cff/(srho(ised,ng)* &
775 & (1.0_r8-bed(i,j,1,iporo))), &
776 & 0.0_r8)
777# ifdef MASKING
778 bed(i,j,1,ithck)=bed(i,j,1,ithck)*rmask(i,j)
779# endif
780 END DO
781 END DO
782!
783!-----------------------------------------------------------------------
784! Output bedload fluxes.
785!-----------------------------------------------------------------------
786!
787 cff=0.5_r8/dt(ng)
788 DO j=jstrr,jendr
789 DO i=istr,iendr
790 bedldu(i,j,ised)=fx(i,j)*(pn(i-1,j)+pn(i,j))*cff
791 END DO
792 END DO
793 DO j=jstr,jendr
794 DO i=istrr,iendr
795 bedldv(i,j,ised)=fe(i,j)*(pm(i,j-1)+pm(i,j))*cff
796 END DO
797 END DO
798 END DO
799!
800! Update mean surface properties.
801! Sd50 must be positive definite, due to BBL routines.
802! Srho must be >1000, due to (s-1) in BBL routines.
803!
804 DO j=jstr,jend
805 DO i=istr,iend
806 cff3=0.0_r8
807 DO ised=1,nst
808 cff3=cff3+bed_mass(i,j,1,nnew,ised)
809 END DO
810 IF (cff3.eq.0.0_r8) THEN
811 cff3=eps
812 END IF
813 DO ised=1,nst
814 bed_frac(i,j,1,ised)=bed_mass(i,j,1,nnew,ised)/cff3
815 END DO
816!
817 cff1=1.0_r8
818 cff2=1.0_r8
819 cff3=1.0_r8
820 cff4=1.0_r8
821 DO ised=1,nst
822 cff1=cff1*tau_ce(ised,ng)**bed_frac(i,j,1,ised)
823 cff2=cff2*sd50(ised,ng)**bed_frac(i,j,1,ised)
824 cff3=cff3*(wsed(ised,ng)+eps)**bed_frac(i,j,1,ised)
825 cff4=cff4*srho(ised,ng)**bed_frac(i,j,1,ised)
826 END DO
827 bottom(i,j,itauc)=cff1
828 bottom(i,j,isd50)=min(cff2,zob(ng))
829 bottom(i,j,iwsed)=cff3
830 bottom(i,j,idens)=max(cff4,1050.0_r8)
831 END DO
832 END DO
833# endif
834!
835!-----------------------------------------------------------------------
836! Apply periodic or gradient boundary conditions to property arrays.
837!-----------------------------------------------------------------------
838!
839 DO ised=1,nst
840 CALL bc_r3d_tile (ng, tile, &
841 & lbi, ubi, lbj, ubj, 1, nbed, &
842 & bed_frac(:,:,:,ised))
843 CALL bc_r3d_tile (ng, tile, &
844 & lbi, ubi, lbj, ubj, 1, nbed, &
845 & bed_mass(:,:,:,nnew,ised))
846# ifdef BEDLOAD
847 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
848 CALL exchange_u2d_tile (ng, tile, &
849 & lbi, ubi, lbj, ubj, &
850 & bedldu(:,:,ised))
851 CALL exchange_v2d_tile (ng, tile, &
852 & lbi, ubi, lbj, ubj, &
853 & bedldv(:,:,ised))
854 END IF
855# endif
856 END DO
857# ifdef DISTRIBUTE
858 CALL mp_exchange4d (ng, tile, inlm, 2, &
859 & lbi, ubi, lbj, ubj, 1, nbed, 1, nst, &
860 & nghostpoints, &
861 & ewperiodic(ng), nsperiodic(ng), &
862 & bed_frac, &
863 & bed_mass(:,:,:,nnew,:))
864# ifdef BEDLOAD
865 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
866 CALL mp_exchange3d (ng, tile, inlm, 2, &
867 & lbi, ubi, lbj, ubj, 1, nst, &
868 & nghostpoints, &
869 & ewperiodic(ng), nsperiodic(ng), &
870 & bedldu, bedldv)
871 END IF
872# endif
873# endif
874
875 DO i=1,mbedp
876 CALL bc_r3d_tile (ng, tile, &
877 & lbi, ubi, lbj, ubj, 1, nbed, &
878 & bed(:,:,:,i))
879 END DO
880# ifdef DISTRIBUTE
881 CALL mp_exchange4d (ng, tile, inlm, 1, &
882 & lbi, ubi, lbj, ubj, 1, nbed, 1, mbedp, &
883 & nghostpoints, &
884 & ewperiodic(ng), nsperiodic(ng), &
885 & bed)
886# endif
887
888 CALL bc_r3d_tile (ng, tile, &
889 & lbi, ubi, lbj, ubj, 1, mbotp, &
890 & bottom)
891# ifdef DISTRIBUTE
892 CALL mp_exchange3d (ng, tile, inlm, 1, &
893 & lbi, ubi, lbj, ubj, 1, mbotp, &
894 & nghostpoints, &
895 & ewperiodic(ng), nsperiodic(ng), &
896 & bottom)
897# endif
898!
899 RETURN
subroutine bc_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
Definition bc_3d.F:48
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
integer, dimension(:), allocatable istvar
integer nbed
Definition mod_param.F:517
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer ncs
Definition mod_param.F:525
integer nghostpoints
Definition mod_param.F:710
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer nst
Definition mod_param.F:521
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:,:), allocatable compositegrid
integer, parameter isouth
real(r8), dimension(:), allocatable zob
integer, parameter ieast
real(dp) g
real(dp) rho0
integer, parameter inorth
real(dp), parameter pi
integer, parameter iwsed
real(r8), dimension(:), allocatable bedload_coeff
integer, parameter mbotp
real(r8), dimension(:,:), allocatable srho
integer, dimension(:), allocatable idsed
integer, parameter iporo
integer, parameter idens
integer, parameter isd50
real(r8), dimension(:,:), allocatable morph_fac
real(r8), dimension(:,:), allocatable sd50
real(r8), dimension(:,:), allocatable wsed
integer, parameter itauc
integer, parameter ithck
real(r8), dimension(:,:), allocatable tau_ce
integer, parameter mbedp
subroutine mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, a, b, c)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)

References bc_3d_mod::bc_r3d_tile(), mod_sediment::bedload_coeff, mod_scalars::compositegrid, mod_param::domain, mod_scalars::dt, mod_scalars::ewperiodic, exchange_2d_mod::exchange_u2d_tile(), exchange_2d_mod::exchange_v2d_tile(), mod_scalars::g, mod_sediment::idens, mod_sediment::idsed, mod_scalars::ieast, mod_param::inlm, mod_scalars::inorth, mod_sediment::iporo, mod_sediment::isd50, mod_scalars::isouth, mod_ncparam::istvar, mod_sediment::itauc, mod_sediment::ithck, mod_scalars::iwest, mod_sediment::iwsed, mod_param::lbc, mod_sediment::morph_fac, mp_exchange_mod::mp_exchange3d(), mp_exchange_mod::mp_exchange4d(), mod_param::ncs, mod_param::nghostpoints, mod_scalars::nsperiodic, mod_scalars::pi, mod_scalars::rho0, mod_sediment::sd50, mod_sediment::srho, mod_sediment::tau_ce, mod_sediment::wsed, and mod_scalars::zob.

Referenced by sed_bedload().

Here is the call graph for this function:
Here is the caller graph for this function: