Line data Source code
1 : !-----------------------------------------------------------------------
2 : ! $Id$
3 : !===============================================================================
4 : module mean_adv
5 :
6 : ! Description:
7 : ! Module mean_adv computes the mean advection terms for all of the
8 : ! time-tendency (prognostic) equations in the CLUBB parameterization. All of
9 : ! the mean advection terms are solved for completely implicitly, and therefore
10 : ! become part of the left-hand side of their respective equations.
11 : !
12 : ! Function term_ma_zt_lhs handles the mean advection terms for the variables
13 : ! located at thermodynamic grid levels. These variables are: rtm, thlm, wp3,
14 : ! all hydrometeor species, and sclrm.
15 : !
16 : ! Function term_ma_zm_lhs handles the mean advection terms for the variables
17 : ! located at momentum grid levels. The variables are: wprtp, wpthlp, wp2,
18 : ! rtp2, thlp2, rtpthlp, up2, vp2, wpsclrp, sclrprtp, sclrpthlp, and sclrp2.
19 :
20 : implicit none
21 :
22 : private ! Default scope
23 :
24 : public :: term_ma_zt_lhs, &
25 : term_ma_zm_lhs
26 :
27 : integer, parameter :: &
28 : ndiags3 = 3
29 :
30 : contains
31 :
32 : !=============================================================================
33 8935056 : subroutine term_ma_zt_lhs( nz, ngrdcol, wm_zt, weights_zt2zm, & ! Intent(in)
34 8935056 : invrs_dzt, invrs_dzm, & ! Intent(in)
35 : l_upwind_xm_ma, & ! Intent(in)
36 8935056 : lhs_ma ) ! Intent(out)
37 :
38 : ! Description:
39 : ! Mean advection of var_zt: implicit portion of the code.
40 : !
41 : ! The variable "var_zt" stands for a variable that is located at
42 : ! thermodynamic grid levels.
43 : !
44 : ! The d(var_zt)/dt equation contains a mean advection term:
45 : !
46 : ! - w * d(var_zt)/dz.
47 : !
48 : ! This term is solved for completely implicitly, such that:
49 : !
50 : ! - w * d( var_zt(t+1) )/dz.
51 : !
52 : ! Note: When the term is brought over to the left-hand side, the sign
53 : ! is reversed and the leading "-" in front of the term is changed to
54 : ! a "+".
55 : !
56 : ! The timestep index (t+1) means that the value of var_zt being used is from
57 : ! the next timestep, which is being advanced to in solving the d(var_zt)/dt
58 : ! equation.
59 : !
60 : ! This term is discretized as follows:
61 : !
62 : ! The values of var_zt are found on the thermodynamic levels, as are the
63 : ! values of wm_zt (mean vertical velocity on thermodynamic levels). The
64 : ! variable var_zt is interpolated to the intermediate momentum levels. The
65 : ! derivative of the interpolated values is taken over the central
66 : ! thermodynamic level. The derivative is multiplied by wm_zt at the central
67 : ! thermodynamic level to get the desired result.
68 : !
69 : ! -----var_zt---------------------------------------------- t(k+1)
70 : !
71 : ! =================var_zt(interp)========================== m(k)
72 : !
73 : ! -----var_zt---------------------d(var_zt)/dz-----wm_zt--- t(k)
74 : !
75 : ! =================var_zt(interp)========================== m(k-1)
76 : !
77 : ! -----var_zt---------------------------------------------- t(k-1)
78 : !
79 : ! The vertical indices t(k+1), m(k), t(k), m(k-1), and t(k-1) correspond
80 : ! with altitudes zt(k+1), zm(k), zt(k), zm(k-1), and zt(k-1), respectively.
81 : ! The letter "t" is used for thermodynamic levels and the letter "m" is used
82 : ! for momentum levels.
83 : !
84 : ! invrs_dzt(k) = 1 / ( zm(k) - zm(k-1) )
85 : !
86 : !
87 : ! Special discretization for upper boundary level:
88 : !
89 : ! Method 1: Constant derivative method (or "one-sided" method).
90 : !
91 : ! The values of var_zt are found on the thermodynamic levels, as are the
92 : ! values of wm_zt (mean vertical velocity on the thermodynamic levels). The
93 : ! variable var_zt is interpolated to momentum level gr%nz-1, based on
94 : ! the values of var_zt at thermodynamic levels gr%nz and gr%nz-1.
95 : ! However, the variable var_zt cannot be interpolated to momentum level
96 : ! gr%nz. Rather, a linear extension is used to find the value of var_zt
97 : ! at momentum level gr%nz, based on the values of var_zt at thermodynamic
98 : ! levels gr%nz and gr%nz-1. The derivative of the extended and
99 : ! interpolated values, d(var_zt)/dz, is taken over the central thermodynamic
100 : ! level. Of course, this derivative will be the same as the derivative of
101 : ! var_zt between thermodynamic levels gr%nz and gr%nz-1. The derivative
102 : ! is multiplied by wm_zt at the central thermodynamic level to get the
103 : ! desired result.
104 : !
105 : ! For the following diagram, k = gr%nz, which is the uppermost level of
106 : ! the model:
107 : !
108 : ! =================var_zt(extend)========================== m(k) Boundary
109 : !
110 : ! -----var_zt---------------------d(var_zt)/dz-----wm_zt--- t(k)
111 : !
112 : ! =================var_zt(interp)========================== m(k-1)
113 : !
114 : ! -----var_zt---------------------------------------------- t(k-1)
115 : !
116 : !
117 : ! Method 2: Zero derivative method:
118 : ! the derivative d(var_zt)/dz over the model top is set to 0.
119 : !
120 : ! This method corresponds with the "zero-flux" boundary condition option
121 : ! for eddy diffusion, where d(var_zt)/dz is set to 0 across the upper
122 : ! boundary.
123 : !
124 : ! In order to discretize the upper boundary condition, consider a new level
125 : ! outside the model (thermodynamic level gr%nz+1) just above the upper
126 : ! boundary level (thermodynamic level gr%nz). The value of var_zt at the
127 : ! level just outside the model is defined to be the same as the value of
128 : ! var_zt at thermodynamic level gr%nz. Therefore, the value of
129 : ! d(var_zt)/dz between the level just outside the model and the uppermost
130 : ! thermodynamic level is 0, staying consistent with the zero-flux boundary
131 : ! condition option for the eddy diffusion portion of the code. Therefore,
132 : ! the value of var_zt at momentum level gr%nz, which is the upper boundary
133 : ! of the model, would be the same as the value of var_zt at the uppermost
134 : ! thermodynamic level.
135 : !
136 : ! The values of var_zt are found on the thermodynamic levels, as are the
137 : ! values of wm_zt (mean vertical velocity on the thermodynamic levels). The
138 : ! variable var_zt is interpolated to momentum level gr%nz-1, based on
139 : ! the values of var_zt at thermodynamic levels gr%nz and gr%nz-1. The
140 : ! value of var_zt at momentum level gr%nz is set equal to the value of
141 : ! var_zt at thermodynamic level gr%nz, as described above. The derivative
142 : ! of the set and interpolated values, d(var_zt)/dz, is taken over the
143 : ! central thermodynamic level. The derivative is multiplied by wm_zt at the
144 : ! central thermodynamic level to get the desired result.
145 : !
146 : ! For the following diagram, k = gr%nz, which is the uppermost level of
147 : ! the model:
148 : !
149 : ! --[var_zt(k+1) = var_zt(k)]----(level outside model)----- t(k+1)
150 : !
151 : ! ==[var_zt(top) = var_zt(k)]===[d(var_zt)/dz|_(top) = 0]== m(k) Boundary
152 : !
153 : ! -----var_zt(k)------------------d(var_zt)/dz-----wm_zt--- t(k)
154 : !
155 : ! =================var_zt(interp)========================== m(k-1)
156 : !
157 : ! -----var_zt(k-1)----------------------------------------- t(k-1)
158 : !
159 : ! where (top) stands for the grid index of momentum level k = gr%nz, which
160 : ! is the upper boundary of the model.
161 : !
162 : ! This method of boundary discretization is also similar to the method
163 : ! currently employed at the lower boundary for most thermodynamic-level
164 : ! variables. Since thermodynamic level k = 1 is below the model bottom,
165 : ! mean advection is not applied. Thus, thermodynamic level k = 2 becomes
166 : ! the lower boundary level. Now, the mean advection term at thermodynamic
167 : ! level 2 takes into account var_zt from levels 1, 2, and 3. However, in
168 : ! most cases, the value of var_zt(1) is set equal to var_zt(2) after the
169 : ! matrix of equations has been solved. Therefore, the derivative,
170 : ! d(var_zt)/dz, over the model bottom (momentum level k = 1) becomes 0.
171 : ! Thus, the method of setting d(var_zt)/dz to 0 over the model top keeps
172 : ! the way the upper and lower boundaries are handled consistent with each
173 : ! other.
174 :
175 : ! References:
176 : ! None
177 : !-----------------------------------------------------------------------
178 :
179 : use grid_class, only: &
180 : grid ! Type
181 :
182 : use constants_clubb, only: &
183 : one, & ! Constant(s)
184 : zero
185 :
186 : use clubb_precision, only: &
187 : core_rknd ! Variable(s)
188 :
189 : implicit none
190 :
191 : ! Constant parameters
192 : integer, parameter :: &
193 : kp1_tdiag = 1, & ! Thermodynamic superdiagonal index.
194 : k_tdiag = 2, & ! Thermodynamic main diagonal index.
195 : km1_tdiag = 3 ! Thermodynamic subdiagonal index.
196 :
197 : integer, parameter :: &
198 : t_above = 1, & ! Index for upper thermodynamic level grid weight.
199 : t_below = 2 ! Index for lower thermodynamic level grid weight.
200 :
201 : ! -------------------------- Input Variables --------------------------
202 : integer, intent(in) :: &
203 : nz, &
204 : ngrdcol
205 :
206 : real( kind = core_rknd ), dimension(ngrdcol,nz), intent(in) :: &
207 : wm_zt, & ! wm_zt [m/s]
208 : invrs_dzt, & ! Inverse of grid spacing [1/m]
209 : invrs_dzm ! Inverse of grid spacing [1/m]
210 :
211 : real( kind = core_rknd ), dimension(ngrdcol,nz,t_above:t_below), intent(in) :: &
212 : weights_zt2zm
213 :
214 : logical, intent(in) :: &
215 : l_upwind_xm_ma ! This flag determines whether we want to use an upwind
216 : ! differencing approximation rather than a centered
217 : ! differencing for turbulent or mean advection terms.
218 : ! It affects rtm, thlm, sclrm, um and vm.
219 :
220 : ! -------------------------- Return Variable --------------------------
221 : real( kind = core_rknd ), dimension(ndiags3,ngrdcol,nz), intent(out) :: &
222 : lhs_ma ! Mean advection contributions to lhs [1/s]
223 :
224 : ! -------------------------- Local Variables --------------------------
225 :
226 : integer :: i, k, b ! Vertical level index
227 :
228 : !-------------------------- Begin Code --------------------------
229 :
230 : !$acc data copyin( wm_zt, invrs_dzt, invrs_dzm, weights_zt2zm ) &
231 : !$acc copyout( lhs_ma )
232 :
233 : ! Set lower boundary array to 0
234 : !$acc parallel loop gang vector collapse(2) default(present)
235 149194656 : do i = 1, ngrdcol
236 569973456 : do b = 1, ndiags3
237 561038400 : lhs_ma(b,i,1) = 0.0_core_rknd
238 : end do
239 : end do
240 : !$acc end parallel loop
241 :
242 :
243 8935056 : if ( .not. l_upwind_xm_ma ) then ! Use centered differencing
244 :
245 : ! Most of the interior model; normal conditions.
246 : !$acc parallel loop gang vector collapse(2) default(present)
247 0 : do k = 3, nz-1, 1
248 0 : do i = 1, ngrdcol
249 :
250 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
251 0 : lhs_ma(kp1_tdiag,i,k) = + wm_zt(i,k) * invrs_dzt(i,k) * weights_zt2zm(i,k,t_above)
252 :
253 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
254 : lhs_ma(k_tdiag,i,k) = + wm_zt(i,k) * invrs_dzt(i,k) * ( weights_zt2zm(i,k,t_below) &
255 0 : - weights_zt2zm(i,k-1,t_above) )
256 :
257 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
258 0 : lhs_ma(km1_tdiag,i,k) = - wm_zt(i,k) * invrs_dzt(i,k) * weights_zt2zm(i,k-1,t_below)
259 :
260 : end do
261 : end do ! k = 3, nz-1, 1
262 : !$acc end parallel loop
263 :
264 : ! Upper Boundary
265 :
266 : ! Special discretization for zero derivative method, where the
267 : ! derivative d(var_zt)/dz over the model top is set to 0, in order
268 : ! to stay consistent with the zero-flux boundary condition option
269 : ! in the eddy diffusion code.
270 : !$acc parallel loop gang vector default(present)
271 0 : do i = 1, ngrdcol
272 :
273 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
274 0 : lhs_ma(kp1_tdiag,i,nz) = zero
275 :
276 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
277 : lhs_ma(k_tdiag,i,nz) = + wm_zt(i,nz) * invrs_dzt(i,nz) &
278 0 : * ( one - weights_zt2zm(i,nz-1,t_above) )
279 :
280 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
281 : lhs_ma(km1_tdiag,i,nz) = - wm_zt(i,nz) * invrs_dzt(i,nz) &
282 0 : * weights_zt2zm(i,nz-1,t_below)
283 :
284 : end do
285 : !$acc end parallel loop
286 :
287 : ! Lower Boundary
288 :
289 : ! Special discretization for zero derivative method, where the
290 : ! derivative d(var_zt)/dz over the model bottom is set to 0.
291 : !$acc parallel loop gang vector default(present)
292 0 : do i = 1, ngrdcol
293 :
294 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
295 0 : lhs_ma(kp1_tdiag,i,2) = + wm_zt(i,2) * invrs_dzt(i,2) &
296 0 : * weights_zt2zm(i,2,t_above)
297 :
298 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
299 : lhs_ma(k_tdiag,i,2) = - wm_zt(i,2) * invrs_dzt(i,2) &
300 0 : * ( one - weights_zt2zm(i,2,t_below) )
301 :
302 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
303 0 : lhs_ma(km1_tdiag,i,2) = zero
304 :
305 : end do
306 : !$acc end parallel loop
307 :
308 : else ! l_upwind_xm_ma == .true.; use "upwind" differencing
309 :
310 : ! Most of the interior model; normal conditions.
311 : !$acc parallel loop gang vector collapse(2) default(present)
312 741609648 : do k = 3, nz-1, 1
313 12242896848 : do i = 1, ngrdcol
314 12233961792 : if ( wm_zt(i,k) >= zero ) then ! Mean wind is in upward direction
315 :
316 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
317 5107411116 : lhs_ma(kp1_tdiag,i,k) = zero
318 :
319 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
320 5107411116 : lhs_ma(k_tdiag,i,k) = + wm_zt(i,k) * invrs_dzm(i,k-1)
321 :
322 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
323 5107411116 : lhs_ma(km1_tdiag,i,k) = - wm_zt(i,k) * invrs_dzm(i,k-1)
324 :
325 : else ! wm_zt < 0; Mean wind is in downward direction
326 :
327 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
328 6393876084 : lhs_ma(kp1_tdiag,i,k) = + wm_zt(i,k) * invrs_dzm(i,k)
329 :
330 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
331 6393876084 : lhs_ma(k_tdiag,i,k) = - wm_zt(i,k) * invrs_dzm(i,k)
332 :
333 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
334 6393876084 : lhs_ma(km1_tdiag,i,k) = zero
335 :
336 : endif ! wm_zt > 0
337 :
338 : end do
339 : end do ! k = 3, nz-1, 1
340 : !$acc end parallel loop
341 :
342 : ! Upper Boundary
343 : !$acc parallel loop gang vector default(present)
344 149194656 : do i = 1, ngrdcol
345 149194656 : if ( wm_zt(i,nz) >= zero ) then ! Mean wind is in upward direction
346 :
347 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
348 65343546 : lhs_ma(kp1_tdiag,i,nz) = zero
349 :
350 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
351 65343546 : lhs_ma(k_tdiag,i,nz) = + wm_zt(i,nz) * invrs_dzm(i,nz-1)
352 :
353 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
354 65343546 : lhs_ma(km1_tdiag,i,nz) = - wm_zt(i,nz) * invrs_dzm(i,nz-1)
355 :
356 : else ! wm_zt < 0; Mean wind is in downward direction
357 :
358 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
359 74916054 : lhs_ma(kp1_tdiag,i,nz) = zero
360 :
361 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
362 74916054 : lhs_ma(k_tdiag,i,nz) = zero
363 :
364 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
365 74916054 : lhs_ma(km1_tdiag,i,nz) = zero
366 :
367 : end if ! wm_zt > 0
368 : end do
369 : !$acc end parallel loop
370 :
371 : ! Lower Boundary
372 : !$acc parallel loop gang vector default(present)
373 149194656 : do i = 1, ngrdcol
374 149194656 : if ( wm_zt(i,2) >= zero ) then ! Mean wind is in upward direction
375 :
376 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
377 140259600 : lhs_ma(kp1_tdiag,i,2) = zero
378 :
379 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
380 140259600 : lhs_ma(k_tdiag,i,2) = zero
381 :
382 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
383 140259600 : lhs_ma(km1_tdiag,i,2) = zero
384 :
385 : else ! wm_zt < 0; Mean wind is in downward direction
386 :
387 : ! Thermodynamic superdiagonal: [ x var_zt(k+1,<t+1>) ]
388 0 : lhs_ma(kp1_tdiag,i,2) = + wm_zt(i,2) * invrs_dzm(i,2)
389 :
390 : ! Thermodynamic main diagonal: [ x var_zt(k,<t+1>) ]
391 0 : lhs_ma(k_tdiag,i,2) = - wm_zt(i,2) * invrs_dzm(i,2)
392 :
393 : ! Thermodynamic subdiagonal: [ x var_zt(k-1,<t+1>) ]
394 0 : lhs_ma(km1_tdiag,i,2) = zero
395 :
396 : endif ! wm_zt > 0
397 : end do
398 : !$acc end parallel loop
399 :
400 : endif ! l_upwind_xm_ma
401 :
402 : !$acc end data
403 :
404 8935056 : return
405 :
406 : end subroutine term_ma_zt_lhs
407 :
408 : !=============================================================================
409 26805168 : subroutine term_ma_zm_lhs( nz, ngrdcol, wm_zm, &
410 26805168 : invrs_dzm, weights_zm2zt, &
411 26805168 : lhs_ma )
412 :
413 : ! Description:
414 : ! Mean advection of var_zm: implicit portion of the code.
415 : !
416 : ! The variable "var_zm" stands for a variable that is located at momentum
417 : ! grid levels.
418 : !
419 : ! The d(var_zm)/dt equation contains a mean advection term:
420 : !
421 : ! - w * d(var_zm)/dz.
422 : !
423 : ! This term is solved for completely implicitly, such that:
424 : !
425 : ! - w * d( var_zm(t+1) )/dz.
426 : !
427 : ! Note: When the term is brought over to the left-hand side, the sign
428 : ! is reversed and the leading "-" in front of the term is changed to
429 : ! a "+".
430 : !
431 : ! The timestep index (t+1) means that the value of var_zm being used is from
432 : ! the next timestep, which is being advanced to in solving the d(var_zm)/dt
433 : ! equation.
434 : !
435 : ! This term is discretized as follows:
436 : !
437 : ! The values of var_zm are found on the momentum levels, as are the values
438 : ! of wm_zm (mean vertical velocity on momentum levels). The variable var_zm
439 : ! is interpolated to the intermediate thermodynamic levels. The derivative
440 : ! of the interpolated values is taken over the central momentum level. The
441 : ! derivative is multiplied by wm_zm at the central momentum level to get the
442 : ! desired result.
443 : !
444 : ! =====var_zm============================================== m(k+1)
445 : !
446 : ! -----------------var_zm(interp)-------------------------- t(k+1)
447 : !
448 : ! =====var_zm=====================d(var_zm)/dz=====wm_zm=== m(k)
449 : !
450 : ! -----------------var_zm(interp)-------------------------- t(k)
451 : !
452 : ! =====var_zm============================================== m(k-1)
453 : !
454 : ! The vertical indices m(k+1), t(k+1), m(k), t(k), and m(k-1) correspond
455 : ! with altitudes zm(k+1), zt(k+1), zm(k), zt(k), and zm(k-1), respectively.
456 : ! The letter "t" is used for thermodynamic levels and the letter "m" is used
457 : ! for momentum levels.
458 : !
459 : ! invrs_dzm(k) = 1 / ( zt(k+1) - zt(k) )
460 :
461 : ! References:
462 : !-----------------------------------------------------------------------
463 :
464 : use grid_class, only: &
465 : grid ! Type
466 :
467 : use constants_clubb, only: &
468 : zero ! Constant(s)
469 :
470 : use clubb_precision, only: &
471 : core_rknd ! Variable(s)
472 :
473 : implicit none
474 :
475 : ! Constant parameters
476 : integer, parameter :: &
477 : kp1_mdiag = 1, & ! Momentum superdiagonal index.
478 : k_mdiag = 2, & ! Momentum main diagonal index.
479 : km1_mdiag = 3 ! Momentum subdiagonal index.
480 :
481 : integer, parameter :: &
482 : m_above = 1, & ! Index for upper momentum level grid weight.
483 : m_below = 2 ! Index for lower momentum level grid weight.
484 :
485 : ! -------------------------- Input Variables --------------------------
486 : integer, intent(in) :: &
487 : nz, &
488 : ngrdcol
489 :
490 : real( kind = core_rknd ), dimension(ngrdcol,nz), intent(in) :: &
491 : wm_zm, & ! wm_zm [m/s]
492 : invrs_dzm ! Inverse of grid spacing [1/m]
493 :
494 : real( kind = core_rknd ), dimension(ngrdcol,nz,m_above:m_below), intent(in) :: &
495 : weights_zm2zt
496 :
497 : ! -------------------------- Return Variable --------------------------
498 : real( kind = core_rknd ), dimension(ndiags3,ngrdcol,nz), intent(out) :: &
499 : lhs_ma ! Mean advection contributions to lhs [1/s]
500 :
501 : ! -------------------------- Local Variables
502 : integer :: i, k, b ! Vertical level index
503 :
504 : ! -------------------------- Begin Code --------------------------
505 :
506 : !$acc data copyin( wm_zm, invrs_dzm, weights_zm2zt ) &
507 : !$acc copyout( lhs_ma )
508 :
509 : ! Set lower boundary array to 0
510 : !$acc parallel loop gang vector collapse(2) default(present)
511 447583968 : do i = 1, ngrdcol
512 1709920368 : do b = 1, ndiags3
513 1683115200 : lhs_ma(b,i,1) = zero
514 : end do
515 : end do
516 : !$acc end parallel loop
517 :
518 : ! Most of the interior model; normal conditions.
519 : !$acc parallel loop gang vector collapse(2) default(present)
520 2251634112 : do k = 2, nz-1, 1
521 37176274512 : do i = 1, ngrdcol
522 :
523 : ! Momentum superdiagonal: [ x var_zm(k+1,<t+1>) ]
524 34924640400 : lhs_ma(kp1_mdiag,i,k) = + wm_zm(i,k) * invrs_dzm(i,k) * weights_zm2zt(i,k+1,m_above)
525 :
526 : ! Momentum main diagonal: [ x var_zm(k,<t+1>) ]
527 : lhs_ma(k_mdiag,i,k) = + wm_zm(i,k) * invrs_dzm(i,k) * ( weights_zm2zt(i,k+1,m_below) &
528 34924640400 : - weights_zm2zt(i,k,m_above) )
529 :
530 : ! Momentum subdiagonal: [ x var_zm(k-1,<t+1>) ]
531 37149469344 : lhs_ma(km1_mdiag,i,k) = - wm_zm(i,k) * invrs_dzm(i,k) * weights_zm2zt(i,k,m_below)
532 :
533 : end do
534 : end do ! k = 2, nz-1, 1
535 : !$acc end parallel loop
536 :
537 : ! Set upper boundary array to 0
538 : !$acc parallel loop gang vector collapse(2) default(present)
539 447583968 : do i = 1, ngrdcol
540 1709920368 : do b = 1, ndiags3
541 1683115200 : lhs_ma(b,i,nz) = zero
542 : end do
543 : end do
544 : !$acc end parallel loop
545 :
546 : !$acc end data
547 :
548 26805168 : return
549 :
550 : end subroutine term_ma_zm_lhs
551 :
552 : !=============================================================================
553 :
554 : end module mean_adv
|