-
Notifications
You must be signed in to change notification settings - Fork 26
/
Copy pathmod_tracers.F90
285 lines (248 loc) · 8.99 KB
/
mod_tracers.F90
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
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
! ------------------------------------------------------------------------------
! Copyright (C) 2007-2024 Mats Bentsen, Jörg Schwinger, Jerry Tjiputra,
! Alok Kumar Gupta, Mariana Vertenstein
!
! This file is part of BLOM.
!
! BLOM is free software: you can redistribute it and/or modify it under the
! terms of the GNU Lesser General Public License as published by the Free
! Software Foundation, either version 3 of the License, or (at your option)
! any later version.
!
! BLOM 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 Lesser General Public License for
! more details.
!
! You should have received a copy of the GNU Lesser General Public License
! along with BLOM. If not, see <https://www.gnu.org/licenses/>.
! ------------------------------------------------------------------------------
module mod_tracers
! ----------------------------------------------------------------------
! This module declares variables related to tracers.
! ----------------------------------------------------------------------
use mod_types, only: r8
use mod_constants, only: spval
#ifdef HAMOCC
use mo_param1_bgc, only: init_indices, nocetra
#endif
use mod_ifdefs, only: use_IDLAGE, use_TKE, use_TRC
use mod_xc
implicit none
private
! Number of ocean tracers (ntrocn) does not include turbulent
! kinetic energy, generic length scale, ideal age and HAMOCC tracers.
integer, protected :: ntrocn ! Number of ocean tracers
integer, protected :: natr ! Number of age tracers.
integer, protected :: ntriag ! Number of ideal age tracer.
integer, protected :: ntrtke ! Number of turbulent kinetic energy tracers
integer, protected :: ntrgls ! Number of generic length scale tracers.
integer, protected :: ntrbgc = -999 ! Number of HAMOCC tracers.
integer, protected :: ntr = -999 ! Total number of tracers. (Return error if default value is used)
integer, protected :: itrtke
integer, protected :: itriag ! Index of first ideal age tracer.
integer, protected :: itrgls
integer, protected :: itrbgc = -999 ! Index of first HAMOCC tracer.
real(r8), allocatable, dimension(:,:,:,:) :: &
trc, & ! Tracer array.
trcold ! Tracer array at old time level.
real(r8), allocatable, dimension(:,:,:) :: &
uflxtr, & ! u-component of tracer flux.
vflxtr, & ! v-component of tracer flux.
trflx ! Surface flux of tracer.
public :: ntrocn, ntrtke, ntrgls, ntriag, ntrbgc, ntr, natr
public :: itrtke, itrgls, itriag, itrbgc
public :: trc, trcold, uflxtr, vflxtr, trflx
! Public routines
public :: inivar_tracers
! Private routines
private :: allocate_tracers
contains
subroutine inivar_tracers()
! ---------------------------------------------------------------------------
! Initialize variables related to tracers.
! ---------------------------------------------------------------------------
! Local variables
real(r8), dimension(1 - nbdy:idm + nbdy, 1 - nbdy:jdm + nbdy) :: uflux
real(r8), dimension(1 - nbdy:idm + nbdy, 1 - nbdy:jdm + nbdy) :: vflux
integer :: i, j, k, l, nt
! Set tracer numbers
ntrocn = 0
natr = 0 ! Number of age tracers.
if (use_TKE) then
ntrtke = 1
ntrgls = 1
itrtke = ntrocn - natr + 1 ! Index of first turbulent kinetic energy
itrgls = ntrocn - natr + ntrtke + 1 ! generic length scale tracers.
else
ntrtke = 0
ntrgls = 0
itrtke = -1
itrgls = -1
end if
if (use_IDLAGE) then
ntriag = 1 ! Number of ideal age tracer.
itriag = ntrocn - natr + ntrtke + ntrgls + 1 ! Index of first ideal age tracer.
else
ntriag = 0
itriag = -1
end if
! Set first indices
if (use_TKE) then
itrtke = ntrocn - natr + 1 ! Index of first turbulent kinetic energy
itrgls = ntrocn - natr + ntrtke + 1 ! generic length scale tracers.
else
itrtke = -1
itrgls = -1
end if
if (use_IDLAGE) then
itriag = ntrocn - natr + ntrtke + ntrgls + 1 ! Index of first ideal age tracer.
else
itriag = -1
end if
#ifdef HAMOCC
call init_indices()
ntrbgc = nocetra
itrbgc = ntrocn - natr + ntrtke + ntrgls + ntriag + 1
#else
ntrbgc = 0
itrbgc = -1
#endif
! Total number of tracers.
if (mnproc.eq.1) then
ntr = ntrocn + ntrtke + ntrgls + ntriag + ntrbgc
write(lp,'(A,1X,I8)') 'Total tracer count: ', ntr
if (ntr < 0) then
write(lp,'(a)') 'Number of tracers must be non-negative.'
call xchalt('(allocate_tracers)')
stop '(allocate_tracers)'
endif
endif
call xcbcst(ntr)
if (use_TRC) then
! Allocate tracer arrays
call allocate_tracers
!$omp parallel do private(i, k, nt)
do j = 1 - nbdy, jj + nbdy
do i = 1 - nbdy, ii + nbdy
uflux(i, j) = spval
vflux(i, j) = spval
do nt = 1, ntr
trflx(nt, i, j) = spval
enddo
do k = 1, kk
do nt = 1, ntr
trc(i ,j ,k , nt) = spval
trc(i ,j ,k + kk, nt) = spval
enddo
enddo
enddo
enddo
!$omp end parallel do
! Initialize uflxtr at points located upstream and downstream (in i
! direction) of p-points.
!$omp parallel do private(l, i, k)
do j = 1, jj
do l = 1, isp(j)
do i = max(1, ifp(j, l)), min(ii, ilp(j, l) + 1)
uflux(i, j) = 0._r8
enddo
enddo
enddo
!$omp end parallel do
call xctilr(uflux, 1, 1, nbdy, nbdy, halo_us)
!$omp parallel do private(i, nt)
do j = 1 - nbdy, jdm + nbdy
do i = 1 - nbdy, idm + nbdy
do nt = 1, ntr
uflxtr(nt, i, j) = uflux(i, j)
enddo
enddo
enddo
!$omp end parallel do
! Initialize vflxtr at points located upstream and downstream (in j
! direction) of p-points.
!$omp parallel do private(l, j, k)
do i = 1, ii
do l = 1, jsp(i)
do j = max(1, jfp(i, l)), min(jj, jlp(i, l) + 1)
vflux(i, j)=0._r8
enddo
enddo
enddo
!$omp end parallel do
call xctilr(vflux, 1, 1, nbdy, nbdy, halo_vs)
!$omp parallel do private(i, nt)
do j = 1 - nbdy, jdm + nbdy
do i = 1 - nbdy, idm + nbdy
do nt = 1, ntr
vflxtr(nt, i, j) = vflux(i, j)
enddo
enddo
enddo
!$omp end parallel do
! Initialize surface flux of tracer.
!$omp parallel do private(l, i, nt)
do j = 1, jj
do l = 1, isp(j)
do i = max(1, ifp(j, l)), min(ii, ilp(j, l))
do nt = 1, ntr
trflx(nt, i, j) = 0._r8
enddo
enddo
enddo
enddo
!$omp end parallel do
end if ! if (use_TRC)
end subroutine inivar_tracers
subroutine allocate_tracers
! ---------------------------------------------------------------------------
! If using HAMOCC, set number of bgc tracers : ntrbgc, itrbgc
! Set total number of tracers : ntr
! Allocate tracer arrays.
! ---------------------------------------------------------------------------
integer :: errstat
integer :: num_bgc_tracers = 0
! Number of HAMOCC tracers.
! Tracer array.
allocate (trc(1 - nbdy:idm + nbdy, 1 - nbdy:jdm + nbdy, 2*kdm, ntr), &
stat = errstat)
if(errstat.ne.0) then
write(lp,'(a)') 'Not enough memory for variable trc.'
call xchalt('(allocate_tracers)')
stop '(allocate_tracers)'
endif
! Tracer array at old time level.
allocate (trcold(1 - nbdy:idm + nbdy, 1 - nbdy:jdm + nbdy, kdm, ntr), &
stat = errstat)
if(errstat.ne.0) then
write(lp,'(a)') 'Not enough memory for variable trcold.'
call xchalt('(allocate_tracers)')
stop '(allocate_tracers)'
endif
! u-component of tracer flux.
allocate (uflxtr(ntr, 1 - nbdy:idm + nbdy, 1 - nbdy:jdm + nbdy), &
stat = errstat)
if(errstat.ne.0) then
write(lp,'(a)') 'Not enough memory for variable uflxtr.'
call xchalt('(allocate_tracers)')
stop '(allocate_tracers)'
endif
! v-component of tracer flux.
allocate (vflxtr(ntr, 1 - nbdy:idm + nbdy, 1 - nbdy:jdm + nbdy), &
stat = errstat)
if(errstat.ne.0) then
write(lp,'(a)') 'Not enough memory for variable vflxtr.'
call xchalt('(allocate_tracers)')
stop '(allocate_tracers)'
endif
! Surface flux of tracer.
allocate (trflx(ntr, 1 - nbdy:idm + nbdy, 1 - nbdy:jdm + nbdy), &
stat = errstat)
if(errstat.ne.0) then
write(lp,'(a)') 'allocate_tracers: Not enough memory for variable trflx.'
call xchalt('(allocate_tracers)')
stop '(allocate_tracers)'
endif
end subroutine allocate_tracers
end module mod_tracers