-
Notifications
You must be signed in to change notification settings - Fork 0
/
SwanTranspAc.ftn90
178 lines (178 loc) · 8.15 KB
/
SwanTranspAc.ftn90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
subroutine SwanTranspAc ( amat , rhs , leakcf, ac2 , ac1 , &
cgo , cax , cay , cad , cas , &
anybin, rdx , rdy , spcsig, spcdir, &
obredf, idcmin, idcmax, iscmin, iscmax, &
iddlow, iddtop, isslow, isstop, anyblk, &
trac0 , trac1 )
!
! --|-----------------------------------------------------------|--
! | Delft University of Technology |
! | Faculty of Civil Engineering and Geosciences |
! | Environmental Fluid Mechanics Section |
! | P.O. Box 5048, 2600 GA Delft, The Netherlands |
! | |
! | Programmer: Marcel Zijlema |
! --|-----------------------------------------------------------|--
!
!
! SWAN (Simulating WAves Nearshore); a third generation wave model
! Copyright (C) 1993-2024 Delft University of Technology
!
! This program is free software: you can redistribute it and/or modify
! it under the terms of the GNU General Public License as published by
! the Free Software Foundation, either version 3 of the License, or
! (at your option) any later version.
!
! This program is distributed in the hope that it will be useful,
! but WITHOUT ANY WARRANTY; without even the implied warranty of
! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
! GNU General Public License for more details.
!
! You should have received a copy of the GNU General Public License
! along with this program. If not, see <http://www.gnu.org/licenses/>.
!
!
! Authors
!
! 40.80: Marcel Zijlema
!
! Updates
!
! 40.80, August 2007: New subroutine
! 40.85, August 2008: add tranport terms for output purposes
! 41.00, February 2009: add GSE correction
! 41.07, July 2009: add explicit scheme for sigma space
!
! Purpose
!
! Computes the transport part of the action balance equation
!
! Modules used
!
use ocpcomm4
use swcomm2
use swcomm3
use swcomm4
use SwanGriddata
!
implicit none
!
! Argument variables
!
integer, intent(in) :: iddlow ! minimum direction bin that is propagated within a sweep
integer, intent(in) :: iddtop ! maximum direction bin that is propagated within a sweep
integer, intent(in) :: isslow ! minimum frequency that is propagated within a sweep
integer, intent(in) :: isstop ! maximum frequency that is propagated within a sweep
!
integer, dimension(MSC), intent(in) :: idcmax ! maximum frequency-dependent counter in directional space
integer, dimension(MSC), intent(in) :: idcmin ! minimum frequency-dependent counter in directional space
integer, dimension(MDC), intent(in) :: iscmax ! maximum direction-dependent counter in frequency space
integer, dimension(MDC), intent(in) :: iscmin ! minimum direction-dependent counter in frequency space
!
real, dimension(MDC,MSC,nverts), intent(in) :: ac1 ! action density at previous time level
real, dimension(MDC,MSC,nverts), intent(in) :: ac2 ! action density at current time level
real, dimension(MDC,MSC,5), intent(out) :: amat ! coefficient matrix of system of equations in (sigma,theta) space
! 1: correspond to point (l ,m )
! 2: correspond to point (l-1,m )
! 3: correspond to point (l+1,m )
! 4: correspond to point (l ,m-1)
! 5: correspond to point (l ,m+1)
real, dimension(MDC,MSC), intent(in) :: cad ! wave transport velocity in theta-direction
real, dimension(MDC,MSC), intent(in) :: cas ! wave transport velocity in sigma-direction
real, dimension(MDC,MSC,ICMAX), intent(in) :: cax ! wave transport velocity in x-direction
real, dimension(MDC,MSC,ICMAX), intent(in) :: cay ! wave transport velocity in y-direction
real, dimension(MSC,ICMAX), intent(in) :: cgo ! group velocity
real, dimension(MDC,MSC), intent(out) :: leakcf ! leak coefficient
real, dimension(MDC,MSC,2), intent(in) :: obredf ! action reduction coefficient based on transmission
real, dimension(2), intent(in) :: rdx ! first component of contravariant base vector rdx(b) = a^(b)_1
real, dimension(2), intent(in) :: rdy ! second component of contravariant base vector rdy(b) = a^(b)_2
real, dimension(MDC,MSC), intent(out) :: rhs ! right-hand side of system of equations in (sigma,theta) space
real, dimension(MDC,6), intent(in) :: spcdir ! (*,1): spectral direction bins (radians)
! (*,2): cosine of spectral directions
! (*,3): sine of spectral directions
! (*,4): cosine^2 of spectral directions
! (*,5): cosine*sine of spectral directions
! (*,6): sine^2 of spectral directions
real, dimension(MSC), intent(in) :: spcsig ! relative frequency bins
real, dimension(MDC,MSC,MTRNP), intent(out) :: trac0 ! explicit part of propagation in present vertex for output purposes
real, dimension(MDC,MSC,MTRNP), intent(out) :: trac1 ! implicit part of propagation in present vertex for output purposes
!
logical, dimension(MDC,MSC), intent(in) :: anybin ! true if bin is active in considered sweep
logical, dimension(MDC,MSC), intent(out) :: anyblk ! true if bin is blocked by a counter current based on a CFL criterion
!
! Local variables
!
integer, save :: ient = 0 ! number of entries in this subroutine
!
! Structure
!
! Description of the pseudo code
!
! Source text
!
if (ltrace) call strace (ient,'SwanTranspAc')
!
! set matrix and right-hand side to zero
!
amat = 0.
rhs = 0.
!
! set leak and transport coefficients to zero
!
leakcf = 0.
trac0 = 0.
trac1 = 0.
!
! compute transport in x-y space
!
!TIMG call SWTSTA(140)
call SwanTranspX ( amat , rhs , ac2 , ac1 , cax , cay , &
rdx , rdy , obredf, idcmin, idcmax, isslow, &
isstop , trac0, trac1 )
!
! add GSE correction, if appropriate
!
if ( WAVAGE > 0. ) call SwanGSECorr ( rhs, ac2, cgo, spcdir, idcmin, idcmax, isslow, isstop, trac0 )
!TIMG call SWTSTO(140)
!
! compute transport in theta space
!
!TIMG call SWTSTA(142)
if ( IREFR /= 0 ) then
!
call STRSD ( DDIR , idcmin , idcmax , cad , &
amat(1,1,4), amat(1,1,1), amat(1,1,5), rhs , &
ac2 , isstop , anybin , leakcf , &
trac0 , trac1 )
!
endif
!TIMG call SWTSTO(142)
!
! compute transport in sigma space
!
!TIMG call SWTSTA(141)
if ( (DYNDEP .OR. ICUR /= 0) .and. ITFRE /= 0 ) then
!
if ( int(PNUMS(8)) == 1 ) then
!
! implicit scheme
!
call STRSSI ( spcsig , cas , amat(1,1,2), amat(1,1,1), &
amat(1,1,3), anybin, rhs , ac2 , &
iscmin , iscmax, iddlow , iddtop , &
trac0 , trac1 )
!
elseif ( int(PNUMS(8)) == 2 ) then
!
! explicit scheme
!
call STRSSB ( iddlow, iddtop, idcmin, idcmax, isstop, &
cax , cay , cas , ac2 , spcsig, &
rhs , anyblk, rdx , rdy , trac0 )
!
endif
!
endif
!TIMG call SWTSTO(141)
!
end subroutine SwanTranspAc