ROMS
Loading...
Searching...
No Matches
sed_surface.F
Go to the documentation of this file.
1#include "cppdefs.h"
2
4
5#if defined NONLINEAR && defined SEDIMENT
6!
7!git $Id$
8!==================================================== John C. Warner ===
9! Copyright (c) 2002-2025 The ROMS Group Hernan G. Arango !
10! Licensed under a MIT/X style license !
11! See License_ROMS.md !
12!=======================================================================
13! !
14! This routine computed sediment surface layer (sediment-water !
15! interface) properties. !
16! !
17! References: !
18! !
19! Warner, J.C., C.R. Sherwood, R.P. Signell, C.K. Harris, and H.G. !
20! Arango, 2008: Development of a three-dimensional, regional, !
21! coupled wave, current, and sediment-transport model, Computers !
22! & Geosciences, 34, 1284-1306. !
23! !
24!=======================================================================
25!
26 implicit none
27
28 PRIVATE
29 PUBLIC :: sed_surface
30
31 CONTAINS
32!
33!***********************************************************************
34 SUBROUTINE sed_surface (ng, tile)
35!***********************************************************************
36!
37 USE mod_param
38 USE mod_ocean
39 USE mod_sedbed
40 USE mod_stepping
41!
42! Imported variable declarations.
43!
44 integer, intent(in) :: ng, tile
45!
46! Local variable declarations.
47!
48 character (len=*), parameter :: myfile = &
49 & __FILE__
50!
51# include "tile.h"
52!
53# ifdef PROFILE
54 CALL wclock_on (ng, inlm, 16, __line__, myfile)
55# endif
56 CALL sed_surface_tile (ng, tile, &
57 & lbi, ubi, lbj, ubj, &
58 & imins, imaxs, jmins, jmaxs, &
59 & nstp(ng), nnew(ng), &
60 & sedbed(ng) % bed_frac, &
61 & sedbed(ng) % bottom)
62# ifdef PROFILE
63 CALL wclock_off (ng, inlm, 16, __line__, myfile)
64# endif
65!
66 RETURN
67 END SUBROUTINE sed_surface
68!
69!***********************************************************************
70 SUBROUTINE sed_surface_tile (ng, tile, &
71 & LBi, UBi, LBj, UBj, &
72 & IminS, ImaxS, JminS, JmaxS, &
73 & nstp, nnew, &
74 & bed_frac, bottom)
75!***********************************************************************
76!
77 USE mod_param
78 USE mod_scalars
79 USE mod_sediment
80!
81 USE bc_3d_mod, ONLY : bc_r3d_tile
82# ifdef DISTRIBUTE
84# endif
85!
86! Imported variable declarations.
87!
88 integer, intent(in) :: ng, tile
89 integer, intent(in) :: LBi, UBi, LBj, UBj
90 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
91 integer, intent(in) :: nstp, nnew
92!
93# ifdef ASSUMED_SHAPE
94 real(r8), intent(inout) :: bed_frac(LBi:,LBj:,:,:)
95 real(r8), intent(inout) :: bottom(LBi:,LBj:,:)
96# else
97 real(r8), intent(inout) :: bed_frac(LBi:UBi,LBj:UBj,Nbed,NST)
98 real(r8), intent(inout) :: bottom(LBi:UBi,LBj:UBj,MBOTP)
99# endif
100!
101! Local variable declarations.
102!
103 integer :: i, ised, j
104
105 real(r8) :: cff1, cff2, cff3, cff4
106 real(r8), parameter :: eps = 1.0e-14_r8
107
108# include "set_bounds.h"
109
110!
111! Update mean surface properties.
112! Sd50 must be positive definite, due to BBL routines.
113! Srho must be >1000, due to (s-1) in BBL routines
114!
115 j_loop : DO j=jstr,jend
116 DO i=istr,iend
117 cff1=1.0_r8
118 cff2=1.0_r8
119 cff3=1.0_r8
120 cff4=1.0_r8
121 DO ised=1,nst
122 cff1=cff1*tau_ce(ised,ng)**bed_frac(i,j,1,ised)
123 cff2=cff2*sd50(ised,ng)**bed_frac(i,j,1,ised)
124 cff3=cff3*(wsed(ised,ng)+eps)**bed_frac(i,j,1,ised)
125 cff4=cff4*srho(ised,ng)**bed_frac(i,j,1,ised)
126 END DO
127 bottom(i,j,itauc)=cff1
128 bottom(i,j,isd50)=min(cff2,zob(ng))
129 bottom(i,j,iwsed)=cff3
130 bottom(i,j,idens)=max(cff4,1050.0_r8)
131 END DO
132 END DO j_loop
133!!
134!! Move bottom roughness computations from BBL routines to here.
135!!
136
137!
138!-----------------------------------------------------------------------
139! Apply periodic or gradient boundary conditions to property arrays.
140!-----------------------------------------------------------------------
141!
142 CALL bc_r3d_tile (ng, tile, &
143 & lbi, ubi, lbj, ubj, 1, mbotp, &
144 & bottom)
145# ifdef DISTRIBUTE
146 CALL mp_exchange3d (ng, tile, inlm, 1, &
147 & lbi, ubi, lbj, ubj, 1, mbotp, &
148 & nghostpoints, &
149 & ewperiodic(ng), nsperiodic(ng), &
150 & bottom)
151# endif
152!
153 RETURN
154 END SUBROUTINE sed_surface_tile
155#endif
156 END MODULE sed_surface_mod
subroutine bc_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
Definition bc_3d.F:48
integer, parameter inlm
Definition mod_param.F:662
integer nghostpoints
Definition mod_param.F:710
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
real(r8), dimension(:), allocatable zob
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
integer, parameter iwsed
real(r8), dimension(:,:), allocatable srho
integer, parameter idens
integer, parameter isd50
real(r8), dimension(:,:), allocatable sd50
real(r8), dimension(:,:), allocatable wsed
integer, parameter itauc
real(r8), dimension(:,:), allocatable tau_ce
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
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)
subroutine sed_surface_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, bed_frac, bottom)
Definition sed_surface.F:75
subroutine, public sed_surface(ng, tile)
Definition sed_surface.F:35
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