-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathdiaar5.f90
305 lines (305 loc) · 10.6 KB
/
diaar5.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
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
MODULE diaar5
USE oce
USE dom_oce
USE eosbn2
USE phycst
USE in_out_manager
USE zdfddm
USE zdf_oce
USE lib_mpp
USE iom
USE fldread
USE timing
IMPLICIT NONE
PRIVATE
PUBLIC :: dia_ar5
PUBLIC :: dia_ar5_alloc
PUBLIC :: dia_ar5_hst
REAL(KIND = wp) :: vol0
REAL(KIND = wp) :: area_tot
REAL(KIND = wp), ALLOCATABLE, SAVE, DIMENSION(:, :) :: area
REAL(KIND = wp), ALLOCATABLE, SAVE, DIMENSION(:, :) :: thick0
REAL(KIND = wp), ALLOCATABLE, SAVE, DIMENSION(:, :, :) :: sn0
LOGICAL :: l_ar5
CONTAINS
FUNCTION dia_ar5_alloc()
INTEGER :: dia_ar5_alloc
ALLOCATE(area(jpi, jpj), thick0(jpi, jpj), sn0(jpi, jpj, jpk), STAT = dia_ar5_alloc)
IF (lk_mpp) CALL mpp_sum(dia_ar5_alloc)
IF (dia_ar5_alloc /= 0) CALL ctl_warn('dia_ar5_alloc: failed to allocate arrays')
END FUNCTION dia_ar5_alloc
SUBROUTINE dia_ar5(kt)
INTEGER, INTENT( IN ) :: kt
INTEGER :: ji, jj, jk
REAL(KIND = wp) :: zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass
REAL(KIND = wp) :: zaw, zbw, zrw
REAL(KIND = wp), ALLOCATABLE, DIMENSION(:, :) :: zarea_ssh, zbotpres
REAL(KIND = wp), ALLOCATABLE, DIMENSION(:, :) :: zpe
REAL(KIND = wp), ALLOCATABLE, DIMENSION(:, :, :) :: zrhd, zrhop
REAL(KIND = wp), ALLOCATABLE, DIMENSION(:, :, :, :) :: ztsn
IF (ln_timing) CALL timing_start('dia_ar5')
IF (kt == nit000) CALL dia_ar5_init
IF (l_ar5) THEN
ALLOCATE(zarea_ssh(jpi, jpj), zbotpres(jpi, jpj))
ALLOCATE(zrhd(jpi, jpj, jpk), zrhop(jpi, jpj, jpk))
ALLOCATE(ztsn(jpi, jpj, jpk, jpts))
!$ACC KERNELS
zarea_ssh(:, :) = area(:, :) * sshn(:, :)
!$ACC END KERNELS
END IF
IF (iom_use('voltot') .OR. iom_use('sshtot') .OR. iom_use('sshdyn')) THEN
zvolssh = SUM(zarea_ssh(:, :))
IF (lk_mpp) CALL mpp_sum(zvolssh)
zvol = vol0 + zvolssh
CALL iom_put('voltot', zvol)
CALL iom_put('sshtot', zvolssh / area_tot)
CALL iom_put('sshdyn', sshn(:, :) - (zvolssh / area_tot))
END IF
IF (iom_use('botpres') .OR. iom_use('sshthster') .OR. iom_use('sshsteric')) THEN
!$ACC KERNELS
ztsn(:, :, :, jp_tem) = tsn(:, :, :, jp_tem)
ztsn(:, :, :, jp_sal) = sn0(:, :, :)
!$ACC END KERNELS
CALL eos(ztsn, zrhd, gdept_n(:, :, :))
!$ACC KERNELS
zbotpres(:, :) = 0._wp
DO jk = 1, jpkm1
zbotpres(:, :) = zbotpres(:, :) + e3t_n(:, :, jk) * zrhd(:, :, jk)
END DO
!$ACC END KERNELS
IF (ln_linssh) THEN
IF (ln_isfcav) THEN
!$ACC KERNELS
DO ji = 1, jpi
DO jj = 1, jpj
zbotpres(ji, jj) = zbotpres(ji, jj) + sshn(ji, jj) * zrhd(ji, jj, mikt(ji, jj)) + riceload(ji, jj)
END DO
END DO
!$ACC END KERNELS
ELSE
!$ACC KERNELS
zbotpres(:, :) = zbotpres(:, :) + sshn(:, :) * zrhd(:, :, 1)
!$ACC END KERNELS
END IF
END IF
zarho = SUM(area(:, :) * zbotpres(:, :))
IF (lk_mpp) CALL mpp_sum(zarho)
zssh_steric = - zarho / area_tot
CALL iom_put('sshthster', zssh_steric)
CALL eos(tsn, zrhd, zrhop, gdept_n(:, :, :))
!$ACC KERNELS
zrhop(:, :, jpk) = 0._wp
!$ACC END KERNELS
CALL iom_put('rhop', zrhop)
!$ACC KERNELS
zbotpres(:, :) = 0._wp
DO jk = 1, jpkm1
zbotpres(:, :) = zbotpres(:, :) + e3t_n(:, :, jk) * zrhd(:, :, jk)
END DO
!$ACC END KERNELS
IF (ln_linssh) THEN
IF (ln_isfcav) THEN
!$ACC KERNELS
DO ji = 1, jpi
DO jj = 1, jpj
zbotpres(ji, jj) = zbotpres(ji, jj) + sshn(ji, jj) * zrhd(ji, jj, mikt(ji, jj)) + riceload(ji, jj)
END DO
END DO
!$ACC END KERNELS
ELSE
!$ACC KERNELS
zbotpres(:, :) = zbotpres(:, :) + sshn(:, :) * zrhd(:, :, 1)
!$ACC END KERNELS
END IF
END IF
zarho = SUM(area(:, :) * zbotpres(:, :))
IF (lk_mpp) CALL mpp_sum(zarho)
zssh_steric = - zarho / area_tot
CALL iom_put('sshsteric', zssh_steric)
!$ACC KERNELS
zztmp = rau0 * grav * 1.E-4_wp
zbotpres(:, :) = zztmp * (zbotpres(:, :) + sshn(:, :) + thick0(:, :))
!$ACC END KERNELS
CALL iom_put('botpres', zbotpres)
END IF
IF (iom_use('masstot') .OR. iom_use('temptot') .OR. iom_use('saltot')) THEN
!$ACC KERNELS
ztemp = 0._wp
zsal = 0._wp
DO jk = 1, jpkm1
DO jj = 1, jpj
DO ji = 1, jpi
zztmp = area(ji, jj) * e3t_n(ji, jj, jk)
ztemp = ztemp + zztmp * tsn(ji, jj, jk, jp_tem)
zsal = zsal + zztmp * tsn(ji, jj, jk, jp_sal)
END DO
END DO
END DO
!$ACC END KERNELS
IF (ln_linssh) THEN
IF (ln_isfcav) THEN
!$ACC KERNELS
DO ji = 1, jpi
DO jj = 1, jpj
ztemp = ztemp + zarea_ssh(ji, jj) * tsn(ji, jj, mikt(ji, jj), jp_tem)
zsal = zsal + zarea_ssh(ji, jj) * tsn(ji, jj, mikt(ji, jj), jp_sal)
END DO
END DO
!$ACC END KERNELS
ELSE
ztemp = ztemp + SUM(zarea_ssh(:, :) * tsn(:, :, 1, jp_tem))
zsal = zsal + SUM(zarea_ssh(:, :) * tsn(:, :, 1, jp_sal))
END IF
END IF
IF (lk_mpp) THEN
CALL mpp_sum(ztemp)
CALL mpp_sum(zsal)
END IF
zmass = rau0 * (zarho + zvol)
ztemp = ztemp / zvol
zsal = zsal / zvol
CALL iom_put('masstot', zmass)
CALL iom_put('temptot', ztemp)
CALL iom_put('saltot', zsal)
END IF
IF (iom_use('tnpeo')) THEN
ALLOCATE(zpe(jpi, jpj))
!$ACC KERNELS
zpe(:, :) = 0._wp
!$ACC END KERNELS
IF (ln_zdfddm) THEN
!$ACC KERNELS
DO jk = 2, jpk
DO jj = 1, jpj
DO ji = 1, jpi
IF (rn2(ji, jj, jk) > 0._wp) THEN
zrw = (gdepw_n(ji, jj, jk) - gdept_n(ji, jj, jk)) / (gdept_n(ji, jj, jk - 1) - gdept_n(ji, jj, jk))
zaw = rab_n(ji, jj, jk, jp_tem) * (1. - zrw) + rab_n(ji, jj, jk - 1, jp_tem) * zrw
zbw = rab_n(ji, jj, jk, jp_sal) * (1. - zrw) + rab_n(ji, jj, jk - 1, jp_sal) * zrw
zpe(ji, jj) = zpe(ji, jj) - grav * (avt(ji, jj, jk) * zaw * (tsn(ji, jj, jk - 1, jp_tem) - tsn(ji, jj, jk, jp_tem)) - avs(ji, jj, jk) * zbw * (tsn(ji, jj, jk - 1, jp_sal) - tsn(ji, jj, jk, jp_sal)))
END IF
END DO
END DO
END DO
!$ACC END KERNELS
ELSE
!$ACC KERNELS
DO jk = 1, jpk
DO ji = 1, jpi
DO jj = 1, jpj
zpe(ji, jj) = zpe(ji, jj) + avt(ji, jj, jk) * MIN(0._wp, rn2(ji, jj, jk)) * rau0 * e3w_n(ji, jj, jk)
END DO
END DO
END DO
!$ACC END KERNELS
END IF
CALL iom_put('tnpeo', zpe)
DEALLOCATE(zpe)
END IF
IF (l_ar5) THEN
DEALLOCATE(zarea_ssh, zbotpres)
DEALLOCATE(zrhd, zrhop)
DEALLOCATE(ztsn)
END IF
IF (ln_timing) CALL timing_stop('dia_ar5')
END SUBROUTINE dia_ar5
SUBROUTINE dia_ar5_hst(ktra, cptr, pua, pva)
INTEGER, INTENT(IN ) :: ktra
CHARACTER(LEN = 3), INTENT(IN) :: cptr
REAL(KIND = wp), DIMENSION(jpi, jpj, jpk), INTENT(IN) :: pua
REAL(KIND = wp), DIMENSION(jpi, jpj, jpk), INTENT(IN) :: pva
INTEGER :: ji, jj, jk
REAL(KIND = wp), DIMENSION(jpi, jpj) :: z2d
!$ACC KERNELS
z2d(:, :) = pua(:, :, 1)
DO jk = 1, jpkm1
DO jj = 2, jpjm1
DO ji = 2, jpim1
z2d(ji, jj) = z2d(ji, jj) + pua(ji, jj, jk)
END DO
END DO
END DO
!$ACC END KERNELS
CALL lbc_lnk(z2d, 'U', - 1.)
IF (cptr == 'adv') THEN
IF (ktra == jp_tem) CALL iom_put("uadv_heattr", rau0_rcp * z2d)
IF (ktra == jp_sal) CALL iom_put("uadv_salttr", rau0 * z2d)
END IF
IF (cptr == 'ldf') THEN
IF (ktra == jp_tem) CALL iom_put("udiff_heattr", rau0_rcp * z2d)
IF (ktra == jp_sal) CALL iom_put("udiff_salttr", rau0 * z2d)
END IF
!$ACC KERNELS
z2d(:, :) = pva(:, :, 1)
DO jk = 1, jpkm1
DO jj = 2, jpjm1
DO ji = 2, jpim1
z2d(ji, jj) = z2d(ji, jj) + pva(ji, jj, jk)
END DO
END DO
END DO
!$ACC END KERNELS
CALL lbc_lnk(z2d, 'V', - 1.)
IF (cptr == 'adv') THEN
IF (ktra == jp_tem) CALL iom_put("vadv_heattr", rau0_rcp * z2d)
IF (ktra == jp_sal) CALL iom_put("vadv_salttr", rau0 * z2d)
END IF
IF (cptr == 'ldf') THEN
IF (ktra == jp_tem) CALL iom_put("vdiff_heattr", rau0_rcp * z2d)
IF (ktra == jp_sal) CALL iom_put("vdiff_salttr", rau0 * z2d)
END IF
END SUBROUTINE dia_ar5_hst
SUBROUTINE dia_ar5_init
INTEGER :: inum
INTEGER :: ik
INTEGER :: ji, jj, jk
REAL(KIND = wp) :: zztmp
REAL(KIND = wp), ALLOCATABLE, DIMENSION(:, :, :, :) :: zsaldta
l_ar5 = .FALSE.
IF (iom_use('voltot') .OR. iom_use('sshtot') .OR. iom_use('sshdyn') .OR. iom_use('masstot') .OR. iom_use('temptot') .OR. iom_use('saltot') .OR. iom_use('botpres') .OR. iom_use('sshthster') .OR. iom_use('sshsteric')) L_ar5 = .TRUE.
IF (l_ar5) THEN
IF (dia_ar5_alloc() /= 0) CALL ctl_stop('STOP', 'dia_ar5_init : unable to allocate arrays')
!$ACC KERNELS
area(:, :) = e1e2t(:, :) * tmask_i(:, :)
!$ACC END KERNELS
area_tot = SUM(area(:, :))
IF (lk_mpp) CALL mpp_sum(area_tot)
!$ACC KERNELS
vol0 = 0._wp
thick0(:, :) = 0._wp
!$ACC END KERNELS
DO jk = 1, jpkm1
vol0 = vol0 + SUM(area(:, :) * tmask(:, :, jk) * e3t_0(:, :, jk))
!$ACC KERNELS
thick0(:, :) = thick0(:, :) + tmask_i(:, :) * tmask(:, :, jk) * e3t_0(:, :, jk)
!$ACC END KERNELS
END DO
IF (lk_mpp) CALL mpp_sum(vol0)
IF (iom_use('sshthster')) THEN
ALLOCATE(zsaldta(jpi, jpj, jpj, jpts))
CALL iom_open('sali_ref_clim_monthly', inum)
CALL iom_get(inum, jpdom_data, 'vosaline', zsaldta(:, :, :, 1), 1)
CALL iom_get(inum, jpdom_data, 'vosaline', zsaldta(:, :, :, 2), 12)
CALL iom_close(inum)
!$ACC KERNELS
sn0(:, :, :) = 0.5_wp * (zsaldta(:, :, :, 1) + zsaldta(:, :, :, 2))
sn0(:, :, :) = sn0(:, :, :) * tmask(:, :, :)
!$ACC END KERNELS
IF (ln_zps) THEN
!$ACC KERNELS
DO jj = 1, jpj
DO ji = 1, jpi
ik = mbkt(ji, jj)
IF (ik > 1) THEN
zztmp = (gdept_1d(ik) - gdept_0(ji, jj, ik)) / (gdept_1d(ik) - gdept_1d(ik - 1))
sn0(ji, jj, ik) = (1._wp - zztmp) * sn0(ji, jj, ik) + zztmp * sn0(ji, jj, ik - 1)
END IF
END DO
END DO
!$ACC END KERNELS
END IF
DEALLOCATE(zsaldta)
END IF
END IF
END SUBROUTINE dia_ar5_init
END MODULE diaar5