1 | #include "cppdefs.h"
|
---|
2 | #define MPDATA_HOT
|
---|
3 |
|
---|
4 | MODULE mpdata_adiff_mod
|
---|
5 | #if defined NONLINEAR && defined TS_MPDATA && defined SOLVE3D
|
---|
6 | !
|
---|
7 | !svn $Id: mpdata_adiff.F 732 2008-09-07 01:55:51Z jcwarner $
|
---|
8 | !================================================== Hernan G. Arango ===
|
---|
9 | ! Copyright (c) 2002-2014 The ROMS/TOMS Group John C. Warner !
|
---|
10 | ! Licensed under a MIT/X style license !
|
---|
11 | ! See License_ROMS.txt !
|
---|
12 | !=======================================================================
|
---|
13 | ! !
|
---|
14 | ! This routine computes anti-diffusive velocities to correct tracer !
|
---|
15 | ! advection using MPDATA Recursive method. !
|
---|
16 | ! !
|
---|
17 | ! On Output: !
|
---|
18 | ! !
|
---|
19 | ! Ua Andi-diffusive velocity in the XI-direction (m/s). !
|
---|
20 | ! Va Anti-diffusive velocity in the ETA-direction (m/s). !
|
---|
21 | ! Wa Anti-diffusive velocity in the S-direction (m3/s). !
|
---|
22 | ! !
|
---|
23 | ! Reference: !
|
---|
24 | ! !
|
---|
25 | ! Margolin, L. and P.K. Smolarkiewicz, 1998: Antidiffusive !
|
---|
26 | ! velocities for multipass donor cell advection, SIAM J. !
|
---|
27 | ! Sci. Comput., 907-929. !
|
---|
28 | ! !
|
---|
29 | !=======================================================================
|
---|
30 | !
|
---|
31 | implicit none
|
---|
32 |
|
---|
33 | PUBLIC :: mpdata_adiff_tile
|
---|
34 |
|
---|
35 | CONTAINS
|
---|
36 | !
|
---|
37 | !***********************************************************************
|
---|
38 | SUBROUTINE mpdata_adiff_tile (ng, tile, &
|
---|
39 | & LBi, UBi, LBj, UBj, &
|
---|
40 | & IminS, ImaxS, JminS, JmaxS, &
|
---|
41 | # ifdef MASKING
|
---|
42 | & rmask, umask, vmask, &
|
---|
43 | # endif
|
---|
44 | # ifdef WET_DRY
|
---|
45 | & rmask_wet, umask_wet, vmask_wet, &
|
---|
46 | # endif
|
---|
47 | & pm, pn, omn, om_u, on_v, &
|
---|
48 | & z_r, oHz, &
|
---|
49 | & Huon, Hvom, w, &
|
---|
50 | # ifdef WEC_VF
|
---|
51 | & W_stokes, &
|
---|
52 | # endif
|
---|
53 | & t, Ta, Ua, Va, Wa)
|
---|
54 | !***********************************************************************
|
---|
55 | !
|
---|
56 | USE mod_param
|
---|
57 | USE mod_ncparam
|
---|
58 | USE mod_scalars
|
---|
59 | !
|
---|
60 | ! Imported variable declarations.
|
---|
61 | !
|
---|
62 | integer, intent(in) :: ng, tile
|
---|
63 | integer, intent(in) :: LBi, UBi, LBj, UBj
|
---|
64 | integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
|
---|
65 | !
|
---|
66 | # ifdef ASSUMED_SHAPE
|
---|
67 | # ifdef MASKING
|
---|
68 | real(r8), intent(in) :: rmask(LBi:,LBj:)
|
---|
69 | real(r8), intent(in) :: umask(LBi:,LBj:)
|
---|
70 | real(r8), intent(in) :: vmask(LBi:,LBj:)
|
---|
71 | # endif
|
---|
72 | # ifdef WET_DRY
|
---|
73 | real(r8), intent(in) :: rmask_wet(LBi:,LBj:)
|
---|
74 | real(r8), intent(in) :: umask_wet(LBi:,LBj:)
|
---|
75 | real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
|
---|
76 | # endif
|
---|
77 | real(r8), intent(in) :: pm(LBi:,LBj:)
|
---|
78 | real(r8), intent(in) :: pn(LBi:,LBj:)
|
---|
79 | real(r8), intent(in) :: omn(LBi:,LBj:)
|
---|
80 | real(r8), intent(in) :: om_u(LBi:,LBj:)
|
---|
81 | real(r8), intent(in) :: on_v(LBi:,LBj:)
|
---|
82 | real(r8), intent(in) :: z_r(LBi:,LBj:,:)
|
---|
83 | real(r8), intent(in) :: oHz(IminS:,JminS:,:)
|
---|
84 | real(r8), intent(in) :: Huon(LBi:,LBj:,:)
|
---|
85 | real(r8), intent(in) :: Hvom(LBi:,LBj:,:)
|
---|
86 | real(r8), intent(in) :: t(LBi:,LBj:,:)
|
---|
87 | real(r8), intent(in) :: w(LBi:,LBj:,0:)
|
---|
88 | # ifdef WEC_VF
|
---|
89 | real(r8), intent(in) :: W_stokes(LBi:,LBj:,0:)
|
---|
90 | # endif
|
---|
91 | real(r8), intent(inout) :: Ta(IminS:,JminS:,:)
|
---|
92 |
|
---|
93 | real(r8), intent(out) :: Ua(IminS:,JminS:,:)
|
---|
94 | real(r8), intent(out) :: Va(IminS:,JminS:,:)
|
---|
95 | real(r8), intent(out) :: Wa(IminS:,JminS:,0:)
|
---|
96 | # else
|
---|
97 | # ifdef MASKING
|
---|
98 | real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
|
---|
99 | real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
|
---|
100 | real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
|
---|
101 | # endif
|
---|
102 | # ifdef WET_DRY
|
---|
103 | real(r8), intent(in) :: rmask_wet(LBi:UBi,LBj:UBj)
|
---|
104 | real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
|
---|
105 | real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
|
---|
106 | # endif
|
---|
107 | real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
|
---|
108 | real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
|
---|
109 | real(r8), intent(in) :: omn(LBi:UBi,LBj:UBj)
|
---|
110 | real(r8), intent(in) :: om_u(LBi:UBi,LBj:UBj)
|
---|
111 | real(r8), intent(in) :: on_v(LBi:UBi,LBj:UBj)
|
---|
112 | real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
|
---|
113 | real(r8), intent(in) :: oHz(IminS:ImaxS,JminS:JmaxS,N(ng))
|
---|
114 | real(r8), intent(in) :: Huon(LBi:UBi,LBj:UBj,N(ng))
|
---|
115 | real(r8), intent(in) :: Hvom(LBi:UBi,LBj:UBj,N(ng))
|
---|
116 | real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng))
|
---|
117 | real(r8), intent(in) :: w(LBi:UBi,LBj:UBj,0:N(ng))
|
---|
118 | # ifdef WEC_VF
|
---|
119 | real(r8), intent(in) :: W_stokes(LBi:UBi,LBj:UBj,0:N(ng))
|
---|
120 | # endif
|
---|
121 | real(r8), intent(inout) :: Ta(IminS:ImaxS,JminS:JmaxS,N(ng))
|
---|
122 |
|
---|
123 | real(r8), intent(out) :: Ua(IminS:ImaxS,JminS:JmaxS,N(ng))
|
---|
124 | real(r8), intent(out) :: Va(IminS:ImaxS,JminS:JmaxS,N(ng))
|
---|
125 | real(r8), intent(out) :: Wa(IminS:ImaxS,JminS:JmaxS,0:N(ng))
|
---|
126 | # endif
|
---|
127 | !
|
---|
128 | ! Local variable declarations.
|
---|
129 | !
|
---|
130 | integer :: i, is, j, k
|
---|
131 |
|
---|
132 | real(r8), parameter :: Large = 1.0E+20_r8
|
---|
133 | real(r8), parameter :: eps = 1.0E-18_r8
|
---|
134 | real(r8), parameter :: eps2= 1.0E-10_r8
|
---|
135 |
|
---|
136 | real(r8) :: A, B, Tmax, Tmin, Um, Vm, X, Y, Z
|
---|
137 | real(r8) :: cff, cff1, cff2, sig_alfa, fac
|
---|
138 | # ifdef MPDATA_HOT
|
---|
139 | real(r8) :: AA, BB, CC, AB, AC, BC
|
---|
140 | real(r8) :: XX, YY, ZZ, XY, XZ, YZ
|
---|
141 | real(r8) :: sig_beta, sig_gama
|
---|
142 | real(r8) :: sig_a, sig_b, sig_c, sig_d, sig_e, sig_f
|
---|
143 | # endif
|
---|
144 |
|
---|
145 | real(r8), dimension(IminS:ImaxS,N(ng)) :: C
|
---|
146 | real(r8), dimension(IminS:ImaxS,N(ng)) :: Wm
|
---|
147 |
|
---|
148 | real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: mask_dn
|
---|
149 | real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: mask_up
|
---|
150 |
|
---|
151 | real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: beta_dn
|
---|
152 | real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: beta_up
|
---|
153 | real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: odz
|
---|
154 |
|
---|
155 | # include "set_bounds.h"
|
---|
156 | # ifdef TS_MPDATA_LIMIT
|
---|
157 | fac=0.25_r8
|
---|
158 | # else
|
---|
159 | fac=1.0_r8
|
---|
160 | # endif
|
---|
161 | !
|
---|
162 | !-----------------------------------------------------------------------
|
---|
163 | ! Compute anti-diffusive MPDATA velocities (Ua,Va,Wa).
|
---|
164 | !-----------------------------------------------------------------------
|
---|
165 | !
|
---|
166 | ! Set boundary conditions to advected tracer, Ta.
|
---|
167 | !
|
---|
168 | IF (.not.EWperiodic(ng)) THEN
|
---|
169 | IF (DOMAIN(ng)%Western_Edge(tile)) THEN
|
---|
170 | DO k=1,N(ng)
|
---|
171 | DO j=JstrVm2,Jendp2i
|
---|
172 | Ta(Istr-1,j,k)=Ta(Istr,j,k)
|
---|
173 | END DO
|
---|
174 | END DO
|
---|
175 | END IF
|
---|
176 | IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN
|
---|
177 | DO k=1,N(ng)
|
---|
178 | DO j=JstrVm2,Jendp2i
|
---|
179 | Ta(Iend+1,j,k)=Ta(Iend,j,k)
|
---|
180 | END DO
|
---|
181 | END DO
|
---|
182 | END IF
|
---|
183 | END IF
|
---|
184 |
|
---|
185 | IF (.not.NSperiodic(ng)) THEN
|
---|
186 | IF (DOMAIN(ng)%Southern_Edge(tile)) THEN
|
---|
187 | DO k=1,N(ng)
|
---|
188 | DO i=IstrUm2,Iendp2i
|
---|
189 | Ta(i,Jstr-1,k)=Ta(i,Jstr,k)
|
---|
190 | END DO
|
---|
191 | END DO
|
---|
192 | END IF
|
---|
193 | IF (DOMAIN(ng)%Northern_Edge(tile)) THEN
|
---|
194 | DO k=1,N(ng)
|
---|
195 | DO i=IstrUm2,Iendp2i
|
---|
196 | Ta(i,Jend+1,k)=Ta(i,Jend,k)
|
---|
197 | END DO
|
---|
198 | END DO
|
---|
199 | END IF
|
---|
200 | END IF
|
---|
201 |
|
---|
202 | IF (.not.(EWperiodic(ng).or.NSperiodic(ng))) THEN
|
---|
203 | IF (DOMAIN(ng)%SouthWest_Corner(tile)) THEN
|
---|
204 | DO k=1,N(ng)
|
---|
205 | Ta(Istr-1,Jstr-1,k)=0.5_r8*(Ta(Istr ,Jstr-1,k)+ &
|
---|
206 | & Ta(Istr-1,Jstr ,k))
|
---|
207 | END DO
|
---|
208 | END IF
|
---|
209 | IF (DOMAIN(ng)%SouthEast_Corner(tile)) THEN
|
---|
210 | DO k=1,N(ng)
|
---|
211 | Ta(Iend+1,Jstr-1,k)=0.5_r8*(Ta(Iend+1,Jstr ,k)+ &
|
---|
212 | & Ta(Iend ,Jstr-1,k))
|
---|
213 | END DO
|
---|
214 | END IF
|
---|
215 | IF (DOMAIN(ng)%NorthWest_Corner(tile)) THEN
|
---|
216 | DO k=1,N(ng)
|
---|
217 | Ta(Istr-1,Jend+1,k)=0.5_r8*(Ta(Istr-1,Jend ,k)+ &
|
---|
218 | & Ta(Istr ,Jend+1,k))
|
---|
219 | END DO
|
---|
220 | END IF
|
---|
221 | IF (DOMAIN(ng)%NorthEast_Corner(tile)) THEN
|
---|
222 | DO k=1,N(ng)
|
---|
223 | Ta(Iend+1,Jend+1,k)=0.5_r8*(Ta(Iend+1,Jend ,k)+ &
|
---|
224 | & Ta(Iend ,Jend+1,k))
|
---|
225 | END DO
|
---|
226 | END IF
|
---|
227 | END IF
|
---|
228 | !
|
---|
229 | ! Compute inverse vertical grid spacing at W-points.
|
---|
230 | !
|
---|
231 | DO k=1,N(ng)-1
|
---|
232 | DO j=Jstrm2,Jendp2
|
---|
233 | DO i=Istrm2,Iendp2
|
---|
234 | odz(i,j,k)=1.0_r8/(z_r(i,j,k+1)-z_r(i,j,k))
|
---|
235 | END DO
|
---|
236 | END DO
|
---|
237 | END DO
|
---|
238 | cff=1.0_r8/dt(ng)
|
---|
239 | !
|
---|
240 | ! Compute nondimensional U-antidiffusive velocities, Ua. If applicable,
|
---|
241 | ! retain up to third-order terms of the power series.
|
---|
242 | !
|
---|
243 | DO j=JstrV-1,Jendp1
|
---|
244 | k=1
|
---|
245 | ! DO i=IstrU-1,Iendp2
|
---|
246 | DO i=IstrUm1,Iendp2
|
---|
247 | C(i,k)=0.25_r8* &
|
---|
248 | & ((Ta(i ,j,k+1)-Ta(i ,j,k ))*odz(i ,j,k )+ &
|
---|
249 | & (Ta(i-1,j,k+1)-Ta(i-1,j,k ))*odz(i-1,j,k ))* &
|
---|
250 | & (z_r(i ,j,k+1)-z_r(i ,j,k)+ &
|
---|
251 | & z_r(i-1,j,k+1)-z_r(i-1,j,k))/ &
|
---|
252 | & (Ta(i-1,j,k)+Ta(i,j,k)+eps)
|
---|
253 | Wm(i,k)=0.25_r8*dt(ng)* &
|
---|
254 | & (w(i-1,j,k )*odz(i-1,j,k)*pm(i-1,j)*pn(i-1,j)+ &
|
---|
255 | & w(i ,j,k )*odz(i ,j,k)*pm(i ,j)*pn(i ,j))
|
---|
256 | # ifdef WEC_VF
|
---|
257 | Wm(i,k)=Wm(i,k)+ &
|
---|
258 | & 0.25_r8*dt(ng)* &
|
---|
259 | & (W_stokes(i-1,j,k )*odz(i-1,j,k)*pm(i-1,j)*pn(i-1,j)+&
|
---|
260 | & W_stokes(i ,j,k )*odz(i ,j,k)*pm(i ,j)*pn(i ,j))
|
---|
261 | # endif
|
---|
262 | END DO
|
---|
263 | DO k=2,N(ng)-1
|
---|
264 | DO i=IstrU-1,Iendp2
|
---|
265 | C(i,k)=0.0625_r8* &
|
---|
266 | & ((Ta(i ,j,k+1)-Ta(i ,j,k ))*odz(i ,j,k )+ &
|
---|
267 | & (Ta(i ,j,k )-Ta(i ,j,k-1))*odz(i ,j,k-1)+ &
|
---|
268 | & (Ta(i-1,j,k+1)-Ta(i-1,j,k ))*odz(i-1,j,k )+ &
|
---|
269 | & (Ta(i-1,j,k )-Ta(i-1,j,k-1))*odz(i-1,j,k-1))* &
|
---|
270 | & (z_r(i ,j,k+1)-z_r(i ,j,k-1)+ &
|
---|
271 | & z_r(i-1,j,k+1)-z_r(i-1,j,k-1))/ &
|
---|
272 | & (Ta(i-1,j,k)+Ta(i,j,k)+eps)
|
---|
273 | Wm(i,k)=0.25_r8*dt(ng)* &
|
---|
274 | & ((w(i-1,j,k-1)*odz(i-1,j,k-1)+ &
|
---|
275 | & w(i-1,j,k )*odz(i-1,j,k ))*pm(i-1,j)*pn(i-1,j)+ &
|
---|
276 | & (w(i ,j,k )*odz(i ,j,k )+ &
|
---|
277 | & w(i ,j,k-1)*odz(i ,j,k-1))*pm(i ,j)*pn(i ,j))
|
---|
278 | # ifdef WEC_VF
|
---|
279 | Wm(i,k)=Wm(i,k)+ &
|
---|
280 | & 0.25_r8*dt(ng)* &
|
---|
281 | & ((W_stokes(i-1,j,k-1)*odz(i-1,j,k-1)+ &
|
---|
282 | & W_stokes(i-1,j,k )*odz(i-1,j,k ))* &
|
---|
283 | & pm(i-1,j)*pn(i-1,j)+ &
|
---|
284 | & (W_stokes(i ,j,k )*odz(i ,j,k )+ &
|
---|
285 | & W_stokes(i ,j,k-1)*odz(i ,j,k-1))* &
|
---|
286 | & pm(i ,j)*pn(i ,j))
|
---|
287 | # endif
|
---|
288 | END DO
|
---|
289 | END DO
|
---|
290 | k=N(ng)
|
---|
291 | DO i=IstrU-1,Iendp2
|
---|
292 | C(i,k)=0.25_r8* &
|
---|
293 | & ((Ta(i ,j,k )-Ta(i ,j,k-1))*odz(i ,j,k-1)+ &
|
---|
294 | & (Ta(i-1,j,k )-Ta(i-1,j,k-1))*odz(i-1,j,k-1))* &
|
---|
295 | & (z_r(i ,j,k )-z_r(i ,j,k-1)+ &
|
---|
296 | & z_r(i-1,j,k )-z_r(i-1,j,k-1))/ &
|
---|
297 | & (Ta(i-1,j,k)+Ta(i,j,k)+eps)
|
---|
298 | Wm(i,k)=0.25_r8*dt(ng)* &
|
---|
299 | & (w(i-1,j,k-1)*odz(i-1,j,k-1)*pm(i-1,j)*pn(i-1,j)+ &
|
---|
300 | & w(i ,j,k-1)*odz(i ,j,k-1)*pm(i ,j)*pn(i ,j))
|
---|
301 | # ifdef WEC_VF
|
---|
302 | Wm(i,k)=Wm(i,k)+ &
|
---|
303 | & 0.25_r8*dt(ng)* &
|
---|
304 | & (W_stokes(i-1,j,k-1)*odz(i-1,j,k-1)* &
|
---|
305 | & pm(i-1,j)*pn(i-1,j)+ &
|
---|
306 | & W_stokes(i ,j,k-1)*odz(i ,j,k-1)* &
|
---|
307 | & pm(i ,j)*pn(i ,j))
|
---|
308 | # endif
|
---|
309 | END DO
|
---|
310 | DO k=1,N(ng)
|
---|
311 | DO i=IstrU-1,Iendp2
|
---|
312 | IF ((Ta(i-1,j,k).le.0.0_r8).or. &
|
---|
313 | & (Ta(i ,j,k).le.0.0_r8).or. &
|
---|
314 | & (ABS(Ta(i-1,j,k)-Ta(i,j,k)).le.eps2)) THEN
|
---|
315 | Ua(i,j,k)=0.0_r8
|
---|
316 | ELSE
|
---|
317 | A=(Ta(i,j,k)-Ta(i-1,j,k))/ &
|
---|
318 | & (Ta(i,j,k)+Ta(i-1,j,k)+eps)
|
---|
319 | # ifdef MASKING
|
---|
320 | B=0.03125_r8* &
|
---|
321 | & ((Ta(i ,j+1,k)-Ta(i ,j ,k))* &
|
---|
322 | & (pn(i ,j )+pn(i ,j+1))*vmask(i ,j+1)+ &
|
---|
323 | & (Ta(i ,j ,k)-Ta(i ,j-1,k))* &
|
---|
324 | & (pn(i ,j-1)+pn(i ,j ))*vmask(i ,j )+ &
|
---|
325 | & (Ta(i-1,j+1,k)-Ta(i-1,j ,k))* &
|
---|
326 | & (pn(i-1,j )+pn(i-1,j+1))*vmask(i-1,j+1)+ &
|
---|
327 | & (Ta(i-1,j ,k)-Ta(i-1,j-1,k))* &
|
---|
328 | & (pn(i-1,j-1)+pn(i-1,j ))*vmask(i-1,j ))
|
---|
329 | # else
|
---|
330 | B=0.03125_r8* &
|
---|
331 | & ((Ta(i ,j+1,k)-Ta(i ,j ,k))* &
|
---|
332 | & (pn(i ,j )+pn(i ,j+1))+ &
|
---|
333 | & (Ta(i ,j ,k)-Ta(i ,j-1,k))* &
|
---|
334 | & (pn(i ,j-1)+pn(i ,j ))+ &
|
---|
335 | & (Ta(i-1,j+1,k)-Ta(i-1,j ,k))* &
|
---|
336 | & (pn(i-1,j )+pn(i-1,j+1))+ &
|
---|
337 | & (Ta(i-1,j ,k)-Ta(i-1,j-1,k))* &
|
---|
338 | & (pn(i-1,j-1)+pn(i-1,j )))
|
---|
339 | # endif
|
---|
340 | B=B*(on_v(i ,j )+on_v(i ,j+1)+ &
|
---|
341 | & on_v(i-1,j )+on_v(i-1,j+1))/ &
|
---|
342 | & (Ta(i-1,j,k)+Ta(i,j,k)+eps)
|
---|
343 | !
|
---|
344 | Um=0.125_r8*Huon(i,j,k)* &
|
---|
345 | & dt(ng)*(pm(i,j)+pm(i-1,j))*(pn(i,j)+pn(i-1,j))* &
|
---|
346 | & (oHz(i-1,j,k)+oHz(i,j,k))
|
---|
347 | Vm=0.03125_r8*dt(ng)* &
|
---|
348 | & (Hvom(i-1,j ,k)*(pm(i-1,j)+pm(i-1,j-1))* &
|
---|
349 | & (pn(i-1,j)+pn(i-1,j-1))* &
|
---|
350 | & (oHz(i-1,j, k)+oHz(i-1,j-1,k))+ &
|
---|
351 | & Hvom(i-1,j+1,k)*(pm(i-1,j+1)+pm(i-1,j))* &
|
---|
352 | & (pn(i-1,j+1)+pn(i-1,j))* &
|
---|
353 | & (oHz(i-1,j+1,k)+oHz(i-1,j ,k))+ &
|
---|
354 | & Hvom(i ,j ,k)*(pm(i ,j)+pm(i ,j-1))* &
|
---|
355 | & (pn(i ,j)+pn(i ,j-1))* &
|
---|
356 | & (oHz(i ,j ,k)+oHz(i ,j-1,k))+ &
|
---|
357 | & Hvom(i ,j+1,k)*(pm(i ,j+1)+pm(i ,j))* &
|
---|
358 | & (pn(i ,j+1)+pn(i ,j))* &
|
---|
359 | & (oHz(i ,j+1,k)+oHz(i ,j ,k)))
|
---|
360 | !
|
---|
361 | X=(ABS(Um)-Um*Um)*A-B*Um*Vm-C(i,k)*Um*Wm(i,k)
|
---|
362 | Y=(ABS(Vm)-Vm*Vm)*B-A*Um*Vm-C(i,k)*Vm*Wm(i,k)
|
---|
363 | Z=(ABS(Wm(i,k))-Wm(i,k)*Wm(i,k))*C(i,k)- &
|
---|
364 | & A*Um*Wm(i,k)-B*Vm*Wm(i,k)
|
---|
365 |
|
---|
366 | # ifdef MPDATA_HOT
|
---|
367 | !
|
---|
368 | AA=A*A
|
---|
369 | BB=B*B
|
---|
370 | CC=C(i,k)*C(i,k)
|
---|
371 | AB=A*B
|
---|
372 | AC=A*C(i,k)
|
---|
373 | BC=B*C(i,k)
|
---|
374 |
|
---|
375 | XX=X*X
|
---|
376 | YY=Y*Y
|
---|
377 | ZZ=Z*Z
|
---|
378 | XY=X*Y
|
---|
379 | XZ=X*Z
|
---|
380 | YZ=Y*Z
|
---|
381 | # endif
|
---|
382 | !
|
---|
383 | sig_alfa=1.0_r8/(1.0_r8-ABS(A)+eps)
|
---|
384 | # ifdef MPDATA_HOT
|
---|
385 | sig_beta=-A/((1.0_r8-ABS(A))* &
|
---|
386 | & (1.0_r8-AA)+eps)
|
---|
387 | sig_gama=2.0_r8*ABS(AA*A)/((1.0_r8-ABS(A))* &
|
---|
388 | & (1.0_r8-AA)* &
|
---|
389 | & (1.0_r8-ABS(AA*A))+eps)
|
---|
390 | sig_a=-B/((1.0_r8-ABS(A))* &
|
---|
391 | & (1.0_r8-ABS(AB))+eps)
|
---|
392 | sig_b=AB/((1.0_r8-ABS(A))* &
|
---|
393 | & (1.0_r8-AA*ABS(B))+eps)* &
|
---|
394 | & (ABS(B)/(1.0_r8-ABS(AB)+eps)+ &
|
---|
395 | & 2.0_r8*A/(1.0_r8-AA+eps))
|
---|
396 | sig_c=ABS(A)*BB/((1.0_r8-ABS(A))* &
|
---|
397 | & (1.0_r8-BB*ABS(A))* &
|
---|
398 | & (1.0_r8-ABS(AB))+eps)
|
---|
399 | sig_d=-C(i,k)/((1.0_r8-ABS(A))* &
|
---|
400 | & (1.0_r8-ABS(AC))+eps)
|
---|
401 | sig_e=AC/((1.0_r8-ABS(A))* &
|
---|
402 | & (1.0_r8-AA*ABS(C(i,k)))+eps)* &
|
---|
403 | & (ABS(C(i,k))/(1.0_r8-ABS(AC)+eps)+ &
|
---|
404 | & 2.0_r8*A/(1.0_r8-AA+eps))
|
---|
405 | sig_f=ABS(A)*CC/((1.0_r8-ABS(A))* &
|
---|
406 | & (1.0_r8-CC*ABS(A))* &
|
---|
407 | & (1.0_r8-ABS(AC))+eps)
|
---|
408 |
|
---|
409 | Ua(i,j,k)=sig_alfa*X+ &
|
---|
410 | & sig_beta*XX+ &
|
---|
411 | & sig_gama*XX*X+ &
|
---|
412 | & sig_a*XY+ &
|
---|
413 | & sig_b*XX*Y+ &
|
---|
414 | & sig_c*X*YY+ &
|
---|
415 | & sig_d*XZ+ &
|
---|
416 | & sig_e*XX*Z+ &
|
---|
417 | & sig_f*X*ZZ
|
---|
418 | # else
|
---|
419 | Ua(i,j,k)=sig_alfa*X
|
---|
420 | # endif
|
---|
421 | !
|
---|
422 | ! Limit by physical velocity.
|
---|
423 | !
|
---|
424 | Ua(i,j,k)=MIN(ABS(Ua(i,j,k)),fac*ABS(Um))* &
|
---|
425 | & SIGN(1.0_r8,Ua(i,j,k))
|
---|
426 | # ifdef MASKING
|
---|
427 | Ua(i,j,k)=Ua(i,j,k)*umask(i,j)
|
---|
428 | # endif
|
---|
429 | # ifdef WET_DRY
|
---|
430 | Ua(i,j,k)=Ua(i,j,k)*umask_wet(i,j)
|
---|
431 | # endif
|
---|
432 | END IF
|
---|
433 | END DO
|
---|
434 | END DO
|
---|
435 | END DO
|
---|
436 | !
|
---|
437 | ! Compute nondimensional V-antidiffusive velocities, Va. If applicable,
|
---|
438 | ! retain up to third-order terms of the power series.
|
---|
439 | !
|
---|
440 | ! DO j=JstrV-1,Jendp2
|
---|
441 | DO j=JstrVm1,Jendp2
|
---|
442 | k=1
|
---|
443 | DO i=IstrU-1,Iendp1
|
---|
444 | C(i,k)=0.25_r8* &
|
---|
445 | & ((Ta(i,j ,k+1)-Ta(i,j ,k ))*odz(i,j ,k )+ &
|
---|
446 | & (Ta(i,j-1,k+1)-Ta(i,j-1,k ))*odz(i,j-1,k ))* &
|
---|
447 | & (z_r(i,j ,k+1)-z_r(i,j ,k)+ &
|
---|
448 | & z_r(i,j-1,k+1)-z_r(i,j-1,k))/ &
|
---|
449 | & (Ta(i,j-1,k)+Ta(i,j,k)+eps)
|
---|
450 | Wm(i,k)=0.25_r8*dt(ng)* &
|
---|
451 | & (w(i,j-1,k )*odz(i,j-1,k )*pm(i,j-1)*pn(i,j-1)+ &
|
---|
452 | & w(i,j ,k )*odz(i,j ,k )*pm(i,j )*pn(i,j ))
|
---|
453 | # ifdef WEC_VF
|
---|
454 | Wm(i,k)=Wm(i,k)+ &
|
---|
455 | & 0.25_r8*dt(ng)* &
|
---|
456 | & (W_stokes(i,j-1,k )*odz(i,j-1,k )*pm(i,j-1)*pn(i,j-1)+&
|
---|
457 | & W_stokes(i,j ,k )*odz(i,j ,k )*pm(i,j )*pn(i,j ))
|
---|
458 | # endif
|
---|
459 | END DO
|
---|
460 | DO k=2,N(ng)-1
|
---|
461 | DO i=IstrU-1,Iendp1
|
---|
462 | C(i,k)=0.0625_r8* &
|
---|
463 | & ((Ta(i,j ,k+1)-Ta(i,j ,k ))*odz(i,j ,k )+ &
|
---|
464 | & (Ta(i,j ,k )-Ta(i,j ,k-1))*odz(i,j ,k-1)+ &
|
---|
465 | & (Ta(i,j-1,k+1)-Ta(i,j-1,k ))*odz(i,j-1,k )+ &
|
---|
466 | & (Ta(i,j-1,k )-Ta(i,j-1,k-1))*odz(i,j-1,k-1))* &
|
---|
467 | & (z_r(i,j ,k+1)-z_r(i,j ,k-1)+ &
|
---|
468 | & z_r(i,j-1,k+1)-z_r(i,j-1,k-1))/ &
|
---|
469 | & (Ta(i,j-1,k)+Ta(i,j,k)+eps)
|
---|
470 | Wm(i,k)=0.25_r8*dt(ng)* &
|
---|
471 | & ((w(i,j-1,k-1)*odz(i,j-1,k-1)+ &
|
---|
472 | & w(i,j-1,k )*odz(i,j-1,k ))*pm(i,j-1)*pn(i,j-1)+ &
|
---|
473 | & (w(i,j ,k )*odz(i,j ,k )+ &
|
---|
474 | & w(i,j ,k-1)*odz(i,j ,k-1))*pm(i,j )*pn(i,j ))
|
---|
475 | # ifdef WEC_VF
|
---|
476 | Wm(i,k)=Wm(i,k)+ &
|
---|
477 | & 0.25_r8*dt(ng)* &
|
---|
478 | & ((W_stokes(i,j-1,k-1)*odz(i,j-1,k-1)+ &
|
---|
479 | & W_stokes(i,j-1,k )*odz(i,j-1,k ))* &
|
---|
480 | & pm(i,j-1)*pn(i,j-1)+ &
|
---|
481 | & (W_stokes(i,j ,k )*odz(i,j ,k )+ &
|
---|
482 | & W_stokes(i,j ,k-1)*odz(i,j ,k-1))* &
|
---|
483 | & pm(i,j )*pn(i,j ))
|
---|
484 | # endif
|
---|
485 | END DO
|
---|
486 | END DO
|
---|
487 | k=N(ng)
|
---|
488 | DO i=IstrU-1,Iendp1
|
---|
489 | C(i,k)=0.25_r8* &
|
---|
490 | & ((Ta(i,j ,k )-Ta(i,j ,k-1))*odz(i,j ,k-1)+ &
|
---|
491 | & (Ta(i,j-1,k )-Ta(i,j-1,k-1))*odz(i,j-1,k-1))* &
|
---|
492 | & (z_r(i,j ,k )-z_r(i,j ,k-1)+ &
|
---|
493 | & z_r(i,j-1,k )-z_r(i,j-1,k-1))/ &
|
---|
494 | & (Ta(i,j-1,k)+Ta(i,j,k)+eps)
|
---|
495 | Wm(i,k)=0.25_r8*dt(ng)* &
|
---|
496 | & (w(i,j-1,k-1)*odz(i,j-1,k-1)*pm(i,j-1)*pn(i,j-1)+ &
|
---|
497 | & w(i,j ,k-1)*odz(i,j ,k-1)*pm(i,j )*pn(i,j ))
|
---|
498 | # ifdef WEC_VF
|
---|
499 | Wm(i,k)=Wm(i,k)+ &
|
---|
500 | & 0.25_r8*dt(ng)* &
|
---|
501 | & (W_stokes(i,j-1,k-1)*odz(i,j-1,k-1)* &
|
---|
502 | & pm(i,j-1)*pn(i,j-1)+ &
|
---|
503 | & W_stokes(i,j ,k-1)*odz(i,j ,k-1)* &
|
---|
504 | & pm(i,j )*pn(i,j ))
|
---|
505 | # endif
|
---|
506 | END DO
|
---|
507 | DO k=1,N(ng)
|
---|
508 | DO i=IstrU-1,Iendp1
|
---|
509 | IF ((Ta(i,j-1,k).le.0.0_r8).or. &
|
---|
510 | & (Ta(i,j ,k).le.0.0_r8).or. &
|
---|
511 | & (ABS(Ta(i,j-1,k)-Ta(i,j,k)).le.eps2)) THEN
|
---|
512 | Va(i,j,k)=0.0_r8
|
---|
513 | ELSE
|
---|
514 | # ifdef MASKING
|
---|
515 | A=0.03125_r8* &
|
---|
516 | & ((Ta(i+1,j ,k)-Ta(i ,j ,k))* &
|
---|
517 | & (pm(i+1,j )+pm(i ,j ))*umask(i+1,j )+ &
|
---|
518 | & (Ta(i ,j ,k)-Ta(i-1,j ,k))* &
|
---|
519 | & (pm(i-1,j )+pm(i ,j ))*umask(i ,j )+ &
|
---|
520 | & (Ta(i+1,j-1,k)-Ta(i ,j-1,k))* &
|
---|
521 | & (pm(i+1,j-1)+pm(i ,j-1))*umask(i+1,j-1)+ &
|
---|
522 | & (Ta(i ,j-1,k)-Ta(i-1,j-1,k))* &
|
---|
523 | & (pm(i-1,j-1)+pm(i ,j-1))*umask(i ,j-1))
|
---|
524 | # else
|
---|
525 | A=0.03125_r8* &
|
---|
526 | & ((Ta(i+1,j ,k)-Ta(i ,j ,k))* &
|
---|
527 | & (pm(i+1,j )+pm(i ,j ))+ &
|
---|
528 | & (Ta(i ,j ,k)-Ta(i-1,j ,k))* &
|
---|
529 | & (pm(i-1,j )+pm(i ,j ))+ &
|
---|
530 | & (Ta(i+1,j-1,k)-Ta(i ,j-1,k))* &
|
---|
531 | & (pm(i+1,j-1)+pm(i ,j-1))+ &
|
---|
532 | & (Ta(i ,j-1,k)-Ta(i-1,j-1,k))* &
|
---|
533 | & (pm(i-1,j-1)+pm(i ,j-1)))
|
---|
534 | # endif
|
---|
535 | A=A*(om_u(i ,j )+om_u(i+1,j )+ &
|
---|
536 | & om_u(i ,j-1)+om_u(i+1,j-1))/ &
|
---|
537 | & (Ta(i ,j-1,k)+Ta(i,j,k)+eps)
|
---|
538 |
|
---|
539 | B=(Ta(i,j,k)-Ta(i,j-1,k))/ &
|
---|
540 | & (Ta(i,j,k)+Ta(i,j-1,k)+eps)
|
---|
541 | !
|
---|
542 | Um=0.03125_r8*dt(ng)* &
|
---|
543 | & (Huon(i+1,j ,k)*(pm(i+1,j)+pm(i,j))* &
|
---|
544 | & (pn(i+1,j)+pn(i,j))* &
|
---|
545 | & (oHz(i+1,j ,k)+oHz(i,j ,k))+ &
|
---|
546 | & Huon(i+1,j-1,k)*(pm(i+1,j-1)+pm(i,j-1))* &
|
---|
547 | & (pn(i+1,j-1)+pn(i,j-1))* &
|
---|
548 | & (oHz(i+1,j-1,k)+oHz(i,j-1,k))+ &
|
---|
549 | & Huon(i ,j ,k)*(pm(i-1,j)+pm(i,j))* &
|
---|
550 | & (pn(i-1,j)+pn(i,j))* &
|
---|
551 | & (oHz(i-1,j ,k)+oHz(i,j ,k))+ &
|
---|
552 | & Huon(i ,j-1,k)*(pm(i-1,j-1)+pm(i,j-1))* &
|
---|
553 | & (pn(i-1,j-1)+pn(i,j-1))* &
|
---|
554 | & (oHz(i-1,j-1,k)+oHz(i,j-1,k)))
|
---|
555 | Vm=0.125_r8*Hvom(i,j,k)* &
|
---|
556 | & dt(ng)*(pn(i,j-1)+pn(i,j))*(pm(i,j-1)+pm(i,j))* &
|
---|
557 | & (oHz(i,j-1,k)+oHz(i,j,k))
|
---|
558 | !
|
---|
559 | X=(ABS(Um)-Um*Um)*A-B*Um*Vm-C(i,k)*Um*Wm(i,k)
|
---|
560 | Y=(ABS(Vm)-Vm*Vm)*B-A*Um*Vm-C(i,k)*Vm*Wm(i,k)
|
---|
561 | Z=(ABS(Wm(i,k))-Wm(i,k)*Wm(i,k))*C(i,k)- &
|
---|
562 | & A*Um*Wm(i,k)-B*Vm*Wm(i,k)
|
---|
563 |
|
---|
564 | # ifdef MPDATA_HOT
|
---|
565 | !
|
---|
566 | AA=A*A
|
---|
567 | BB=B*B
|
---|
568 | CC=C(i,k)*C(i,k)
|
---|
569 | AB=A*B
|
---|
570 | AC=A*C(i,k)
|
---|
571 | BC=B*C(i,k)
|
---|
572 |
|
---|
573 | XX=X*X
|
---|
574 | YY=Y*Y
|
---|
575 | ZZ=Z*Z
|
---|
576 | XY=X*Y
|
---|
577 | XZ=X*Z
|
---|
578 | YZ=Y*Z
|
---|
579 | # endif
|
---|
580 | !
|
---|
581 | sig_alfa=1.0_r8/(1.0_r8-ABS(B)+eps)
|
---|
582 | # ifdef MPDATA_HOT
|
---|
583 | sig_beta=-B/((1.0_r8-ABS(B))* &
|
---|
584 | & (1.0_r8-BB)+eps)
|
---|
585 | sig_gama=2.0_r8*ABS(BB*B)/((1.0_r8-ABS(B))* &
|
---|
586 | & (1.0_r8-BB)* &
|
---|
587 | & (1.0_r8-ABS(BB*B))+eps)
|
---|
588 | sig_a=-A/((1.0_r8-ABS(B))* &
|
---|
589 | & (1.0_r8-ABS(AB))+eps)
|
---|
590 | sig_b=AB/((1.0_r8-ABS(B))* &
|
---|
591 | & (1.0_r8-BB*ABS(A))+eps)* &
|
---|
592 | & (ABS(A)/(1.0_r8-ABS(AB)+eps)+ &
|
---|
593 | & 2.0_r8*B/(1.0_r8-BB+eps))
|
---|
594 | sig_c=ABS(B)*AA/((1.0_r8-ABS(B))* &
|
---|
595 | & (1.0_r8-AA*ABS(B))* &
|
---|
596 | & (1.0_r8-ABS(AB))+eps)
|
---|
597 | sig_d=-C(i,k)/((1.0_r8-ABS(B))* &
|
---|
598 | & (1.0_r8-ABS(BC))+eps)
|
---|
599 | sig_e=BC/((1.0_r8-ABS(B))* &
|
---|
600 | & (1.0_r8-BB*ABS(C(i,k)))+eps)* &
|
---|
601 | & (ABS(C(i,k))/(1.0_r8-ABS(BC)+eps)+ &
|
---|
602 | & 2.0_r8*B/(1.0_r8-BB+eps))
|
---|
603 | sig_f=ABS(B)*CC/((1.0_r8-ABS(B))* &
|
---|
604 | & (1.0_r8-CC*ABS(B))* &
|
---|
605 | & (1.0_r8-ABS(BC))+eps)
|
---|
606 |
|
---|
607 | Va(i,j,k)=sig_alfa*Y+ &
|
---|
608 | & sig_beta*YY+ &
|
---|
609 | & sig_gama*YY*Y+ &
|
---|
610 | & sig_a*XY+ &
|
---|
611 | & sig_b*Y*XX+ &
|
---|
612 | & sig_c*YY*X+ &
|
---|
613 | & sig_d*YZ+ &
|
---|
614 | & sig_e*YY*Z+ &
|
---|
615 | & sig_f*Y*ZZ
|
---|
616 | # else
|
---|
617 | Va(i,j,k)=sig_alfa*Y
|
---|
618 | # endif
|
---|
619 | !
|
---|
620 | ! Limit by physical velocity.
|
---|
621 | !
|
---|
622 | Va(i,j,k)=MIN(ABS(Va(i,j,k)),fac*ABS(Vm))* &
|
---|
623 | & SIGN(1.0_r8,Va(i,j,k))
|
---|
624 | # ifdef MASKING
|
---|
625 | Va(i,j,k)=Va(i,j,k)*vmask(i,j)
|
---|
626 | # endif
|
---|
627 | # ifdef WET_DRY
|
---|
628 | Va(i,j,k)=Va(i,j,k)*vmask_wet(i,j)
|
---|
629 | # endif
|
---|
630 | END IF
|
---|
631 | END DO
|
---|
632 | END DO
|
---|
633 | END DO
|
---|
634 | !
|
---|
635 | ! Apply boundary conditions to anti-diffusive velocities.
|
---|
636 | !
|
---|
637 | IF (.not.EWperiodic(ng)) THEN
|
---|
638 | IF (DOMAIN(ng)%Western_Edge(tile)) THEN
|
---|
639 | IF (LBC(iwest,isBu3d,ng)%closed) THEN
|
---|
640 | DO k=1,N(ng)
|
---|
641 | ! DO j=Jstr,Jend
|
---|
642 | DO j=Jstrm1,Jendp1
|
---|
643 | Ua(Istr,j,k)=0.0_r8
|
---|
644 | END DO
|
---|
645 | END DO
|
---|
646 | ELSE
|
---|
647 | DO k=1,N(ng)
|
---|
648 | ! DO j=Jstr,Jend
|
---|
649 | DO j=Jstrm1,Jendp1
|
---|
650 | Ua(Istr,j,k)=Ua(Istr+1,j,k)
|
---|
651 | END DO
|
---|
652 | END DO
|
---|
653 | END IF
|
---|
654 | END IF
|
---|
655 |
|
---|
656 | IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN
|
---|
657 | IF (LBC(ieast,isBu3d,ng)%closed) THEN
|
---|
658 | DO k=1,N(ng)
|
---|
659 | ! DO j=Jstr,Jend
|
---|
660 | DO j=Jstrm1,Jendp1
|
---|
661 | Ua(Iend+1,j,k)=0.0_r8
|
---|
662 | END DO
|
---|
663 | END DO
|
---|
664 | ELSE
|
---|
665 | DO k=1,N(ng)
|
---|
666 | ! DO j=Jstr,Jend
|
---|
667 | DO j=Jstrm1,Jendp1
|
---|
668 | Ua(Iend+1,j,k)=Ua(Iend,j,k)
|
---|
669 | END DO
|
---|
670 | END DO
|
---|
671 | END IF
|
---|
672 | END IF
|
---|
673 | END IF
|
---|
674 |
|
---|
675 | IF (.not.NSperiodic(ng)) THEN
|
---|
676 | IF (DOMAIN(ng)%Southern_Edge(tile)) THEN
|
---|
677 | IF (LBC(isouth,isBv3d,ng)%closed) THEN
|
---|
678 | DO k=1,N(ng)
|
---|
679 | ! DO i=Istr,Iend
|
---|
680 | DO i=Istrm1,Iendp1
|
---|
681 | Va(i,Jstr,k)=0.0_r8
|
---|
682 | END DO
|
---|
683 | END DO
|
---|
684 | ELSE
|
---|
685 | DO k=1,N(ng)
|
---|
686 | ! DO i=Istr,Iend
|
---|
687 | DO i=Istrm1,Iendp1
|
---|
688 | Va(i,Jstr,k)=Va(i,Jstr+1,k)
|
---|
689 | END DO
|
---|
690 | END DO
|
---|
691 | END IF
|
---|
692 | END IF
|
---|
693 |
|
---|
694 | IF (DOMAIN(ng)%Northern_Edge(tile)) THEN
|
---|
695 | IF (LBC(inorth,isBv3d,ng)%closed) THEN
|
---|
696 | DO k=1,N(ng)
|
---|
697 | ! DO i=Istr,Iend
|
---|
698 | DO i=Istrm1,Iendp1
|
---|
699 | Va(i,Jend+1,k)=0.0_r8
|
---|
700 | END DO
|
---|
701 | END DO
|
---|
702 | ELSE
|
---|
703 | DO k=1,N(ng)
|
---|
704 | ! DO i=Istr,Iend
|
---|
705 | DO i=Istrm1,Iendp1
|
---|
706 | Va(i,Jend+1,k)=Va(i,Jend,k)
|
---|
707 | END DO
|
---|
708 | END DO
|
---|
709 | END IF
|
---|
710 | END IF
|
---|
711 | END IF
|
---|
712 |
|
---|
713 | !
|
---|
714 | ! Compute nondimensional W-antidiffusive velocities, Wa. If applicable,
|
---|
715 | ! retain up to third-order terms of the power series.
|
---|
716 | !
|
---|
717 | DO j=JstrV-1,Jendp1
|
---|
718 | DO k=1,N(ng)-1
|
---|
719 | DO i=IstrU-1,Iendp1
|
---|
720 | IF ((Ta(i,j,k ).le.0.0_r8).or. &
|
---|
721 | & (Ta(i,j,k+1).le.0.0_r8).or. &
|
---|
722 | & (ABS(Ta(i,j,k)-Ta(i,j,k+1)).le.eps2)) THEN
|
---|
723 | Wa(i,j,k)=0.0_r8
|
---|
724 | ELSE
|
---|
725 | C(i,k)=(Ta(i,j,k+1)-Ta(i,j,k))/ &
|
---|
726 | & (Ta(i,j,k+1)+Ta(i,j,k)+eps)
|
---|
727 | # ifdef MASKING
|
---|
728 | A=0.0625_r8* &
|
---|
729 | & ((Ta(i+1,j,k+1)-Ta(i ,j,k+1))* &
|
---|
730 | & (pm(i+1,j )+pm(i ,j ))*umask(i+1,j)+ &
|
---|
731 | & (Ta(i ,j,k+1)-Ta(i-1,j,k+1))* &
|
---|
732 | & (pm(i ,j )+pm(i-1,j ))*umask(i ,j)+ &
|
---|
733 | & (Ta(i+1,j,k )-Ta(i ,j,k ))* &
|
---|
734 | & (pm(i+1,j )+pm(i ,j ))*umask(i+1,j)+ &
|
---|
735 | & (Ta(i ,j,k )-Ta(i-1,j,k ))* &
|
---|
736 | & (pm(i ,j )+pm(i-1,j ))*umask(i ,j))
|
---|
737 | B=0.0625_r8* &
|
---|
738 | & ((Ta(i,j+1,k+1)-Ta(i,j ,k+1))* &
|
---|
739 | & (pn(i ,j+1)+pn(i ,j ))*vmask(i,j+1)+ &
|
---|
740 | & (Ta(i,j ,k+1)-Ta(i,j-1,k+1))* &
|
---|
741 | & (pn(i ,j )+pn(i ,j-1))*vmask(i,j )+ &
|
---|
742 | & (Ta(i,j+1,k )-Ta(i,j ,k ))* &
|
---|
743 | & (pn(i ,j+1)+pn(i ,j ))*vmask(i,j+1)+ &
|
---|
744 | & (Ta(i,j ,k )-Ta(i,j-1,k ))* &
|
---|
745 | & (pn(i ,j )+pn(i ,j-1))*vmask(i,j ))
|
---|
746 | # else
|
---|
747 | A=0.0625_r8* &
|
---|
748 | & ((Ta(i+1,j,k+1)-Ta(i ,j,k+1))* &
|
---|
749 | & (pm(i+1,j )+pm(i ,j ))+ &
|
---|
750 | & (Ta(i ,j,k+1)-Ta(i-1,j,k+1))* &
|
---|
751 | & (pm(i ,j )+pm(i-1,j ))+ &
|
---|
752 | & (Ta(i+1,j,k )-Ta(i ,j,k ))* &
|
---|
753 | & (pm(i+1,j )+pm(i ,j ))+ &
|
---|
754 | & (Ta(i ,j,k )-Ta(i-1,j,k ))* &
|
---|
755 | & (pm(i ,j )+pm(i-1,j )))
|
---|
756 | B=0.0625_r8* &
|
---|
757 | & ((Ta(i,j+1,k+1)-Ta(i,j ,k+1))* &
|
---|
758 | & (pn(i ,j+1)+pn(i ,j ))+ &
|
---|
759 | & (Ta(i,j ,k+1)-Ta(i,j-1,k+1))* &
|
---|
760 | & (pn(i ,j )+pn(i ,j-1))+ &
|
---|
761 | & (Ta(i,j+1,k )-Ta(i,j ,k ))* &
|
---|
762 | & (pn(i ,j+1)+pn(i ,j ))+ &
|
---|
763 | & (Ta(i,j ,k )-Ta(i,j-1,k ))* &
|
---|
764 | & (pn(i ,j )+pn(i ,j-1)))
|
---|
765 | # endif
|
---|
766 | A=A*(om_u(i+1,j)+om_u(i ,j))/ &
|
---|
767 | & (Ta(i,j,k+1)+Ta(i,j,k)+eps)
|
---|
768 | B=B*(on_v(i,j+1)+on_v(i,j ))/ &
|
---|
769 | & (Ta(i,j,k+1)+Ta(i,j,k)+eps)
|
---|
770 | !
|
---|
771 | Um=0.03125_r8*dt(ng)* &
|
---|
772 | & (Huon(i ,j,k )*(pm(i,j)+pm(i-1,j))* &
|
---|
773 | & (pn(i,j)+pn(i-1,j))* &
|
---|
774 | & (oHz(i,j,k )+oHz(i-1,j,k ))+ &
|
---|
775 | & Huon(i ,j,k+1)*(pm(i,j)+pm(i-1,j))* &
|
---|
776 | & (pn(i,j)+pn(i-1,j))* &
|
---|
777 | & (oHz(i,j,k+1)+oHz(i-1,j,k+1))+ &
|
---|
778 | & Huon(i+1,j,k )*(pm(i,j)+pm(i+1,j))* &
|
---|
779 | & (pn(i,j)+pn(i+1,j))* &
|
---|
780 | & (oHz(i,j,k )+oHz(i+1,j,k ))+ &
|
---|
781 | & Huon(i+1,j,k+1)*(pm(i,j)+pm(i+1,j))* &
|
---|
782 | & (pn(i,j)+pn(i+1,j))* &
|
---|
783 | & (oHz(i,j,k+1)+oHz(i+1,j,k+1)))
|
---|
784 | Vm=0.03125_r8*dt(ng)* &
|
---|
785 | & (Hvom(i,j ,k )*(pm(i,j)+pm(i,j-1))* &
|
---|
786 | & (pn(i,j)+pn(i,j-1))* &
|
---|
787 | & (oHz(i,j,k )+oHz(i,j-1,k ))+ &
|
---|
788 | & Hvom(i,j ,k+1)*(pm(i,j)+pm(i,j-1))* &
|
---|
789 | & (pn(i,j)+pn(i,j-1))* &
|
---|
790 | & (oHz(i,j,k+1)+oHz(i,j-1,k+1))+ &
|
---|
791 | & Hvom(i,j+1,k )*(pm(i,j)+pm(i,j+1))* &
|
---|
792 | & (pn(i,j)+pn(i,j+1))* &
|
---|
793 | & (oHz(i,j,k )+oHz(i,j+1,k ))+ &
|
---|
794 | & Hvom(i,j+1,k+1)*(pm(i,j)+pm(i,j+1))* &
|
---|
795 | & (pn(i,j)+pn(i,j+1))* &
|
---|
796 | & (oHz(i,j,k+1)+oHz(i,j+1,k+1)))
|
---|
797 |
|
---|
798 | Wm(i,k)=w(i,j,k)*odz(i,j,k)*pm(i,j)*pn(i,j)*dt(ng)
|
---|
799 | # ifdef WEC_VF
|
---|
800 | Wm(i,k)=Wm(i,k)+W_stokes(i,j,k)*odz(i,j,k)* &
|
---|
801 | & pm(i,j)*pn(i,j)*dt(ng)
|
---|
802 | # endif
|
---|
803 | !
|
---|
804 | X=(ABS(Um)-Um*Um)*A-B*Um*Vm-C(i,k)*Um*Wm(i,k)
|
---|
805 | Y=(ABS(Vm)-Vm*Vm)*B-A*Um*Vm-C(i,k)*Vm*Wm(i,k)
|
---|
806 | Z=(ABS(Wm(i,k))-Wm(i,k)*Wm(i,k))*C(i,k)- &
|
---|
807 | & A*Um*Wm(i,k)-B*Vm*Wm(i,k)
|
---|
808 |
|
---|
809 | # ifdef MPDATA_HOT
|
---|
810 | !
|
---|
811 | AA=A*A
|
---|
812 | BB=B*B
|
---|
813 | CC=C(i,k)*C(i,k)
|
---|
814 | AB=A*B
|
---|
815 | AC=A*C(i,k)
|
---|
816 | BC=B*C(i,k)
|
---|
817 |
|
---|
818 | XX=X*X
|
---|
819 | YY=Y*Y
|
---|
820 | ZZ=Z*Z
|
---|
821 | XY=X*Y
|
---|
822 | XZ=X*Z
|
---|
823 | YZ=Y*Z
|
---|
824 | # endif
|
---|
825 | !
|
---|
826 | sig_alfa=1.0_r8/(1.0_r8-ABS(C(i,k))+eps)
|
---|
827 | # ifdef MPDATA_HOT
|
---|
828 | sig_beta=-C(i,k)/((1.0_r8-ABS(C(i,k)))* &
|
---|
829 | & (1.0_r8-CC)+eps)
|
---|
830 | sig_gama=2.0_r8*ABS(CC*C(i,k))/ &
|
---|
831 | & ((1.0_r8-ABS(C(i,k)))* &
|
---|
832 | & (1.0_r8-CC)* &
|
---|
833 | & (1.0_r8-ABS(CC*C(i,k)))+eps)
|
---|
834 | sig_a=-B/((1.0_r8-ABS(C(i,k)))* &
|
---|
835 | & (1.0_r8-ABS(BC))+eps)
|
---|
836 | sig_b=BC/((1.0_r8-ABS(C(i,k)))* &
|
---|
837 | & (1.0_r8-CC*ABS(B))+eps)* &
|
---|
838 | & (ABS(B)/(1.0_r8-ABS(BC)+eps)+ &
|
---|
839 | & 2.0_r8*C(i,k)/(1.0_r8-CC+eps))
|
---|
840 | sig_c=ABS(C(i,k))*BB/((1.0_r8-ABS(C(i,k)))* &
|
---|
841 | & (1.0_r8-B*B*ABS(C(i,k)))* &
|
---|
842 | & (1.0_r8-ABS(BC))+eps)
|
---|
843 | sig_d=-A/((1.0_r8-ABS(C(i,k)))* &
|
---|
844 | & (1.0_r8-ABS(AC))+eps)
|
---|
845 | sig_e=AC/((1.0_r8-ABS(C(i,k)))* &
|
---|
846 | & (1.0_r8-CC*ABS(A))+eps)* &
|
---|
847 | & (ABS(A)/(1.0_r8-ABS(AC)+eps)+ &
|
---|
848 | & 2.0_r8*C(i,k)/(1.0_r8-CC+eps))
|
---|
849 | sig_f=ABS(C(i,k))*AA/((1.0_r8-ABS(C(i,k)))* &
|
---|
850 | & (1.0_r8-AA*ABS(C(i,k)))* &
|
---|
851 | & (1.0_r8-ABS(AC))+eps)
|
---|
852 |
|
---|
853 | Wa(i,j,k)=sig_alfa*Z+ &
|
---|
854 | & sig_beta*ZZ+ &
|
---|
855 | & sig_gama*ZZ*Z+ &
|
---|
856 | & sig_a*YZ+ &
|
---|
857 | & sig_b*ZZ*Y+ &
|
---|
858 | & sig_c*Z*YY+ &
|
---|
859 | & sig_d*XZ+ &
|
---|
860 | & sig_e*ZZ*X+ &
|
---|
861 | & sig_f*Z*XX
|
---|
862 | # else
|
---|
863 | Wa(i,j,k)=sig_alfa*Z
|
---|
864 | # endif
|
---|
865 | !
|
---|
866 | ! Limit by physical velocity.
|
---|
867 | !
|
---|
868 | Wa(i,j,k)=MIN(ABS(Wa(i,j,k)),fac*ABS(Wm(i,k)))* &
|
---|
869 | & SIGN(1.0_r8,Wa(i,j,k))
|
---|
870 | # ifdef MASKING
|
---|
871 | Wa(i,j,k)=Wa(i,j,k)*rmask(i,j)
|
---|
872 | # endif
|
---|
873 | # ifdef WET_DRY
|
---|
874 | Wa(i,j,k)=Wa(i,j,k)*rmask_wet(i,j)
|
---|
875 | # endif
|
---|
876 | END IF
|
---|
877 | END DO
|
---|
878 | END DO
|
---|
879 | DO i=IstrU-1,Iendp1
|
---|
880 | Wa(i,j,0)=0.0_r8
|
---|
881 | Wa(i,j,N(ng))=0.0_r8
|
---|
882 | END DO
|
---|
883 | END DO
|
---|
884 | !
|
---|
885 | !-----------------------------------------------------------------------
|
---|
886 | ! Supress false oscillations in the solution by imposing appropriate
|
---|
887 | ! limits on the transport fluxes. Compute the UP and DOWN beta-ratios
|
---|
888 | ! described in Smolarkiewicz and Grabowski (1990).
|
---|
889 | !-----------------------------------------------------------------------
|
---|
890 | !
|
---|
891 | ! Build special mask array used to limit the UP and DOWN beta-ratios
|
---|
892 | ! to avoid including land/sea masking when computing limiting Tmin
|
---|
893 | ! and Tmax values. Notice that a zero Tmin due to land mask is not
|
---|
894 | ! considered.
|
---|
895 | !
|
---|
896 | DO j=Jstrm2,Jendp2
|
---|
897 | DO i=Istrm2,Iendp2
|
---|
898 | # ifdef MASKING
|
---|
899 | mask_up(i,j)=rmask(i,j)
|
---|
900 | mask_dn(i,j)=MAX(1.0_r8,MIN(Large,(1.0_r8-rmask(i,j))*Large))
|
---|
901 | # else
|
---|
902 | mask_up(i,j)=1.0_r8
|
---|
903 | mask_dn(i,j)=1.0_r8
|
---|
904 | # endif
|
---|
905 | END DO
|
---|
906 | END DO
|
---|
907 | !
|
---|
908 | ! Compute UP and DOWN beta-ratios.
|
---|
909 | !
|
---|
910 | DO j=JstrV-1,Jendp1
|
---|
911 | k=1
|
---|
912 | DO i=IstrU-1,Iendp1
|
---|
913 | Tmax=MAX(Ta(i-1,j ,k )*mask_up(i-1,j ), &
|
---|
914 | & t (i-1,j ,k )*mask_up(i-1,j ), &
|
---|
915 | & Ta(i ,j ,k )*mask_up(i ,j ), &
|
---|
916 | & t (i ,j ,k )*mask_up(i ,j ), &
|
---|
917 | & Ta(i+1,j ,k )*mask_up(i+1,j ), &
|
---|
918 | & t (i+1,j ,k )*mask_up(i+1,j ), &
|
---|
919 | & Ta(i ,j-1,k )*mask_up(i ,j-1), &
|
---|
920 | & t (i ,j-1,k )*mask_up(i ,j-1), &
|
---|
921 | & Ta(i ,j+1,k )*mask_up(i ,j+1), &
|
---|
922 | & t (i ,j+1,k )*mask_up(i ,j+1), &
|
---|
923 | & Ta(i ,j ,k+1)*mask_up(i ,j ), &
|
---|
924 | & t (i ,j ,k+1)*mask_up(i ,j ))
|
---|
925 | cff1=Ta(i-1,j ,k )*MAX(0.0_r8,Ua(i ,j ,k ))- &
|
---|
926 | & Ta(i+1,j ,k )*MIN(0.0_r8,Ua(i+1,j ,k ))+ &
|
---|
927 | & Ta(i ,j-1,k )*MAX(0.0_r8,Va(i ,j ,k ))- &
|
---|
928 | & Ta(i ,j+1,k )*MIN(0.0_r8,Va(i ,j+1,k ))- &
|
---|
929 | & Ta(i ,j ,k+1)*MIN(0.0_r8,Wa(i ,j ,k ))
|
---|
930 | beta_up(i,j,k)=(Tmax-Ta(i,j,k))/(cff1+eps)
|
---|
931 | !
|
---|
932 | Tmin=MIN(Ta(i-1,j ,k )*mask_dn(i-1,j ), &
|
---|
933 | & t (i-1,j ,k )*mask_dn(i-1,j ), &
|
---|
934 | & Ta(i ,j ,k )*mask_dn(i ,j ), &
|
---|
935 | & t (i ,j ,k )*mask_dn(i ,j ), &
|
---|
936 | & Ta(i+1,j ,k )*mask_dn(i+1,j ), &
|
---|
937 | & t (i+1,j ,k )*mask_dn(i+1,j ), &
|
---|
938 | & Ta(i ,j-1,k )*mask_dn(i ,j-1), &
|
---|
939 | & t (i ,j-1,k )*mask_dn(i ,j-1), &
|
---|
940 | & Ta(i ,j+1,k )*mask_dn(i ,j+1), &
|
---|
941 | & t (i ,j+1,k )*mask_dn(i ,j+1), &
|
---|
942 | & Ta(i ,j ,k+1)*mask_dn(i ,j ), &
|
---|
943 | & t (i ,j ,k+1)*mask_dn(i ,j ))
|
---|
944 | cff2=Ta(i ,j ,k )*MAX(0.0_r8,Ua(i+1,j ,k ))- &
|
---|
945 | & Ta(i ,j ,k )*MIN(0.0_r8,Ua(i ,j ,k ))+ &
|
---|
946 | & Ta(i ,j ,k )*MAX(0.0_r8,Va(i ,j+1,k ))- &
|
---|
947 | & Ta(i ,j ,k )*MIN(0.0_r8,Va(i ,j ,k ))+ &
|
---|
948 | & Ta(i ,j ,k )*MAX(0.0_r8,Wa(i ,j ,k ))
|
---|
949 | beta_dn(i,j,k)=(Ta(i,j,k)-Tmin)/(cff2+eps)
|
---|
950 | END DO
|
---|
951 | DO k=2,N(ng)-1
|
---|
952 | DO i=IstrU-1,Iendp1
|
---|
953 | Tmax=MAX(Ta(i-1,j ,k )*mask_up(i-1,j ), &
|
---|
954 | & t (i-1,j ,k )*mask_up(i-1,j ), &
|
---|
955 | & Ta(i ,j ,k )*mask_up(i ,j ), &
|
---|
956 | & t (i ,j ,k )*mask_up(i ,j ), &
|
---|
957 | & Ta(i+1,j ,k )*mask_up(i+1,j ), &
|
---|
958 | & t (i+1,j ,k )*mask_up(i+1,j ), &
|
---|
959 | & Ta(i ,j-1,k )*mask_up(i ,j-1), &
|
---|
960 | & t (i ,j-1,k )*mask_up(i ,j-1), &
|
---|
961 | & Ta(i ,j+1,k )*mask_up(i ,j+1), &
|
---|
962 | & t (i ,j+1,k )*mask_up(i ,j+1), &
|
---|
963 | & Ta(i ,j ,k-1)*mask_up(i ,j ), &
|
---|
964 | & t (i ,j ,k-1)*mask_up(i ,j ), &
|
---|
965 | & Ta(i ,j ,k+1)*mask_up(i ,j ), &
|
---|
966 | & t (i ,j ,k+1)*mask_up(i ,j ))
|
---|
967 | cff1=Ta(i-1,j ,k )*MAX(0.0_r8,Ua(i ,j ,k ))- &
|
---|
968 | & Ta(i+1,j ,k )*MIN(0.0_r8,Ua(i+1,j ,k ))+ &
|
---|
969 | & Ta(i ,j-1,k )*MAX(0.0_r8,Va(i ,j ,k ))- &
|
---|
970 | & Ta(i ,j+1,k )*MIN(0.0_r8,Va(i ,j+1,k ))+ &
|
---|
971 | & Ta(i ,j ,k-1)*MAX(0.0_r8,Wa(i ,j ,k-1))- &
|
---|
972 | & Ta(i ,j ,k+1)*MIN(0.0_r8,Wa(i ,j ,k ))
|
---|
973 | beta_up(i,j,k)=(Tmax-Ta(i,j,k))/(cff1+eps)
|
---|
974 | !
|
---|
975 | Tmin=MIN(Ta(i-1,j ,k )*mask_dn(i-1,j ), &
|
---|
976 | & t (i-1,j ,k )*mask_dn(i-1,j ), &
|
---|
977 | & Ta(i ,j ,k )*mask_dn(i ,j ), &
|
---|
978 | & t (i ,j ,k )*mask_dn(i ,j ), &
|
---|
979 | & Ta(i+1,j ,k )*mask_dn(i+1,j ), &
|
---|
980 | & t (i+1,j ,k )*mask_dn(i+1,j ), &
|
---|
981 | & Ta(i ,j-1,k )*mask_dn(i ,j-1), &
|
---|
982 | & t (i ,j-1,k )*mask_dn(i ,j-1), &
|
---|
983 | & Ta(i ,j+1,k )*mask_dn(i ,j+1), &
|
---|
984 | & t (i ,j+1,k )*mask_dn(i ,j+1), &
|
---|
985 | & Ta(i ,j ,k-1)*mask_dn(i ,j ), &
|
---|
986 | & t (i ,j ,k-1)*mask_dn(i ,j ), &
|
---|
987 | & Ta(i ,j ,k+1)*mask_dn(i ,j ), &
|
---|
988 | & t (i ,j ,k+1)*mask_dn(i ,j ))
|
---|
989 | cff2=Ta(i ,j ,k )*MAX(0.0_r8,Ua(i+1,j ,k ))- &
|
---|
990 | & Ta(i ,j ,k )*MIN(0.0_r8,Ua(i ,j ,k ))+ &
|
---|
991 | & Ta(i ,j ,k )*MAX(0.0_r8,Va(i ,j+1,k ))- &
|
---|
992 | & Ta(i ,j ,k )*MIN(0.0_r8,Va(i ,j ,k ))+ &
|
---|
993 | & Ta(i ,j ,k )*MAX(0.0_r8,Wa(i ,j ,k ))- &
|
---|
994 | & Ta(i ,j ,k )*MIN(0.0_r8,Wa(i ,j ,k-1))
|
---|
995 | beta_dn(i,j,k)=(Ta(i,j,k)-Tmin)/(cff2+eps)
|
---|
996 | END DO
|
---|
997 | END DO
|
---|
998 | k=N(ng)
|
---|
999 | DO i=IstrU-1,Iendp1
|
---|
1000 | Tmax=MAX(Ta(i-1,j ,k )*mask_up(i-1,j ), &
|
---|
1001 | & t (i-1,j ,k )*mask_up(i-1,j ), &
|
---|
1002 | & Ta(i ,j ,k )*mask_up(i ,j ), &
|
---|
1003 | & t (i ,j ,k )*mask_up(i ,j ), &
|
---|
1004 | & Ta(i+1,j ,k )*mask_up(i+1,j ), &
|
---|
1005 | & t (i+1,j ,k )*mask_up(i+1,j ), &
|
---|
1006 | & Ta(i ,j-1,k )*mask_up(i ,j-1), &
|
---|
1007 | & t (i ,j-1,k )*mask_up(i ,j-1), &
|
---|
1008 | & Ta(i ,j+1,k )*mask_up(i ,j+1), &
|
---|
1009 | & t (i ,j+1,k )*mask_up(i ,j+1), &
|
---|
1010 | & Ta(i ,j ,k-1)*mask_up(i ,j ), &
|
---|
1011 | & t (i ,j ,k-1)*mask_up(i ,j ))
|
---|
1012 | cff1=Ta(i-1,j ,k )*MAX(0.0_r8,Ua(i ,j ,k ))- &
|
---|
1013 | & Ta(i+1,j ,k )*MIN(0.0_r8,Ua(i+1,j ,k ))+ &
|
---|
1014 | & Ta(i ,j-1,k )*MAX(0.0_r8,Va(i ,j ,k ))- &
|
---|
1015 | & Ta(i ,j+1,k )*MIN(0.0_r8,Va(i ,j+1,k ))+ &
|
---|
1016 | & Ta(i ,j ,k-1)*MAX(0.0_r8,Wa(i ,j ,k-1))
|
---|
1017 | beta_up(i,j,k)=(Tmax-Ta(i,j,k))/(cff1+eps)
|
---|
1018 | !
|
---|
1019 | Tmin=MIN(Ta(i-1,j ,k )*mask_dn(i-1,j ), &
|
---|
1020 | & t (i-1,j ,k )*mask_dn(i-1,j ), &
|
---|
1021 | & Ta(i ,j ,k )*mask_dn(i ,j ), &
|
---|
1022 | & t (i ,j ,k )*mask_dn(i ,j ), &
|
---|
1023 | & Ta(i+1,j ,k )*mask_dn(i+1,j ), &
|
---|
1024 | & t (i+1,j ,k )*mask_dn(i+1,j ), &
|
---|
1025 | & Ta(i ,j-1,k )*mask_dn(i ,j-1), &
|
---|
1026 | & t (i ,j-1,k )*mask_dn(i ,j-1), &
|
---|
1027 | & Ta(i ,j+1,k )*mask_dn(i ,j+1), &
|
---|
1028 | & t (i ,j+1,k )*mask_dn(i ,j+1), &
|
---|
1029 | & Ta(i ,j ,k-1)*mask_dn(i ,j ), &
|
---|
1030 | & t (i ,j ,k-1)*mask_dn(i ,j ))
|
---|
1031 | cff2=Ta(i ,j ,k )*MAX(0.0_r8,Ua(i+1,j ,k ))- &
|
---|
1032 | & Ta(i ,j ,k )*MIN(0.0_r8,Ua(i ,j ,k ))+ &
|
---|
1033 | & Ta(i ,j ,k )*MAX(0.0_r8,Va(i ,j+1,k ))- &
|
---|
1034 | & Ta(i ,j ,k )*MIN(0.0_r8,Va(i ,j ,k ))- &
|
---|
1035 | & Ta(i ,j ,k )*MIN(0.0_r8,Wa(i ,j ,k-1))
|
---|
1036 | beta_dn(i,j,k)=(Ta(i,j,k)-Tmin)/(cff2+eps)
|
---|
1037 | END DO
|
---|
1038 | END DO
|
---|
1039 | DO k=1,N(ng)
|
---|
1040 | DO j=JstrV-1,Jendp1
|
---|
1041 | DO i=IstrU-1,Iendp1
|
---|
1042 | IF (mask_up(i,j).eq.0.0_r8) THEN
|
---|
1043 | beta_up(i,j,k)=2.0_r8
|
---|
1044 | beta_dn(i,j,k)=2.0_r8
|
---|
1045 | END IF
|
---|
1046 | END DO
|
---|
1047 | END DO
|
---|
1048 | END DO
|
---|
1049 | !
|
---|
1050 | ! Calculate monotonic velocities. Scale back to dimensional units.
|
---|
1051 | !
|
---|
1052 | DO k=1,N(ng)
|
---|
1053 | DO j=Jstr,Jend
|
---|
1054 | DO i=IstrU,Iendp1
|
---|
1055 | cff1=MIN(beta_dn(i-1,j,k),beta_up(i,j,k),1.0_r8)
|
---|
1056 | cff2=MIN(beta_up(i-1,j,k),beta_dn(i,j,k),1.0_r8)
|
---|
1057 | Ua(i,j,k)=(cff1*MAX(0.0_r8,Ua(i,j,k))+ &
|
---|
1058 | & cff2*MIN(0.0_r8,Ua(i,j,k)))* &
|
---|
1059 | & cff*om_u(i,j)
|
---|
1060 | # ifdef MASKING
|
---|
1061 | Ua(i,j,k)=Ua(i,j,k)*umask(i,j)
|
---|
1062 | # endif
|
---|
1063 | # ifdef WET_DRY
|
---|
1064 | Ua(i,j,k)=Ua(i,j,k)*umask_wet(i,j)
|
---|
1065 | # endif
|
---|
1066 | END DO
|
---|
1067 | END DO
|
---|
1068 | DO j=JstrV,Jendp1
|
---|
1069 | DO i=Istr,Iend
|
---|
1070 | cff1=MIN(beta_dn(i,j-1,k),beta_up(i,j,k),1.0_r8)
|
---|
1071 | cff2=MIN(beta_up(i,j-1,k),beta_dn(i,j,k),1.0_r8)
|
---|
1072 | Va(i,j,k)=(cff1*MAX(0.0_r8,Va(i,j,k))+ &
|
---|
1073 | & cff2*MIN(0.0_r8,Va(i,j,k)))* &
|
---|
1074 | & cff*on_v(i,j)
|
---|
1075 | # ifdef MASKING
|
---|
1076 | Va(i,j,k)=Va(i,j,k)*vmask(i,j)
|
---|
1077 | # endif
|
---|
1078 | # ifdef WET_DRY
|
---|
1079 | Va(i,j,k)=Va(i,j,k)*vmask_wet(i,j)
|
---|
1080 | # endif
|
---|
1081 | END DO
|
---|
1082 | END DO
|
---|
1083 | IF (k.lt.N(ng)) THEN
|
---|
1084 | DO j=Jstr,Jend
|
---|
1085 | DO i=Istr,Iend
|
---|
1086 | cff1=MIN(beta_dn(i,j,k),beta_up(i,j,k+1),1.0_r8)
|
---|
1087 | cff2=MIN(beta_up(i,j,k),beta_dn(i,j,k+1),1.0_r8)
|
---|
1088 | Wa(i,j,k)=(cff1*MAX(0.0_r8,Wa(i,j,k))+ &
|
---|
1089 | & cff2*MIN(0.0_r8,Wa(i,j,k)))* &
|
---|
1090 | & cff*omn(i,j)*(z_r(i,j,k+1)-z_r(i,j,k))
|
---|
1091 | # ifdef MASKING
|
---|
1092 | Wa(i,j,k)=Wa(i,j,k)*rmask(i,j)
|
---|
1093 | # endif
|
---|
1094 | # ifdef WET_DRY
|
---|
1095 | Wa(i,j,k)=Wa(i,j,k)*rmask_wet(i,j)
|
---|
1096 | # endif
|
---|
1097 | END DO
|
---|
1098 | END DO
|
---|
1099 | END IF
|
---|
1100 | END DO
|
---|
1101 | !
|
---|
1102 | ! Apply boundary conditions to anti-diffusive velocities.
|
---|
1103 | !
|
---|
1104 | IF (.not.EWperiodic(ng)) THEN
|
---|
1105 | IF (DOMAIN(ng)%Western_Edge(tile)) THEN
|
---|
1106 | IF (LBC(iwest,isBu3d,ng)%closed) THEN
|
---|
1107 | DO k=1,N(ng)
|
---|
1108 | DO j=Jstr,Jend
|
---|
1109 | Ua(Istr,j,k)=0.0_r8
|
---|
1110 | END DO
|
---|
1111 | END DO
|
---|
1112 | ELSE
|
---|
1113 | DO k=1,N(ng)
|
---|
1114 | DO j=Jstr,Jend
|
---|
1115 | Ua(Istr,j,k)=Ua(Istr+1,j,k)
|
---|
1116 | END DO
|
---|
1117 | END DO
|
---|
1118 | END IF
|
---|
1119 | END IF
|
---|
1120 |
|
---|
1121 | IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN
|
---|
1122 | IF (LBC(ieast,isBu3d,ng)%closed) THEN
|
---|
1123 | DO k=1,N(ng)
|
---|
1124 | DO j=Jstr,Jend
|
---|
1125 | Ua(Iend+1,j,k)=0.0_r8
|
---|
1126 | END DO
|
---|
1127 | END DO
|
---|
1128 | ELSE
|
---|
1129 | DO k=1,N(ng)
|
---|
1130 | DO j=Jstr,Jend
|
---|
1131 | Ua(Iend+1,j,k)=Ua(Iend,j,k)
|
---|
1132 | END DO
|
---|
1133 | END DO
|
---|
1134 | END IF
|
---|
1135 | END IF
|
---|
1136 | END IF
|
---|
1137 |
|
---|
1138 | IF (.not.NSperiodic(ng)) THEN
|
---|
1139 | IF (DOMAIN(ng)%Southern_Edge(tile)) THEN
|
---|
1140 | IF (LBC(isouth,isBv3d,ng)%closed) THEN
|
---|
1141 | DO k=1,N(ng)
|
---|
1142 | DO i=Istr,Iend
|
---|
1143 | Va(i,Jstr,k)=0.0_r8
|
---|
1144 | END DO
|
---|
1145 | END DO
|
---|
1146 | ELSE
|
---|
1147 | DO k=1,N(ng)
|
---|
1148 | DO i=Istr,Iend
|
---|
1149 | Va(i,Jstr,k)=Va(i,Jstr+1,k)
|
---|
1150 | END DO
|
---|
1151 | END DO
|
---|
1152 | END IF
|
---|
1153 | END IF
|
---|
1154 |
|
---|
1155 | IF (DOMAIN(ng)%Northern_Edge(tile)) THEN
|
---|
1156 | IF (LBC(inorth,isBv3d,ng)%closed) THEN
|
---|
1157 | DO k=1,N(ng)
|
---|
1158 | DO i=Istr,Iend
|
---|
1159 | Va(i,Jend+1,k)=0.0_r8
|
---|
1160 | END DO
|
---|
1161 | END DO
|
---|
1162 | ELSE
|
---|
1163 | DO k=1,N(ng)
|
---|
1164 | DO i=Istr,Iend
|
---|
1165 | Va(i,Jend+1,k)=Va(i,Jend,k)
|
---|
1166 | END DO
|
---|
1167 | END DO
|
---|
1168 | END IF
|
---|
1169 | END IF
|
---|
1170 | END IF
|
---|
1171 |
|
---|
1172 | RETURN
|
---|
1173 | END SUBROUTINE mpdata_adiff_tile
|
---|
1174 | #endif
|
---|
1175 | END MODULE mpdata_adiff_mod
|
---|