SCALE-RM
scale_atmos_dyn_fvm_flux_cd8.F90
Go to the documentation of this file.
1 
2 !-------------------------------------------------------------------------------
11 !-------------------------------------------------------------------------------
12 ! Warning: This file was generated from atmos-rm/dynamics/scale_atmos_dyn_fvm_flux_udcd.F90.erb.
13 ! Do not edit this file.
14 !-------------------------------------------------------------------------------
15 #include "scalelib.h"
17  !-----------------------------------------------------------------------------
18  !
19  !++ used modules
20  !
21  use scale_precision
22  use scale_io
23  use scale_prof
25  use scale_index
26  use scale_tracer
27  use scale_prc
28 #ifdef DEBUG
29  use scale_debug, only: &
30  check
31  use scale_const, only: &
32  undef => const_undef, &
33  iundef => const_undef2
34 #endif
35  !-----------------------------------------------------------------------------
36  implicit none
37  private
38  !-----------------------------------------------------------------------------
39  !
40  !++ Public procedure
41  !
43 
47 
53 
59 
65 
66  !-----------------------------------------------------------------------------
67  !
68  !++ Public parameters & variables
69  !
70  !-----------------------------------------------------------------------------
71  !
72  !++ Private procedure
73  !
74 #if 1
75 #define F2H(k,p,q) (CDZ(k+p-1)*GSQRT(k+p-1,i,j)/(CDZ(k)*GSQRT(k,i,j)+CDZ(k+1)*GSQRT(k+1,i,j)))
76 #else
77 #define F2H(k,p,q) 0.5_RP
78 #endif
79  !-----------------------------------------------------------------------------
80  !
81  !++ Private parameters & variables
82  !
83 
84  real(RP), parameter :: F2 = 0.5_rp ! F2 is always used to calculate flux near boundary.
85 
86 
87  real(RP), parameter :: F41 = 7.0_rp/12.0_rp
88  real(RP), parameter :: F42 = -1.0_rp/12.0_rp
89 
90 
91  real(RP), parameter :: F61 = 37.0_rp/60.0_rp
92  real(RP), parameter :: F62 = -8.0_rp/60.0_rp
93  real(RP), parameter :: F63 = 1.0_rp/60.0_rp
94 
95 
96  real(RP), parameter :: F81 = 533.0_rp/840.0_rp
97  real(RP), parameter :: F82 = -139.0_rp/840.0_rp
98  real(RP), parameter :: F83 = 29.0_rp/840.0_rp
99  real(RP), parameter :: F84 = -3.0_rp/840.0_rp
100 
101 
102 
103 contains
104  !-----------------------------------------------------------------------------
106 !OCL SERIAL
107  subroutine atmos_dyn_fvm_flux_valuew_z_cd8( &
108  valW, &
109  mflx, val, GSQRT, &
110  CDZ )
111  implicit none
112 
113  real(rp), intent(out) :: valw (ka)
114  real(rp), intent(in) :: mflx (ka)
115  real(rp), intent(in) :: val (ka)
116  real(rp), intent(in) :: gsqrt(ka)
117  real(rp), intent(in) :: cdz (ka)
118 
119  integer :: k
120  !---------------------------------------------------------------------------
121 
122  do k = ks+3, ke-4
123 #ifdef DEBUG
124  call check( __line__, mflx(k) )
125 
126  call check( __line__, val(k) )
127  call check( __line__, val(k+1) )
128 
129  call check( __line__, val(k-1) )
130  call check( __line__, val(k+2) )
131 
132  call check( __line__, val(k-2) )
133  call check( __line__, val(k+3) )
134 
135  call check( __line__, val(k-3) )
136  call check( __line__, val(k+4) )
137 
138 #endif
139  valw(k) = f81 * ( val(k+1)+val(k) ) &
140  + f82 * ( val(k+2)+val(k-1) ) &
141  + f83 * ( val(k+3)+val(k-2) ) &
142  + f84 * ( val(k+4)+val(k-3) )
143  enddo
144 #ifdef DEBUG
145  k = iundef
146 #endif
147 
148 #ifdef DEBUG
149 
150  call check( __line__, mflx(ks) )
151  call check( __line__, val(ks ) )
152  call check( __line__, val(ks+1) )
153  call check( __line__, mflx(ke-1) )
154  call check( __line__, val(ke ) )
155  call check( __line__, val(ke-1) )
156 
157  call check( __line__, mflx(ks+1) )
158  call check( __line__, val(ks+2 ) )
159  call check( __line__, val(ks+3) )
160  call check( __line__, mflx(ke-2) )
161  call check( __line__, val(ke-2 ) )
162  call check( __line__, val(ke-3) )
163 
164  call check( __line__, mflx(ks+2) )
165  call check( __line__, val(ks+4 ) )
166  call check( __line__, val(ks+5) )
167  call check( __line__, mflx(ke-3) )
168  call check( __line__, val(ke-4 ) )
169  call check( __line__, val(ke-5) )
170 
171 #endif
172 
173  valw(ks) = f2 * ( val(ks+1)+val(ks) )
174  valw(ke-1) = f2 * ( val(ke)+val(ke-1) )
175 
176  valw(ks+1) = f41 * ( val(ks+2)+val(ks+1) ) &
177  + f42 * ( val(ks+3)+val(ks) )
178  valw(ke-2) = f41 * ( val(ke-1)+val(ke-2) ) &
179  + f42 * ( val(ke)+val(ke-3) )
180 
181  valw(ks+2) = f61 * ( val(ks+3)+val(ks+2) ) &
182  + f62 * ( val(ks+4)+val(ks+1) ) &
183  + f63 * ( val(ks+5)+val(ks) )
184  valw(ke-3) = f61 * ( val(ke-2)+val(ke-3) ) &
185  + f62 * ( val(ke-1)+val(ke-4) ) &
186  + f63 * ( val(ke)+val(ke-5) )
187 
188 
189  return
190  end subroutine atmos_dyn_fvm_flux_valuew_z_cd8
191 
192  !-----------------------------------------------------------------------------
194  subroutine atmos_dyn_fvm_fluxz_xyz_cd8( &
195  flux, &
196  mflx, val, GSQRT, &
197  num_diff, &
198  CDZ, &
199  IIS, IIE, JJS, JJE )
200  use scale_const, only: &
201  eps => const_eps
202  implicit none
203 
204  real(rp), intent(inout) :: flux (ka,ia,ja)
205  real(rp), intent(in) :: mflx (ka,ia,ja)
206  real(rp), intent(in) :: val (ka,ia,ja)
207  real(rp), intent(in) :: gsqrt (ka,ia,ja)
208  real(rp), intent(in) :: num_diff(ka,ia,ja)
209  real(rp), intent(in) :: cdz (ka)
210  integer, intent(in) :: iis, iie, jjs, jje
211 
212  real(rp) :: vel
213  integer :: k, i, j
214  !---------------------------------------------------------------------------
215 
216  !$omp parallel default(none) private(i,j,k, vel) &
217  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mflx,val,flux,GSQRT,num_diff,EPS)
218 
219  !$omp do OMP_SCHEDULE_ collapse(2)
220  do j = jjs, jje
221  do i = iis, iie
222  do k = ks+3, ke-4
223 #ifdef DEBUG
224  call check( __line__, mflx(k,i,j) )
225 
226  call check( __line__, val(k,i,j) )
227  call check( __line__, val(k+1,i,j) )
228 
229  call check( __line__, val(k-1,i,j) )
230  call check( __line__, val(k+2,i,j) )
231 
232  call check( __line__, val(k-2,i,j) )
233  call check( __line__, val(k+3,i,j) )
234 
235  call check( __line__, val(k-3,i,j) )
236  call check( __line__, val(k+4,i,j) )
237 
238 #endif
239  vel = mflx(k,i,j)
240  flux(k,i,j) = vel &
241  * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
242  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
243  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
244  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) ) &
245  + gsqrt(k,i,j) * num_diff(k,i,j)
246  enddo
247  enddo
248  enddo
249  !$omp end do nowait
250 #ifdef DEBUG
251  k = iundef; i = iundef; j = iundef
252 #endif
253 
254  !$omp do OMP_SCHEDULE_ collapse(2)
255  do j = jjs, jje
256  do i = iis, iie
257 #ifdef DEBUG
258 
259  call check( __line__, mflx(ks,i,j) )
260  call check( __line__, val(ks ,i,j) )
261  call check( __line__, val(ks+1,i,j) )
262  call check( __line__, mflx(ke-1,i,j) )
263  call check( __line__, val(ke ,i,j) )
264  call check( __line__, val(ke-1,i,j) )
265 
266  call check( __line__, mflx(ks+1,i,j) )
267  call check( __line__, val(ks+2 ,i,j) )
268  call check( __line__, val(ks+3,i,j) )
269  call check( __line__, mflx(ke-2,i,j) )
270  call check( __line__, val(ke-2 ,i,j) )
271  call check( __line__, val(ke-3,i,j) )
272 
273  call check( __line__, mflx(ks+2,i,j) )
274  call check( __line__, val(ks+4 ,i,j) )
275  call check( __line__, val(ks+5,i,j) )
276  call check( __line__, mflx(ke-3,i,j) )
277  call check( __line__, val(ke-4 ,i,j) )
278  call check( __line__, val(ke-5,i,j) )
279 
280 #endif
281  flux(ks-1,i,j) = 0.0_rp
282 
283  vel = mflx(ks,i,j)
284  flux(ks,i,j) = vel &
285  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) &
286  + gsqrt(ks,i,j) * num_diff(ks,i,j)
287  vel = mflx(ke-1,i,j)
288  flux(ke-1,i,j) = vel &
289  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) ) &
290  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
291 
292  vel = mflx(ks+1,i,j)
293  flux(ks+1,i,j) = vel &
294  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
295  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) ) &
296  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j)
297  vel = mflx(ke-2,i,j)
298  flux(ke-2,i,j) = vel &
299  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
300  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) ) &
301  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j)
302 
303  vel = mflx(ks+2,i,j)
304  flux(ks+2,i,j) = vel &
305  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
306  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
307  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) ) &
308  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j)
309  vel = mflx(ke-3,i,j)
310  flux(ke-3,i,j) = vel &
311  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
312  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
313  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) ) &
314  + gsqrt(ke-3,i,j) * num_diff(ke-3,i,j)
315 
316  flux(ke ,i,j) = 0.0_rp
317  enddo
318  enddo
319  !$omp end do nowait
320 
321  !$omp end parallel
322 #ifdef DEBUG
323  k = iundef; i = iundef; j = iundef
324 #endif
325 
326  return
327  end subroutine atmos_dyn_fvm_fluxz_xyz_cd8
328 
329  !-----------------------------------------------------------------------------
331  subroutine atmos_dyn_fvm_fluxx_xyz_cd8( &
332  flux, &
333  mflx, val, GSQRT, &
334  num_diff, &
335  CDZ, &
336  IIS, IIE, JJS, JJE )
337  implicit none
338 
339  real(rp), intent(inout) :: flux (ka,ia,ja)
340  real(rp), intent(in) :: mflx (ka,ia,ja)
341  real(rp), intent(in) :: val (ka,ia,ja)
342  real(rp), intent(in) :: gsqrt (ka,ia,ja)
343  real(rp), intent(in) :: num_diff(ka,ia,ja)
344  real(rp), intent(in) :: cdz(ka)
345  integer, intent(in) :: iis, iie, jjs, jje
346 
347  real(rp) :: vel
348  integer :: k, i, j
349  !---------------------------------------------------------------------------
350 
351  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
352  !$omp private(vel) &
353  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mflx,val,flux,GSQRT,num_diff)
354  do j = jjs, jje
355  do i = iis-1, iie
356  do k = ks, ke
357 #ifdef DEBUG
358  call check( __line__, mflx(k,i,j) )
359 
360  call check( __line__, val(k,i,j) )
361  call check( __line__, val(k,i+1,j) )
362 
363  call check( __line__, val(k,i-1,j) )
364  call check( __line__, val(k,i+2,j) )
365 
366  call check( __line__, val(k,i-2,j) )
367  call check( __line__, val(k,i+3,j) )
368 
369  call check( __line__, val(k,i-3,j) )
370  call check( __line__, val(k,i+4,j) )
371 
372 #endif
373  vel = mflx(k,i,j)
374  flux(k,i,j) = vel &
375  * ( f81 * ( val(k,i+1,j)+val(k,i,j) ) &
376  + f82 * ( val(k,i+2,j)+val(k,i-1,j) ) &
377  + f83 * ( val(k,i+3,j)+val(k,i-2,j) ) &
378  + f84 * ( val(k,i+4,j)+val(k,i-3,j) ) ) &
379  + gsqrt(k,i,j) * num_diff(k,i,j)
380  enddo
381  enddo
382  enddo
383 #ifdef DEBUG
384  k = iundef; i = iundef; j = iundef
385 #endif
386 
387  return
388  end subroutine atmos_dyn_fvm_fluxx_xyz_cd8
389 
390  !-----------------------------------------------------------------------------
392  subroutine atmos_dyn_fvm_fluxy_xyz_cd8( &
393  flux, &
394  mflx, val, GSQRT, &
395  num_diff, &
396  CDZ, &
397  IIS, IIE, JJS, JJE )
398  implicit none
399 
400  real(rp), intent(inout) :: flux (ka,ia,ja)
401  real(rp), intent(in) :: mflx (ka,ia,ja)
402  real(rp), intent(in) :: val (ka,ia,ja)
403  real(rp), intent(in) :: gsqrt (ka,ia,ja)
404  real(rp), intent(in) :: num_diff(ka,ia,ja)
405  real(rp), intent(in) :: cdz(ka)
406  integer, intent(in) :: iis, iie, jjs, jje
407 
408  real(rp) :: vel
409  integer :: k, i, j
410  !---------------------------------------------------------------------------
411 
412  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
413  !$omp private(vel) &
414  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mflx,val,flux,GSQRT,num_diff)
415  do j = jjs-1, jje
416  do i = iis, iie
417  do k = ks, ke
418 #ifdef DEBUG
419  call check( __line__, mflx(k,i,j) )
420 
421  call check( __line__, val(k,i,j) )
422  call check( __line__, val(k,i,j+1) )
423 
424  call check( __line__, val(k,i,j-1) )
425  call check( __line__, val(k,i,j+2) )
426 
427  call check( __line__, val(k,i,j-2) )
428  call check( __line__, val(k,i,j+3) )
429 
430  call check( __line__, val(k,i,j-3) )
431  call check( __line__, val(k,i,j+4) )
432 
433 #endif
434  vel = mflx(k,i,j)
435  flux(k,i,j) = vel &
436  * ( f81 * ( val(k,i,j+1)+val(k,i,j) ) &
437  + f82 * ( val(k,i,j+2)+val(k,i,j-1) ) &
438  + f83 * ( val(k,i,j+3)+val(k,i,j-2) ) &
439  + f84 * ( val(k,i,j+4)+val(k,i,j-3) ) ) &
440  + gsqrt(k,i,j) * num_diff(k,i,j)
441  enddo
442  enddo
443  enddo
444 #ifdef DEBUG
445  k = iundef; i = iundef; j = iundef
446 #endif
447 
448  return
449  end subroutine atmos_dyn_fvm_fluxy_xyz_cd8
450 
451 
452  !-----------------------------------------------------------------------------
454  subroutine atmos_dyn_fvm_fluxz_xyw_cd8( &
455  flux, &
456  mom, val, DENS, &
457  GSQRT, J33G, &
458  num_diff, &
459  CDZ, FDZ, &
460  dtrk, &
461  IIS, IIE, JJS, JJE )
462  implicit none
463 
464  real(rp), intent(inout) :: flux (ka,ia,ja)
465  real(rp), intent(in) :: mom (ka,ia,ja)
466  real(rp), intent(in) :: val (ka,ia,ja)
467  real(rp), intent(in) :: dens (ka,ia,ja)
468  real(rp), intent(in) :: gsqrt (ka,ia,ja)
469  real(rp), intent(in) :: j33g
470  real(rp), intent(in) :: num_diff(ka,ia,ja)
471  real(rp), intent(in) :: cdz (ka)
472  real(rp), intent(in) :: fdz (ka-1)
473  real(rp), intent(in) :: dtrk
474  integer, intent(in) :: iis, iie, jjs, jje
475 
476  real(rp) :: vel
477  integer :: k, i, j
478  !---------------------------------------------------------------------------
479 
480  ! note than z-index is added by -1
481 
482  !$omp parallel default(none) private(i,j,k,vel) &
483  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,flux,J33G,GSQRT,num_diff,DENS,FDZ,dtrk)
484 
485  !$omp do OMP_SCHEDULE_ collapse(2)
486  do j = jjs, jje
487  do i = iis, iie
488  do k = ks+4, ke-3
489 #ifdef DEBUG
490  call check( __line__, mom(k-1,i,j) )
491  call check( __line__, mom(k ,i,j) )
492 
493  call check( __line__, val(k-1,i,j) )
494  call check( __line__, val(k,i,j) )
495 
496  call check( __line__, val(k-2,i,j) )
497  call check( __line__, val(k+1,i,j) )
498 
499  call check( __line__, val(k-3,i,j) )
500  call check( __line__, val(k+2,i,j) )
501 
502  call check( __line__, val(k-4,i,j) )
503  call check( __line__, val(k+3,i,j) )
504 
505 #endif
506  vel = ( 0.5_rp * ( mom(k-1,i,j) &
507  + mom(k,i,j) ) ) &
508  / dens(k,i,j)
509  flux(k-1,i,j) = j33g * vel &
510  * ( f81 * ( val(k,i,j)+val(k-1,i,j) ) &
511  + f82 * ( val(k+1,i,j)+val(k-2,i,j) ) &
512  + f83 * ( val(k+2,i,j)+val(k-3,i,j) ) &
513  + f84 * ( val(k+3,i,j)+val(k-4,i,j) ) ) &
514  + gsqrt(k,i,j) * num_diff(k,i,j)
515  enddo
516  enddo
517  enddo
518  !$omp end do nowait
519 #ifdef DEBUG
520  k = iundef; i = iundef; j = iundef
521 #endif
522 
523  !$omp do OMP_SCHEDULE_ collapse(2)
524  do j = jjs, jje
525  do i = iis, iie
526 #ifdef DEBUG
527 
528  call check( __line__, val(ks,i,j) )
529  call check( __line__, val(ks+1,i,j) )
530  call check( __line__, val(ks+2,i,j) )
531  call check( __line__, val(ks+3,i,j) )
532  call check( __line__, val(ks+4,i,j) )
533  call check( __line__, val(ks+5,i,j) )
534 
535 
536  call check( __line__, val(ke-5,i,j) )
537  call check( __line__, val(ke-4,i,j) )
538  call check( __line__, val(ke-3,i,j) )
539  call check( __line__, val(ke-2,i,j) )
540  call check( __line__, val(ke-1,i,j) )
541 
542 #endif
543  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS.
544  ! The flux at KS can be non-zero.
545  ! To reduce calculations, all the fluxes are set to zero.
546  flux(ks-1,i,j) = 0.0_rp ! k = KS
547 
548  vel = ( 0.5_rp * ( mom(ks,i,j) &
549  + mom(ks+1,i,j) ) ) &
550  / dens(ks+1,i,j)
551  flux(ks,i,j) = j33g * vel &
552  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) &
553  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j) ! k = KS+1
554 
555  vel = ( 0.5_rp * ( mom(ks+1,i,j) &
556  + mom(ks+2,i,j) ) ) &
557  / dens(ks+2,i,j)
558  flux(ks+1,i,j) = j33g * vel &
559  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
560  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) ) &
561  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j) ! k = KS+2
562 
563  vel = ( 0.5_rp * ( mom(ks+2,i,j) &
564  + mom(ks+3,i,j) ) ) &
565  / dens(ks+3,i,j)
566  flux(ks+2,i,j) = j33g * vel &
567  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
568  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
569  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) ) &
570  + gsqrt(ks+3,i,j) * num_diff(ks+3,i,j) ! k = KS+3
571 
572 
573 
574  vel = ( 0.5_rp * ( mom(ke-3,i,j) &
575  + mom(ke-2,i,j) ) ) &
576  / dens(ke-1,i,j)
577  flux(ke-3,i,j) = j33g * vel &
578  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
579  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
580  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) ) &
581  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j) ! k = KE-2
582 
583  vel = ( 0.5_rp * ( mom(ke-2,i,j) &
584  + mom(ke-1,i,j) ) ) &
585  / dens(ke-1,i,j)
586  flux(ke-2,i,j) = j33g * vel &
587  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
588  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) ) &
589  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j) ! k = KE-1
590 
591  flux(ke-1,i,j) = 0.0_rp ! k = KE
592  flux(ke ,i,j) = 0.0_rp ! k = KE+1
593  enddo
594  enddo
595  !$omp end do nowait
596 
597  !$omp end parallel
598 
599  return
600  end subroutine atmos_dyn_fvm_fluxz_xyw_cd8
601 
602 
603  !-----------------------------------------------------------------------------
605  subroutine atmos_dyn_fvm_fluxj13_xyw_cd8( &
606  flux, &
607  mom, val, DENS, &
608  GSQRT, J13G, MAPF, &
609  CDZ, TwoD, &
610  IIS, IIE, JJS, JJE )
611  implicit none
612 
613  real(rp), intent(inout) :: flux (ka,ia,ja)
614  real(rp), intent(in) :: mom (ka,ia,ja)
615  real(rp), intent(in) :: val (ka,ia,ja)
616  real(rp), intent(in) :: dens (ka,ia,ja)
617  real(rp), intent(in) :: gsqrt (ka,ia,ja)
618  real(rp), intent(in) :: j13g (ka,ia,ja)
619  real(rp), intent(in) :: mapf ( ia,ja,2)
620  real(rp), intent(in) :: cdz (ka)
621  logical, intent(in) :: twod
622  integer, intent(in) :: iis, iie, jjs, jje
623 
624  real(rp) :: vel
625  integer :: k, i, j
626  !---------------------------------------------------------------------------
627 
628  !$omp parallel default(none) private(i,j,k,vel) &
629  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J13G,MAPF)
630 
631  !$omp do OMP_SCHEDULE_ collapse(2)
632  do j = jjs, jje
633  do i = iis, iie
634  do k = ks+3, ke-3
635  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i-1,j) ) ) &
636  / dens(k,i,j)
637  vel = vel * j13g(k,i,j)
638  flux(k-1,i,j) = vel / mapf(i,j,+2) &
639  * ( f81 * ( val(k,i,j)+val(k-1,i,j) ) &
640  + f82 * ( val(k+1,i,j)+val(k-2,i,j) ) &
641  + f83 * ( val(k+2,i,j)+val(k-3,i,j) ) &
642  + f84 * ( val(k+3,i,j)+val(k-4,i,j) ) )
643  enddo
644  enddo
645  enddo
646  !$omp end do nowait
647 
648  !$omp do OMP_SCHEDULE_ collapse(2)
649  do j = jjs, jje
650  do i = iis, iie
651  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS.
652  ! The flux at KS can be non-zero.
653  ! To reduce calculations, all the fluxes are set to zero.
654  flux(ks-1,i,j) = 0.0_rp ! k = KS
655 
656  ! physically incorrect but for numerical stability
657  vel = ( ( 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i-1,j) ) ) / dens(ks+1,i,j) &
658  + ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i-1,j) ) ) / dens(ks ,i,j) ) * 0.5_rp
659 ! vel = ( 0.5_RP * ( mom(KS+1,i,j)+mom(KS+1,i-1,j) ) ) &
660 ! / DENS(KS+1,i,j)
661  vel = vel * j13g(ks+1,i,j)
662  flux(ks,i,j) = vel / mapf(i,j,+2) &
663  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) ! k = KS+1
664 
665  vel = ( 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i-1,j) ) ) &
666  / dens(ks+2,i,j)
667  vel = vel * j13g(ks,i,j)
668  flux(ks+1,i,j) = vel / mapf(i,j,+2) &
669  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
670  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) ) ! k = KS+2
671 
672 
673  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i-1,j) ) ) &
674  / dens(ke-1,i,j)
675  vel = vel * j13g(ke-1,i,j)
676  flux(ke-2,i,j) = vel / mapf(i,j,+2) &
677  * ( f2 * ( val(ke-1,i,j)+val(ke-2,i,j) ) ) ! k = KE-3
678 
679  vel = ( 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i-1,j) ) ) &
680  / dens(ke-2,i,j)
681  vel = vel * j13g(ke-2,i,j)
682  flux(ke-3,i,j) = vel / mapf(i,j,+2) &
683  * ( f41 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
684  + f42 * ( val(ke-1,i,j)+val(ke-4,i,j) ) ) ! k = KE-4
685 
686  flux(ke-1,i,j) = 0.0_rp
687  enddo
688  enddo
689  !$omp end do nowait
690 
691  !$omp end parallel
692 
693  return
694  end subroutine atmos_dyn_fvm_fluxj13_xyw_cd8
695 
696  !-----------------------------------------------------------------------------
698  subroutine atmos_dyn_fvm_fluxj23_xyw_cd8( &
699  flux, &
700  mom, val, DENS, &
701  GSQRT, J23G, MAPF, &
702  CDZ, TwoD, &
703  IIS, IIE, JJS, JJE )
704  implicit none
705 
706  real(rp), intent(inout) :: flux (ka,ia,ja)
707  real(rp), intent(in) :: mom (ka,ia,ja)
708  real(rp), intent(in) :: val (ka,ia,ja)
709  real(rp), intent(in) :: dens (ka,ia,ja)
710  real(rp), intent(in) :: gsqrt (ka,ia,ja)
711  real(rp), intent(in) :: j23g (ka,ia,ja)
712  real(rp), intent(in) :: mapf ( ia,ja,2)
713  real(rp), intent(in) :: cdz (ka)
714  logical, intent(in) :: twod
715  integer, intent(in) :: iis, iie, jjs, jje
716 
717  real(rp) :: vel
718  integer :: k, i, j
719  !---------------------------------------------------------------------------
720 
721  !$omp parallel default(none) private(i,j,k,vel) &
722  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J23G,MAPF)
723 
724  !$omp do OMP_SCHEDULE_ collapse(2)
725  do j = jjs, jje
726  do i = iis, iie
727  do k = ks+3, ke-3
728  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j-1) ) ) &
729  / dens(k,i,j)
730  vel = vel * j23g(k,i,j)
731  flux(k-1,i,j) = vel / mapf(i,j,+1) &
732  * ( f81 * ( val(k,i,j)+val(k-1,i,j) ) &
733  + f82 * ( val(k+1,i,j)+val(k-2,i,j) ) &
734  + f83 * ( val(k+2,i,j)+val(k-3,i,j) ) &
735  + f84 * ( val(k+3,i,j)+val(k-4,i,j) ) )
736  enddo
737  enddo
738  enddo
739  !$omp end do nowait
740 
741  !$omp do OMP_SCHEDULE_ collapse(2)
742  do j = jjs, jje
743  do i = iis, iie
744  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS.
745  ! The flux at KS can be non-zero.
746  ! To reduce calculations, all the fluxes are set to zero.
747  flux(ks-1,i,j) = 0.0_rp ! k = KS
748 
749  ! physically incorrect but for numerical stability
750  vel = ( ( 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i,j-1) ) ) / dens(ks+1,i,j) &
751  + ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i,j-1) ) ) / dens(ks ,i,j) ) * 0.5_rp
752 ! vel = ( 0.5_RP * ( mom(KS+1,i,j)+mom(KS+1,i,j-1) ) ) &
753 ! / DENS(KS+1,i,j)
754  vel = vel * j23g(ks+1,i,j)
755  flux(ks,i,j) = vel / mapf(i,j,+1) &
756  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) ! k = KS+1
757 
758  vel = ( 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i,j-1) ) ) &
759  / dens(ks+2,i,j)
760  vel = vel * j23g(ks,i,j)
761  flux(ks+1,i,j) = vel / mapf(i,j,+1) &
762  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
763  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) ) ! k = KS+2
764 
765 
766  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j-1) ) ) &
767  / dens(ke-1,i,j)
768  vel = vel * j23g(ke-1,i,j)
769  flux(ke-2,i,j) = vel / mapf(i,j,+1) &
770  * ( f2 * ( val(ke-1,i,j)+val(ke-2,i,j) ) ) ! k = KE-3
771 
772  vel = ( 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i,j-1) ) ) &
773  / dens(ke-2,i,j)
774  vel = vel * j23g(ke-2,i,j)
775  flux(ke-3,i,j) = vel / mapf(i,j,+1) &
776  * ( f41 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
777  + f42 * ( val(ke-1,i,j)+val(ke-4,i,j) ) ) ! k = KE-4
778 
779  flux(ke-1,i,j) = 0.0_rp
780  enddo
781  enddo
782  !$omp end do nowait
783 
784  !$omp end parallel
785 
786  return
787  end subroutine atmos_dyn_fvm_fluxj23_xyw_cd8
788 
789 
790  !-----------------------------------------------------------------------------
792  subroutine atmos_dyn_fvm_fluxx_xyw_cd8( &
793  flux, &
794  mom, val, DENS, &
795  GSQRT, MAPF, &
796  num_diff, &
797  CDZ, TwoD, &
798  IIS, IIE, JJS, JJE )
799  implicit none
800 
801  real(rp), intent(inout) :: flux (ka,ia,ja)
802  real(rp), intent(in) :: mom (ka,ia,ja)
803  real(rp), intent(in) :: val (ka,ia,ja)
804  real(rp), intent(in) :: dens (ka,ia,ja)
805  real(rp), intent(in) :: gsqrt (ka,ia,ja)
806  real(rp), intent(in) :: mapf ( ia,ja,2)
807  real(rp), intent(in) :: num_diff(ka,ia,ja)
808  real(rp), intent(in) :: cdz (ka)
809  logical, intent(in) :: twod
810  integer, intent(in) :: iis, iie, jjs, jje
811 
812  real(rp) :: vel
813  integer :: k, i, j
814  !---------------------------------------------------------------------------
815 
816  !$omp parallel default(none) private(i,j,k,vel) &
817  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff) &
818  !$omp shared(CDZ)
819 
820  !$omp do OMP_SCHEDULE_ collapse(2)
821  do j = jjs, jje
822  do i = iis-1, iie
823  do k = ks, ke-1
824 #ifdef DEBUG
825  call check( __line__, mom(k ,i,j) )
826  call check( __line__, mom(k+1,i,j) )
827 
828  call check( __line__, val(k,i,j) )
829  call check( __line__, val(k,i+1,j) )
830 
831  call check( __line__, val(k,i-1,j) )
832  call check( __line__, val(k,i+2,j) )
833 
834  call check( __line__, val(k,i-2,j) )
835  call check( __line__, val(k,i+3,j) )
836 
837  call check( __line__, val(k,i-3,j) )
838  call check( __line__, val(k,i+4,j) )
839 
840 #endif
841  vel = ( f2h(k,1,i_uyz) &
842  * mom(k+1,i,j) &
843  + f2h(k,2,i_uyz) &
844  * mom(k,i,j) ) &
845  / ( f2h(k,1,i_uyz) &
846  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
847  + f2h(k,2,i_uyz) &
848  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
849  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
850  * ( f81 * ( val(k,i+1,j)+val(k,i,j) ) &
851  + f82 * ( val(k,i+2,j)+val(k,i-1,j) ) &
852  + f83 * ( val(k,i+3,j)+val(k,i-2,j) ) &
853  + f84 * ( val(k,i+4,j)+val(k,i-3,j) ) ) &
854  + gsqrt(k,i,j) * num_diff(k,i,j)
855  enddo
856  enddo
857  enddo
858  !$omp end do nowait
859 #ifdef DEBUG
860  k = iundef; i = iundef; j = iundef
861 #endif
862 
863  !$omp do OMP_SCHEDULE_ collapse(2)
864  do j = jjs, jje
865  do i = iis-1, iie
866  flux(ke,i,j) = 0.0_rp
867  enddo
868  enddo
869  !$omp end do nowait
870 
871  !$omp end parallel
872 #ifdef DEBUG
873  k = iundef; i = iundef; j = iundef
874 #endif
875 
876  return
877  end subroutine atmos_dyn_fvm_fluxx_xyw_cd8
878 
879  !-----------------------------------------------------------------------------
881  subroutine atmos_dyn_fvm_fluxy_xyw_cd8( &
882  flux, &
883  mom, val, DENS, &
884  GSQRT, MAPF, &
885  num_diff, &
886  CDZ, TwoD, &
887  IIS, IIE, JJS, JJE )
888  implicit none
889 
890  real(rp), intent(inout) :: flux (ka,ia,ja)
891  real(rp), intent(in) :: mom (ka,ia,ja)
892  real(rp), intent(in) :: val (ka,ia,ja)
893  real(rp), intent(in) :: dens (ka,ia,ja)
894  real(rp), intent(in) :: gsqrt (ka,ia,ja)
895  real(rp), intent(in) :: mapf ( ia,ja,2)
896  real(rp), intent(in) :: num_diff(ka,ia,ja)
897  real(rp), intent(in) :: cdz (ka)
898  logical, intent(in) :: twod
899  integer, intent(in) :: iis, iie, jjs, jje
900 
901  real(rp) :: vel
902  integer :: k, i, j
903  !---------------------------------------------------------------------------
904 
905  !$omp parallel default(none) private(i,j,k,vel) &
906  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff) &
907  !$omp shared(CDZ)
908 
909  !$omp do OMP_SCHEDULE_ collapse(2)
910  do j = jjs-1, jje
911  do i = iis, iie
912  do k = ks, ke-1
913 #ifdef DEBUG
914  call check( __line__, mom(k ,i,j) )
915  call check( __line__, mom(k+1,i,j) )
916 
917  call check( __line__, val(k,i,j) )
918  call check( __line__, val(k,i,j+1) )
919 
920  call check( __line__, val(k,i,j-1) )
921  call check( __line__, val(k,i,j+2) )
922 
923  call check( __line__, val(k,i,j-2) )
924  call check( __line__, val(k,i,j+3) )
925 
926  call check( __line__, val(k,i,j-3) )
927  call check( __line__, val(k,i,j+4) )
928 
929 #endif
930  vel = ( f2h(k,1,i_xvz) &
931  * mom(k+1,i,j) &
932  + f2h(k,2,i_xvz) &
933  * mom(k,i,j) ) &
934  / ( f2h(k,1,i_xvz) &
935  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
936  + f2h(k,2,i_xvz) &
937  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
938  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
939  * ( f81 * ( val(k,i,j+1)+val(k,i,j) ) &
940  + f82 * ( val(k,i,j+2)+val(k,i,j-1) ) &
941  + f83 * ( val(k,i,j+3)+val(k,i,j-2) ) &
942  + f84 * ( val(k,i,j+4)+val(k,i,j-3) ) ) &
943  + gsqrt(k,i,j) * num_diff(k,i,j)
944  enddo
945  enddo
946  enddo
947  !$omp end do nowait
948 #ifdef DEBUG
949  k = iundef; i = iundef; j = iundef
950 #endif
951 
952  !$omp do OMP_SCHEDULE_ collapse(2)
953  do j = jjs-1, jje
954  do i = iis, iie
955  flux(ke,i,j) = 0.0_rp
956  enddo
957  enddo
958  !$omp end do nowait
959 
960  !$omp end parallel
961 #ifdef DEBUG
962  k = iundef; i = iundef; j = iundef
963 #endif
964 
965  return
966  end subroutine atmos_dyn_fvm_fluxy_xyw_cd8
967 
968 
969  !-----------------------------------------------------------------------------
971  subroutine atmos_dyn_fvm_fluxz_uyz_cd8( &
972  flux, &
973  mom, val, DENS, &
974  GSQRT, J33G, &
975  num_diff, &
976  CDZ, TwoD, &
977  IIS, IIE, JJS, JJE )
978  implicit none
979 
980  real(rp), intent(inout) :: flux (ka,ia,ja)
981  real(rp), intent(in) :: mom (ka,ia,ja)
982  real(rp), intent(in) :: val (ka,ia,ja)
983  real(rp), intent(in) :: dens (ka,ia,ja)
984  real(rp), intent(in) :: gsqrt (ka,ia,ja)
985  real(rp), intent(in) :: j33g
986  real(rp), intent(in) :: num_diff(ka,ia,ja)
987  real(rp), intent(in) :: cdz (ka)
988  logical, intent(in) :: twod
989  integer, intent(in) :: iis, iie, jjs, jje
990 
991  real(rp) :: vel
992  integer :: k, i, j
993  !---------------------------------------------------------------------------
994 
995  !$omp parallel default(none) private(i,j,k,vel) &
996  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J33G,GSQRT,num_diff) &
997  !$omp shared(CDZ,TwoD)
998 
999 
1000  if ( twod ) then
1001 
1002  !$omp do OMP_SCHEDULE_
1003  do j = jjs, jje
1004  do k = ks+3, ke-4
1005 #ifdef DEBUG
1006  call check( __line__, mom(k,i,j) )
1007 
1008  call check( __line__, val(k,i,j) )
1009  call check( __line__, val(k+1,i,j) )
1010 
1011  call check( __line__, val(k-1,i,j) )
1012  call check( __line__, val(k+2,i,j) )
1013 
1014  call check( __line__, val(k-2,i,j) )
1015  call check( __line__, val(k+3,i,j) )
1016 
1017  call check( __line__, val(k-3,i,j) )
1018  call check( __line__, val(k+4,i,j) )
1019 
1020 #endif
1021  i = iis
1022  vel = ( mom(k,i,j) ) &
1023  / ( f2h(k,1,i_xyz) &
1024  * dens(k+1,i,j) &
1025  + f2h(k,2,i_xyz) &
1026  * dens(k,i,j) )
1027  flux(k,i,j) = j33g * vel &
1028  * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
1029  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1030  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1031  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) ) &
1032  + gsqrt(k,i,j) * num_diff(k,i,j)
1033  enddo
1034  enddo
1035  !$omp end do nowait
1036 #ifdef DEBUG
1037  k = iundef; i = iundef; j = iundef
1038 #endif
1039 
1040  !$omp do OMP_SCHEDULE_
1041  do j = jjs, jje
1042 #ifdef DEBUG
1043 
1044  call check( __line__, mom(ks,i ,j) )
1045  call check( __line__, val(ks+1,i,j) )
1046  call check( __line__, val(ks,i,j) )
1047 
1048  call check( __line__, mom(ks+1,i ,j) )
1049  call check( __line__, val(ks+3,i,j) )
1050  call check( __line__, val(ks+2,i,j) )
1051 
1052  call check( __line__, mom(ks+2,i ,j) )
1053  call check( __line__, val(ks+5,i,j) )
1054  call check( __line__, val(ks+4,i,j) )
1055 
1056 #endif
1057  i = iis
1058  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1059  ! The flux at KS-1 can be non-zero.
1060  ! To reduce calculations, all the fluxes are set to zero.
1061  flux(ks-1,i,j) = 0.0_rp
1062 
1063  vel = ( mom(ks,i,j) ) &
1064  / ( f2h(ks,1,i_xyz) &
1065  * dens(ks+1,i,j) &
1066  + f2h(ks,2,i_xyz) &
1067  * dens(ks,i,j) )
1068  flux(ks,i,j) = j33g * vel &
1069  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) &
1070  + gsqrt(ks,i,j) * num_diff(ks,i,j)
1071  vel = ( mom(ke-1,i,j) ) &
1072  / ( f2h(ke-1,1,i_xyz) &
1073  * dens(ke,i,j) &
1074  + f2h(ke-1,2,i_xyz) &
1075  * dens(ke-1,i,j) )
1076  flux(ke-1,i,j) = j33g * vel &
1077  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) ) &
1078  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
1079 
1080  vel = ( mom(ks+1,i,j) ) &
1081  / ( f2h(ks+1,1,i_xyz) &
1082  * dens(ks+2,i,j) &
1083  + f2h(ks+1,2,i_xyz) &
1084  * dens(ks+1,i,j) )
1085  flux(ks+1,i,j) = j33g * vel &
1086  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
1087  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) ) &
1088  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j)
1089  vel = ( mom(ke-2,i,j) ) &
1090  / ( f2h(ke-2,1,i_xyz) &
1091  * dens(ke-1,i,j) &
1092  + f2h(ke-2,2,i_xyz) &
1093  * dens(ke-2,i,j) )
1094  flux(ke-2,i,j) = j33g * vel &
1095  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
1096  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) ) &
1097  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j)
1098 
1099  vel = ( mom(ks+2,i,j) ) &
1100  / ( f2h(ks+2,1,i_xyz) &
1101  * dens(ks+3,i,j) &
1102  + f2h(ks+2,2,i_xyz) &
1103  * dens(ks+2,i,j) )
1104  flux(ks+2,i,j) = j33g * vel &
1105  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
1106  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
1107  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) ) &
1108  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j)
1109  vel = ( mom(ke-3,i,j) ) &
1110  / ( f2h(ke-3,1,i_xyz) &
1111  * dens(ke-2,i,j) &
1112  + f2h(ke-3,2,i_xyz) &
1113  * dens(ke-3,i,j) )
1114  flux(ke-3,i,j) = j33g * vel &
1115  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
1116  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
1117  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) ) &
1118  + gsqrt(ke-3,i,j) * num_diff(ke-3,i,j)
1119 
1120  flux(ke,i,j) = 0.0_rp
1121  enddo
1122  !$omp end do nowait
1123 
1124  else
1125 
1126 
1127  !$omp do OMP_SCHEDULE_ collapse(2)
1128  do j = jjs, jje
1129  do i = iis, iie
1130  do k = ks+3, ke-4
1131 #ifdef DEBUG
1132  call check( __line__, mom(k,i,j) )
1133  call check( __line__, mom(k,i+1,j) )
1134 
1135  call check( __line__, val(k,i,j) )
1136  call check( __line__, val(k+1,i,j) )
1137 
1138  call check( __line__, val(k-1,i,j) )
1139  call check( __line__, val(k+2,i,j) )
1140 
1141  call check( __line__, val(k-2,i,j) )
1142  call check( __line__, val(k+3,i,j) )
1143 
1144  call check( __line__, val(k-3,i,j) )
1145  call check( __line__, val(k+4,i,j) )
1146 
1147 #endif
1148  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i+1,j) ) ) &
1149  / ( f2h(k,1,i_uyz) &
1150  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
1151  + f2h(k,2,i_uyz) &
1152  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
1153  flux(k,i,j) = j33g * vel &
1154  * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
1155  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1156  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1157  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) ) &
1158  + gsqrt(k,i,j) * num_diff(k,i,j)
1159  enddo
1160  enddo
1161  enddo
1162  !$omp end do nowait
1163 #ifdef DEBUG
1164  k = iundef; i = iundef; j = iundef
1165 #endif
1166 
1167  !$omp do OMP_SCHEDULE_ collapse(2)
1168  do j = jjs, jje
1169  do i = iis, iie
1170 #ifdef DEBUG
1171 
1172  call check( __line__, mom(ks,i ,j) )
1173  call check( __line__, mom(ks,i+1,j) )
1174  call check( __line__, val(ks+1,i,j) )
1175  call check( __line__, val(ks,i,j) )
1176 
1177  call check( __line__, mom(ks+1,i ,j) )
1178  call check( __line__, mom(ks+1,i+1,j) )
1179  call check( __line__, val(ks+3,i,j) )
1180  call check( __line__, val(ks+2,i,j) )
1181 
1182  call check( __line__, mom(ks+2,i ,j) )
1183  call check( __line__, mom(ks+2,i+1,j) )
1184  call check( __line__, val(ks+5,i,j) )
1185  call check( __line__, val(ks+4,i,j) )
1186 
1187 #endif
1188  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1189  ! The flux at KS-1 can be non-zero.
1190  ! To reduce calculations, all the fluxes are set to zero.
1191  flux(ks-1,i,j) = 0.0_rp
1192 
1193  vel = ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i+1,j) ) ) &
1194  / ( f2h(ks,1,i_uyz) &
1195  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
1196  + f2h(ks,2,i_uyz) &
1197  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
1198  flux(ks,i,j) = j33g * vel &
1199  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) &
1200  + gsqrt(ks,i,j) * num_diff(ks,i,j)
1201  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i+1,j) ) ) &
1202  / ( f2h(ke-1,1,i_uyz) &
1203  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
1204  + f2h(ke-1,2,i_uyz) &
1205  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
1206  flux(ke-1,i,j) = j33g * vel &
1207  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) ) &
1208  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
1209 
1210  vel = ( 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i+1,j) ) ) &
1211  / ( f2h(ks+1,1,i_uyz) &
1212  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) &
1213  + f2h(ks+1,2,i_uyz) &
1214  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) )
1215  flux(ks+1,i,j) = j33g * vel &
1216  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
1217  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) ) &
1218  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j)
1219  vel = ( 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i+1,j) ) ) &
1220  / ( f2h(ke-2,1,i_uyz) &
1221  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) &
1222  + f2h(ke-2,2,i_uyz) &
1223  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) )
1224  flux(ke-2,i,j) = j33g * vel &
1225  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
1226  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) ) &
1227  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j)
1228 
1229  vel = ( 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i+1,j) ) ) &
1230  / ( f2h(ks+2,1,i_uyz) &
1231  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i+1,j) ) &
1232  + f2h(ks+2,2,i_uyz) &
1233  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) )
1234  flux(ks+2,i,j) = j33g * vel &
1235  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
1236  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
1237  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) ) &
1238  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j)
1239  vel = ( 0.5_rp * ( mom(ke-3,i,j)+mom(ke-3,i+1,j) ) ) &
1240  / ( f2h(ke-3,1,i_uyz) &
1241  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) &
1242  + f2h(ke-3,2,i_uyz) &
1243  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i+1,j) ) )
1244  flux(ke-3,i,j) = j33g * vel &
1245  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
1246  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
1247  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) ) &
1248  + gsqrt(ke-3,i,j) * num_diff(ke-3,i,j)
1249 
1250  flux(ke,i,j) = 0.0_rp
1251  enddo
1252  enddo
1253  !$omp end do nowait
1254 
1255  end if
1256 
1257 
1258  !$omp end parallel
1259 #ifdef DEBUG
1260  k = iundef; i = iundef; j = iundef
1261 #endif
1262 
1263  return
1264  end subroutine atmos_dyn_fvm_fluxz_uyz_cd8
1265 
1266  !-----------------------------------------------------------------------------
1268  subroutine atmos_dyn_fvm_fluxj13_uyz_cd8( &
1269  flux, &
1270  mom, val, DENS, &
1271  GSQRT, J13G, MAPF, &
1272  CDZ, TwoD, &
1273  IIS, IIE, JJS, JJE )
1274  implicit none
1275 
1276  real(rp), intent(inout) :: flux (ka,ia,ja)
1277  real(rp), intent(in) :: mom (ka,ia,ja)
1278  real(rp), intent(in) :: val (ka,ia,ja)
1279  real(rp), intent(in) :: dens (ka,ia,ja)
1280  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1281  real(rp), intent(in) :: j13g (ka,ia,ja)
1282  real(rp), intent(in) :: mapf ( ia,ja,2)
1283  real(rp), intent(in) :: cdz (ka)
1284  logical, intent(in) :: twod
1285  integer, intent(in) :: iis, iie, jjs, jje
1286 
1287  real(rp) :: vel
1288  integer :: k, i, j
1289  !---------------------------------------------------------------------------
1290 
1291  !$omp parallel default(none) private(i,j,k,vel) &
1292  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J13G,MAPF) &
1293  !$omp shared(GSQRT,CDZ,TwoD)
1294 
1295 
1296 
1297  !$omp do OMP_SCHEDULE_ collapse(2)
1298  do j = jjs, jje
1299  do i = iis, iie
1300  do k = ks+3, ke-4
1301  vel = ( f2h(k,1,i_uyz) &
1302  * mom(k+1,i,j) &
1303  + f2h(k,2,i_uyz) &
1304  * mom(k,i,j) ) &
1305  / ( f2h(k,1,i_uyz) &
1306  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
1307  + f2h(k,2,i_uyz) &
1308  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
1309  vel = vel * j13g(k,i,j)
1310  flux(k,i,j) = vel / mapf(i,j,+2) &
1311  * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
1312  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1313  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1314  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) )
1315  enddo
1316  enddo
1317  enddo
1318  !$omp end do nowait
1319 
1320  !$omp do OMP_SCHEDULE_ collapse(2)
1321  do j = jjs, jje
1322  do i = iis, iie
1323  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1324  ! The flux at KS-1 can be non-zero.
1325  ! To reduce calculations, all the fluxes are set to zero.
1326  flux(ks-1,i,j) = 0.0_rp
1327 
1328  vel = ( f2h(ks,1,i_uyz) &
1329  * mom(ks+1,i,j) &
1330  + f2h(ks,2,i_uyz) &
1331  * mom(ks,i,j) ) &
1332  / ( f2h(ks,1,i_uyz) &
1333  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
1334  + f2h(ks,2,i_uyz) &
1335  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
1336  vel = vel * j13g(ks,i,j)
1337  flux(ks,i,j) = vel / mapf(i,j,+2) &
1338  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
1339 
1340  vel = ( f2h(ke-1,1,i_uyz) &
1341  * mom(ke,i,j) &
1342  + f2h(ke-1,2,i_uyz) &
1343  * mom(ke-1,i,j) ) &
1344  / ( f2h(ke-1,1,i_uyz) &
1345  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
1346  + f2h(ke-1,2,i_uyz) &
1347  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
1348  vel = vel * j13g(ke-1,i,j)
1349  flux(ke-1,i,j) = vel / mapf(i,j,+2) &
1350  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
1351 
1352  vel = ( f2h(ks+1,1,i_uyz) &
1353  * mom(ks+2,i,j) &
1354  + f2h(ks+1,2,i_uyz) &
1355  * mom(ks+1,i,j) ) &
1356  / ( f2h(ks+1,1,i_uyz) &
1357  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) &
1358  + f2h(ks+1,2,i_uyz) &
1359  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) )
1360  vel = vel * j13g(ks+1,i,j)
1361  flux(ks+1,i,j) = vel / mapf(i,j,+2) &
1362  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
1363  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) )
1364 
1365  vel = ( f2h(ke-2,1,i_uyz) &
1366  * mom(ke-1,i,j) &
1367  + f2h(ke-2,2,i_uyz) &
1368  * mom(ke-2,i,j) ) &
1369  / ( f2h(ke-2,1,i_uyz) &
1370  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) &
1371  + f2h(ke-2,2,i_uyz) &
1372  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) )
1373  vel = vel * j13g(ke-2,i,j)
1374  flux(ke-2,i,j) = vel / mapf(i,j,+2) &
1375  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
1376  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) )
1377 
1378  vel = ( f2h(ks+2,1,i_uyz) &
1379  * mom(ks+3,i,j) &
1380  + f2h(ks+2,2,i_uyz) &
1381  * mom(ks+2,i,j) ) &
1382  / ( f2h(ks+2,1,i_uyz) &
1383  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i+1,j) ) &
1384  + f2h(ks+2,2,i_uyz) &
1385  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) )
1386  vel = vel * j13g(ks+2,i,j)
1387  flux(ks+2,i,j) = vel / mapf(i,j,+2) &
1388  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
1389  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
1390  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) )
1391 
1392  vel = ( f2h(ke-3,1,i_uyz) &
1393  * mom(ke-2,i,j) &
1394  + f2h(ke-3,2,i_uyz) &
1395  * mom(ke-3,i,j) ) &
1396  / ( f2h(ke-3,1,i_uyz) &
1397  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) &
1398  + f2h(ke-3,2,i_uyz) &
1399  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i+1,j) ) )
1400  vel = vel * j13g(ke-3,i,j)
1401  flux(ke-3,i,j) = vel / mapf(i,j,+2) &
1402  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
1403  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
1404  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) )
1405 
1406  flux(ke ,i,j) = 0.0_rp
1407  enddo
1408  enddo
1409  !$omp end do nowait
1410 
1411 
1412 
1413  !$omp end parallel
1414  return
1415  end subroutine atmos_dyn_fvm_fluxj13_uyz_cd8
1416 
1417  !-----------------------------------------------------------------------------
1419  subroutine atmos_dyn_fvm_fluxj23_uyz_cd8( &
1420  flux, &
1421  mom, val, DENS, &
1422  GSQRT, J23G, MAPF, &
1423  CDZ, TwoD, &
1424  IIS, IIE, JJS, JJE )
1425  implicit none
1426 
1427  real(rp), intent(inout) :: flux (ka,ia,ja)
1428  real(rp), intent(in) :: mom (ka,ia,ja)
1429  real(rp), intent(in) :: val (ka,ia,ja)
1430  real(rp), intent(in) :: dens (ka,ia,ja)
1431  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1432  real(rp), intent(in) :: j23g (ka,ia,ja)
1433  real(rp), intent(in) :: mapf ( ia,ja,2)
1434  real(rp), intent(in) :: cdz (ka)
1435  logical, intent(in) :: twod
1436  integer, intent(in) :: iis, iie, jjs, jje
1437 
1438  real(rp) :: vel
1439  integer :: k, i, j
1440  !---------------------------------------------------------------------------
1441 
1442  !$omp parallel default(none) private(i,j,k,vel) &
1443  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J23G,MAPF) &
1444  !$omp shared(GSQRT,CDZ,TwoD)
1445 
1446 
1447  if ( twod ) then
1448 
1449  !$omp do OMP_SCHEDULE_
1450  do j = jjs, jje
1451  do k = ks+3, ke-4
1452  i = iis
1453  vel = ( f2h(k,1,i_xyz) &
1454  * 0.5_rp * ( mom(k+1,i,j)+mom(k+1,i,j-1) ) &
1455  + f2h(k,2,i_xyz) &
1456  * 0.5_rp * ( mom(k,i,j)+mom(k,i,j-1) ) ) &
1457  / ( f2h(k,1,i_xyz) &
1458  * dens(k+1,i,j) &
1459  + f2h(k,2,i_xyz) &
1460  * dens(k,i,j) )
1461  vel = vel * j23g(k,i,j)
1462  flux(k,i,j) = vel * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
1463  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1464  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1465  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) )
1466  enddo
1467  enddo
1468  !$omp end do nowait
1469 
1470  !$omp do OMP_SCHEDULE_
1471  do j = jjs, jje
1472  i = iis
1473  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1474  ! The flux at KS-1 can be non-zero.
1475  ! To reduce calculations, all the fluxes are set to zero.
1476  flux(ks-1,i,j) = 0.0_rp
1477 
1478  vel = ( f2h(ks,1,i_xyz) &
1479  * 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i,j-1) ) &
1480  + f2h(ks,2,i_xyz) &
1481  * 0.5_rp * ( mom(ks,i,j)+mom(ks,i,j-1) ) ) &
1482  / ( f2h(ks,1,i_xyz) &
1483  * dens(ks+1,i,j) &
1484  + f2h(ks,2,i_xyz) &
1485  * dens(ks,i,j) )
1486  vel = vel * j23g(ks,i,j)
1487  flux(ks,i,j) = vel / mapf(i,j,+1) &
1488  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
1489 
1490  vel = ( f2h(ke-1,1,i_xyz) &
1491  * 0.5_rp * ( mom(ke,i,j)+mom(ke,i,j-1) ) &
1492  + f2h(ke-1,2,i_xyz) &
1493  * 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j-1) ) ) &
1494  / ( f2h(ke-1,1,i_xyz) &
1495  * dens(ke,i,j) &
1496  + f2h(ke-1,2,i_xyz) &
1497  * dens(ke-1,i,j) )
1498  vel = vel * j23g(ke-1,i,j)
1499  flux(ke-1,i,j) = vel / mapf(i,j,+1) &
1500  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
1501 
1502  vel = ( f2h(ks+1,1,i_xyz) &
1503  * 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i,j-1) ) &
1504  + f2h(ks+1,2,i_xyz) &
1505  * 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i,j-1) ) ) &
1506  / ( f2h(ks+1,1,i_xyz) &
1507  * dens(ks+2,i,j) &
1508  + f2h(ks+1,2,i_xyz) &
1509  * dens(ks+1,i,j) )
1510  vel = vel * j23g(ks+1,i,j)
1511  flux(ks+1,i,j) = vel / mapf(i,j,+1) &
1512  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
1513  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) )
1514 
1515  vel = ( f2h(ke-2,1,i_xyz) &
1516  * 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j-1) ) &
1517  + f2h(ke-2,2,i_xyz) &
1518  * 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i,j-1) ) ) &
1519  / ( f2h(ke-2,1,i_xyz) &
1520  * dens(ke-1,i,j) &
1521  + f2h(ke-2,2,i_xyz) &
1522  * dens(ke-2,i,j) )
1523  vel = vel * j23g(ke-2,i,j)
1524  flux(ke-2,i,j) = vel / mapf(i,j,+1) &
1525  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
1526  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) )
1527 
1528  vel = ( f2h(ks+2,1,i_xyz) &
1529  * 0.5_rp * ( mom(ks+3,i,j)+mom(ks+3,i,j-1) ) &
1530  + f2h(ks+2,2,i_xyz) &
1531  * 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i,j-1) ) ) &
1532  / ( f2h(ks+2,1,i_xyz) &
1533  * dens(ks+3,i,j) &
1534  + f2h(ks+2,2,i_xyz) &
1535  * dens(ks+2,i,j) )
1536  vel = vel * j23g(ks+2,i,j)
1537  flux(ks+2,i,j) = vel / mapf(i,j,+1) &
1538  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
1539  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
1540  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) )
1541 
1542  vel = ( f2h(ke-3,1,i_xyz) &
1543  * 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i,j-1) ) &
1544  + f2h(ke-3,2,i_xyz) &
1545  * 0.5_rp * ( mom(ke-3,i,j)+mom(ke-3,i,j-1) ) ) &
1546  / ( f2h(ke-3,1,i_xyz) &
1547  * dens(ke-2,i,j) &
1548  + f2h(ke-3,2,i_xyz) &
1549  * dens(ke-3,i,j) )
1550  vel = vel * j23g(ke-3,i,j)
1551  flux(ke-3,i,j) = vel / mapf(i,j,+1) &
1552  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
1553  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
1554  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) )
1555 
1556  flux(ke ,i,j) = 0.0_rp
1557  enddo
1558  !$omp end do nowait
1559 
1560  else
1561 
1562 
1563  !$omp do OMP_SCHEDULE_ collapse(2)
1564  do j = jjs, jje
1565  do i = iis, iie
1566  do k = ks+3, ke-4
1567  vel = ( f2h(k,1,i_uyz) &
1568  * 0.25_rp * ( mom(k+1,i,j)+mom(k+1,i+1,j)+mom(k+1,i,j-1)+mom(k+1,i+1,j-1) ) &
1569  + f2h(k,2,i_uyz) &
1570  * 0.25_rp * ( mom(k,i,j)+mom(k,i+1,j)+mom(k,i,j-1)+mom(k,i+1,j-1) ) ) &
1571  / ( f2h(k,1,i_uyz) &
1572  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
1573  + f2h(k,2,i_uyz) &
1574  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
1575  vel = vel * j23g(k,i,j)
1576  flux(k,i,j) = vel / mapf(i,j,+1) &
1577  * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
1578  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1579  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1580  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) )
1581  enddo
1582  enddo
1583  enddo
1584  !$omp end do nowait
1585 
1586  !$omp do OMP_SCHEDULE_ collapse(2)
1587  do j = jjs, jje
1588  do i = iis, iie
1589  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1590  ! The flux at KS-1 can be non-zero.
1591  ! To reduce calculations, all the fluxes are set to zero.
1592  flux(ks-1,i,j) = 0.0_rp
1593 
1594  vel = ( f2h(ks,1,i_uyz) &
1595  * 0.25_rp * ( mom(ks+1,i,j)+mom(ks+1,i+1,j)+mom(ks+1,i,j-1)+mom(ks+1,i+1,j-1) ) &
1596  + f2h(ks,2,i_uyz) &
1597  * 0.25_rp * ( mom(ks,i,j)+mom(ks,i+1,j)+mom(ks,i,j-1)+mom(ks,i+1,j-1) ) ) &
1598  / ( f2h(ks,1,i_uyz) &
1599  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
1600  + f2h(ks,2,i_uyz) &
1601  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
1602  vel = vel * j23g(ks,i,j)
1603  flux(ks,i,j) = vel / mapf(i,j,+1) &
1604  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
1605 
1606  vel = ( f2h(ke-1,1,i_uyz) &
1607  * 0.25_rp * ( mom(ke,i,j)+mom(ke,i+1,j)+mom(ke,i,j-1)+mom(ke,i+1,j-1) ) &
1608  + f2h(ke-1,2,i_uyz) &
1609  * 0.25_rp * ( mom(ke-1,i,j)+mom(ke-1,i+1,j)+mom(ke-1,i,j-1)+mom(ke-1,i+1,j-1) ) ) &
1610  / ( f2h(ke-1,1,i_uyz) &
1611  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
1612  + f2h(ke-1,2,i_uyz) &
1613  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
1614  vel = vel * j23g(ke-1,i,j)
1615  flux(ke-1,i,j) = vel / mapf(i,j,+1) &
1616  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
1617 
1618  vel = ( f2h(ks+1,1,i_uyz) &
1619  * 0.25_rp * ( mom(ks+2,i,j)+mom(ks+2,i+1,j)+mom(ks+2,i,j-1)+mom(ks+2,i+1,j-1) ) &
1620  + f2h(ks+1,2,i_uyz) &
1621  * 0.25_rp * ( mom(ks+1,i,j)+mom(ks+1,i+1,j)+mom(ks+1,i,j-1)+mom(ks+1,i+1,j-1) ) ) &
1622  / ( f2h(ks+1,1,i_uyz) &
1623  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) &
1624  + f2h(ks+1,2,i_uyz) &
1625  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) )
1626  vel = vel * j23g(ks+1,i,j)
1627  flux(ks+1,i,j) = vel / mapf(i,j,+1) &
1628  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
1629  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) )
1630 
1631  vel = ( f2h(ke-2,1,i_uyz) &
1632  * 0.25_rp * ( mom(ke-1,i,j)+mom(ke-1,i+1,j)+mom(ke-1,i,j-1)+mom(ke-1,i+1,j-1) ) &
1633  + f2h(ke-2,2,i_uyz) &
1634  * 0.25_rp * ( mom(ke-2,i,j)+mom(ke-2,i+1,j)+mom(ke-2,i,j-1)+mom(ke-2,i+1,j-1) ) ) &
1635  / ( f2h(ke-2,1,i_uyz) &
1636  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) &
1637  + f2h(ke-2,2,i_uyz) &
1638  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) )
1639  vel = vel * j23g(ke-2,i,j)
1640  flux(ke-2,i,j) = vel / mapf(i,j,+1) &
1641  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
1642  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) )
1643 
1644  vel = ( f2h(ks+2,1,i_uyz) &
1645  * 0.25_rp * ( mom(ks+3,i,j)+mom(ks+3,i+1,j)+mom(ks+3,i,j-1)+mom(ks+3,i+1,j-1) ) &
1646  + f2h(ks+2,2,i_uyz) &
1647  * 0.25_rp * ( mom(ks+2,i,j)+mom(ks+2,i+1,j)+mom(ks+2,i,j-1)+mom(ks+2,i+1,j-1) ) ) &
1648  / ( f2h(ks+2,1,i_uyz) &
1649  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i+1,j) ) &
1650  + f2h(ks+2,2,i_uyz) &
1651  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) )
1652  vel = vel * j23g(ks+2,i,j)
1653  flux(ks+2,i,j) = vel / mapf(i,j,+1) &
1654  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
1655  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
1656  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) )
1657 
1658  vel = ( f2h(ke-3,1,i_uyz) &
1659  * 0.25_rp * ( mom(ke-2,i,j)+mom(ke-2,i+1,j)+mom(ke-2,i,j-1)+mom(ke-2,i+1,j-1) ) &
1660  + f2h(ke-3,2,i_uyz) &
1661  * 0.25_rp * ( mom(ke-3,i,j)+mom(ke-3,i+1,j)+mom(ke-3,i,j-1)+mom(ke-3,i+1,j-1) ) ) &
1662  / ( f2h(ke-3,1,i_uyz) &
1663  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) &
1664  + f2h(ke-3,2,i_uyz) &
1665  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i+1,j) ) )
1666  vel = vel * j23g(ke-3,i,j)
1667  flux(ke-3,i,j) = vel / mapf(i,j,+1) &
1668  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
1669  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
1670  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) )
1671 
1672  flux(ke ,i,j) = 0.0_rp
1673  enddo
1674  enddo
1675  !$omp end do nowait
1676 
1677 
1678  end if
1679 
1680 
1681  !$omp end parallel
1682  return
1683  end subroutine atmos_dyn_fvm_fluxj23_uyz_cd8
1684 
1685  !-----------------------------------------------------------------------------
1687  subroutine atmos_dyn_fvm_fluxx_uyz_cd8( &
1688  flux, &
1689  mom, val, DENS, &
1690  GSQRT, MAPF, &
1691  num_diff, &
1692  CDZ, TwoD, &
1693  IIS, IIE, JJS, JJE )
1694  implicit none
1695 
1696  real(rp), intent(inout) :: flux (ka,ia,ja)
1697  real(rp), intent(in) :: mom (ka,ia,ja)
1698  real(rp), intent(in) :: val (ka,ia,ja)
1699  real(rp), intent(in) :: dens (ka,ia,ja)
1700  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1701  real(rp), intent(in) :: mapf ( ia,ja,2)
1702  real(rp), intent(in) :: num_diff(ka,ia,ja)
1703  real(rp), intent(in) :: cdz (ka)
1704  logical, intent(in) :: twod
1705  integer, intent(in) :: iis, iie, jjs, jje
1706 
1707  real(rp) :: vel
1708  integer :: k, i, j
1709  !---------------------------------------------------------------------------
1710 
1711  ! note that x-index is added by -1
1712 
1713 
1714 
1715  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
1716  !$omp private(vel) &
1717  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff)
1718  do j = jjs, jje
1719  do i = iis, iie+1
1720  do k = ks, ke
1721 #ifdef DEBUG
1722  call check( __line__, mom(k,i ,j) )
1723  call check( __line__, mom(k,i-1,j) )
1724 
1725  call check( __line__, val(k,i-1,j) )
1726  call check( __line__, val(k,i,j) )
1727 
1728  call check( __line__, val(k,i-2,j) )
1729  call check( __line__, val(k,i+1,j) )
1730 
1731  call check( __line__, val(k,i-3,j) )
1732  call check( __line__, val(k,i+2,j) )
1733 
1734  call check( __line__, val(k,i-4,j) )
1735  call check( __line__, val(k,i+3,j) )
1736 
1737 #endif
1738  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i-1,j) ) ) &
1739  / ( dens(k,i,j) )
1740  flux(k,i-1,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
1741  * ( f81 * ( val(k,i,j)+val(k,i-1,j) ) &
1742  + f82 * ( val(k,i+1,j)+val(k,i-2,j) ) &
1743  + f83 * ( val(k,i+2,j)+val(k,i-3,j) ) &
1744  + f84 * ( val(k,i+3,j)+val(k,i-4,j) ) ) &
1745  + gsqrt(k,i,j) * num_diff(k,i,j)
1746  enddo
1747  enddo
1748  enddo
1749 #ifdef DEBUG
1750  k = iundef; i = iundef; j = iundef
1751 #endif
1752 
1753 
1754 
1755  return
1756  end subroutine atmos_dyn_fvm_fluxx_uyz_cd8
1757 
1758  !-----------------------------------------------------------------------------
1760  subroutine atmos_dyn_fvm_fluxy_uyz_cd8( &
1761  flux, &
1762  mom, val, DENS, &
1763  GSQRT, MAPF, &
1764  num_diff, &
1765  CDZ, TwoD, &
1766  IIS, IIE, JJS, JJE )
1767  implicit none
1768 
1769  real(rp), intent(inout) :: flux (ka,ia,ja)
1770  real(rp), intent(in) :: mom (ka,ia,ja)
1771  real(rp), intent(in) :: val (ka,ia,ja)
1772  real(rp), intent(in) :: dens (ka,ia,ja)
1773  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1774  real(rp), intent(in) :: mapf ( ia,ja,2)
1775  real(rp), intent(in) :: num_diff(ka,ia,ja)
1776  real(rp), intent(in) :: cdz (ka)
1777  logical, intent(in) :: twod
1778  integer, intent(in) :: iis, iie, jjs, jje
1779 
1780  real(rp) :: vel
1781  integer :: k, i, j
1782  !---------------------------------------------------------------------------
1783 
1784 
1785 
1786  if ( twod ) then
1787 
1788  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ &
1789  !$omp private(vel) &
1790  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff,TwoD)
1791  do j = jjs-1, jje
1792  do k = ks, ke
1793  i = iis
1794 #ifdef DEBUG
1795  call check( __line__, mom(k,i ,j) )
1796 
1797  call check( __line__, val(k,i,j) )
1798  call check( __line__, val(k,i,j+1) )
1799 
1800  call check( __line__, val(k,i,j-1) )
1801  call check( __line__, val(k,i,j+2) )
1802 
1803  call check( __line__, val(k,i,j-2) )
1804  call check( __line__, val(k,i,j+3) )
1805 
1806  call check( __line__, val(k,i,j-3) )
1807  call check( __line__, val(k,i,j+4) )
1808 
1809 #endif
1810  vel = ( mom(k,i,j) ) &
1811  / ( 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
1812  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
1813  * ( f81 * ( val(k,i,j+1)+val(k,i,j) ) &
1814  + f82 * ( val(k,i,j+2)+val(k,i,j-1) ) &
1815  + f83 * ( val(k,i,j+3)+val(k,i,j-2) ) &
1816  + f84 * ( val(k,i,j+4)+val(k,i,j-3) ) ) &
1817  + gsqrt(k,i,j) * num_diff(k,i,j)
1818  enddo
1819  enddo
1820 #ifdef DEBUG
1821  k = iundef; i = iundef; j = iundef
1822 #endif
1823 
1824  else
1825 
1826 
1827  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
1828  !$omp private(vel) &
1829  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff)
1830  do j = jjs-1, jje
1831  do i = iis, iie
1832  do k = ks, ke
1833 #ifdef DEBUG
1834  call check( __line__, mom(k,i ,j) )
1835  call check( __line__, mom(k,i-1,j) )
1836 
1837  call check( __line__, val(k,i,j) )
1838  call check( __line__, val(k,i,j+1) )
1839 
1840  call check( __line__, val(k,i,j-1) )
1841  call check( __line__, val(k,i,j+2) )
1842 
1843  call check( __line__, val(k,i,j-2) )
1844  call check( __line__, val(k,i,j+3) )
1845 
1846  call check( __line__, val(k,i,j-3) )
1847  call check( __line__, val(k,i,j+4) )
1848 
1849 #endif
1850  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i+1,j) ) ) &
1851  / ( 0.25_rp * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) )
1852  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
1853  * ( f81 * ( val(k,i,j+1)+val(k,i,j) ) &
1854  + f82 * ( val(k,i,j+2)+val(k,i,j-1) ) &
1855  + f83 * ( val(k,i,j+3)+val(k,i,j-2) ) &
1856  + f84 * ( val(k,i,j+4)+val(k,i,j-3) ) ) &
1857  + gsqrt(k,i,j) * num_diff(k,i,j)
1858  enddo
1859  enddo
1860  enddo
1861 #ifdef DEBUG
1862  k = iundef; i = iundef; j = iundef
1863 #endif
1864 
1865 
1866  end if
1867 
1868 
1869  return
1870  end subroutine atmos_dyn_fvm_fluxy_uyz_cd8
1871 
1872 
1873 
1874  !-----------------------------------------------------------------------------
1876  subroutine atmos_dyn_fvm_fluxz_xvz_cd8( &
1877  flux, &
1878  mom, val, DENS, &
1879  GSQRT, J33G, &
1880  num_diff, &
1881  CDZ, TwoD, &
1882  IIS, IIE, JJS, JJE )
1883  implicit none
1884 
1885  real(rp), intent(inout) :: flux (ka,ia,ja)
1886  real(rp), intent(in) :: mom (ka,ia,ja)
1887  real(rp), intent(in) :: val (ka,ia,ja)
1888  real(rp), intent(in) :: dens (ka,ia,ja)
1889  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1890  real(rp), intent(in) :: j33g
1891  real(rp), intent(in) :: num_diff(ka,ia,ja)
1892  real(rp), intent(in) :: cdz (ka)
1893  logical, intent(in) :: twod
1894  integer, intent(in) :: iis, iie, jjs, jje
1895 
1896  real(rp) :: vel
1897  integer :: k, i, j
1898  !---------------------------------------------------------------------------
1899 
1900  !$omp parallel default(none) private(i,j,k,vel) &
1901  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J33G,GSQRT,num_diff) &
1902  !$omp shared(CDZ,TwoD)
1903 
1904 
1905  !$omp do OMP_SCHEDULE_ collapse(2)
1906  do j = jjs, jje
1907  do i = iis, iie
1908  do k = ks+3, ke-4
1909 #ifdef DEBUG
1910  call check( __line__, mom(k,i,j) )
1911  call check( __line__, mom(k,i,j+1) )
1912 
1913  call check( __line__, val(k,i,j) )
1914  call check( __line__, val(k+1,i,j) )
1915 
1916  call check( __line__, val(k-1,i,j) )
1917  call check( __line__, val(k+2,i,j) )
1918 
1919  call check( __line__, val(k-2,i,j) )
1920  call check( __line__, val(k+3,i,j) )
1921 
1922  call check( __line__, val(k-3,i,j) )
1923  call check( __line__, val(k+4,i,j) )
1924 
1925 #endif
1926  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j+1) ) ) &
1927  / ( f2h(k,1,i_xvz) &
1928  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
1929  + f2h(k,2,i_xvz) &
1930  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
1931  flux(k,i,j) = j33g * vel &
1932  * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
1933  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1934  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1935  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) ) &
1936  + gsqrt(k,i,j) * num_diff(k,i,j)
1937  enddo
1938  enddo
1939  enddo
1940  !$omp end do nowait
1941 #ifdef DEBUG
1942  k = iundef; i = iundef; j = iundef
1943 #endif
1944 
1945  !$omp do OMP_SCHEDULE_ collapse(2)
1946  do j = jjs, jje
1947  do i = iis, iie
1948 #ifdef DEBUG
1949 
1950  call check( __line__, mom(ks,i ,j) )
1951  call check( __line__, mom(ks,i,j+1) )
1952  call check( __line__, val(ks+1,i,j) )
1953  call check( __line__, val(ks,i,j) )
1954 
1955  call check( __line__, mom(ks+1,i ,j) )
1956  call check( __line__, mom(ks+1,i,j+1) )
1957  call check( __line__, val(ks+3,i,j) )
1958  call check( __line__, val(ks+2,i,j) )
1959 
1960  call check( __line__, mom(ks+2,i ,j) )
1961  call check( __line__, mom(ks+2,i,j+1) )
1962  call check( __line__, val(ks+5,i,j) )
1963  call check( __line__, val(ks+4,i,j) )
1964 
1965 #endif
1966  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1967  ! The flux at KS-1 can be non-zero.
1968  ! To reduce calculations, all the fluxes are set to zero.
1969  flux(ks-1,i,j) = 0.0_rp
1970 
1971  vel = ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i,j+1) ) ) &
1972  / ( f2h(ks,1,i_xvz) &
1973  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
1974  + f2h(ks,2,i_xvz) &
1975  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
1976  flux(ks,i,j) = j33g * vel &
1977  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) &
1978  + gsqrt(ks,i,j) * num_diff(ks,i,j)
1979  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j+1) ) ) &
1980  / ( f2h(ke-1,1,i_xvz) &
1981  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
1982  + f2h(ke-1,2,i_xvz) &
1983  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
1984  flux(ke-1,i,j) = j33g * vel &
1985  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) ) &
1986  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
1987 
1988  vel = ( 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i,j+1) ) ) &
1989  / ( f2h(ks+1,1,i_xvz) &
1990  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) &
1991  + f2h(ks+1,2,i_xvz) &
1992  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) )
1993  flux(ks+1,i,j) = j33g * vel &
1994  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
1995  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) ) &
1996  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j)
1997  vel = ( 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i,j+1) ) ) &
1998  / ( f2h(ke-2,1,i_xvz) &
1999  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) &
2000  + f2h(ke-2,2,i_xvz) &
2001  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) )
2002  flux(ke-2,i,j) = j33g * vel &
2003  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
2004  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) ) &
2005  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j)
2006 
2007  vel = ( 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i,j+1) ) ) &
2008  / ( f2h(ks+2,1,i_xvz) &
2009  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i,j+1) ) &
2010  + f2h(ks+2,2,i_xvz) &
2011  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) )
2012  flux(ks+2,i,j) = j33g * vel &
2013  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
2014  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
2015  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) ) &
2016  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j)
2017  vel = ( 0.5_rp * ( mom(ke-3,i,j)+mom(ke-3,i,j+1) ) ) &
2018  / ( f2h(ke-3,1,i_xvz) &
2019  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) &
2020  + f2h(ke-3,2,i_xvz) &
2021  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i,j+1) ) )
2022  flux(ke-3,i,j) = j33g * vel &
2023  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
2024  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
2025  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) ) &
2026  + gsqrt(ke-3,i,j) * num_diff(ke-3,i,j)
2027 
2028  flux(ke,i,j) = 0.0_rp
2029  enddo
2030  enddo
2031  !$omp end do nowait
2032 
2033 
2034  !$omp end parallel
2035 #ifdef DEBUG
2036  k = iundef; i = iundef; j = iundef
2037 #endif
2038 
2039  return
2040  end subroutine atmos_dyn_fvm_fluxz_xvz_cd8
2041 
2042  !-----------------------------------------------------------------------------
2044  subroutine atmos_dyn_fvm_fluxj13_xvz_cd8( &
2045  flux, &
2046  mom, val, DENS, &
2047  GSQRT, J13G, MAPF, &
2048  CDZ, TwoD, &
2049  IIS, IIE, JJS, JJE )
2050  implicit none
2051 
2052  real(rp), intent(inout) :: flux (ka,ia,ja)
2053  real(rp), intent(in) :: mom (ka,ia,ja)
2054  real(rp), intent(in) :: val (ka,ia,ja)
2055  real(rp), intent(in) :: dens (ka,ia,ja)
2056  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2057  real(rp), intent(in) :: j13g (ka,ia,ja)
2058  real(rp), intent(in) :: mapf ( ia,ja,2)
2059  real(rp), intent(in) :: cdz (ka)
2060  logical, intent(in) :: twod
2061  integer, intent(in) :: iis, iie, jjs, jje
2062 
2063  real(rp) :: vel
2064  integer :: k, i, j
2065  !---------------------------------------------------------------------------
2066 
2067  !$omp parallel default(none) private(i,j,k,vel) &
2068  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J13G,MAPF) &
2069  !$omp shared(GSQRT,CDZ,TwoD)
2070 
2071 
2072 
2073  !$omp do OMP_SCHEDULE_ collapse(2)
2074  do j = jjs, jje
2075  do i = iis, iie
2076  do k = ks+3, ke-4
2077  vel = ( f2h(k,1,i_xvz) &
2078  * 0.25_rp * ( mom(k+1,i,j)+mom(k+1,i-1,j)+mom(k+1,i,j+1)+mom(k+1,i-1,j+1) ) &
2079  + f2h(k,2,i_xvz) &
2080  * 0.25_rp * ( mom(k,i,j)+mom(k,i-1,j)+mom(k,i,j+1)+mom(k,i-1,j+1) ) ) &
2081  / ( f2h(k,1,i_xvz) &
2082  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
2083  + f2h(k,2,i_xvz) &
2084  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
2085  vel = vel * j13g(k,i,j)
2086  flux(k,i,j) = vel / mapf(i,j,+2) &
2087  * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
2088  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
2089  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
2090  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) )
2091  enddo
2092  enddo
2093  enddo
2094  !$omp end do nowait
2095 
2096  !$omp do OMP_SCHEDULE_ collapse(2)
2097  do j = jjs, jje
2098  do i = iis, iie
2099  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
2100  ! The flux at KS-1 can be non-zero.
2101  ! To reduce calculations, all the fluxes are set to zero.
2102  flux(ks-1,i,j) = 0.0_rp
2103 
2104  vel = ( f2h(ks,1,i_xvz) &
2105  * 0.25_rp * ( mom(ks+1,i,j)+mom(ks+1,i-1,j)+mom(ks+1,i,j+1)+mom(ks+1,i-1,j+1) ) &
2106  + f2h(ks,2,i_xvz) &
2107  * 0.25_rp * ( mom(ks,i,j)+mom(ks,i-1,j)+mom(ks,i,j+1)+mom(ks,i-1,j+1) ) ) &
2108  / ( f2h(ks,1,i_xvz) &
2109  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
2110  + f2h(ks,2,i_xvz) &
2111  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
2112  vel = vel * j13g(ks,i,j)
2113  flux(ks,i,j) = vel / mapf(i,j,+2) &
2114  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
2115 
2116  vel = ( f2h(ke-1,1,i_xvz) &
2117  * 0.25_rp * ( mom(ke,i,j)+mom(ke,i-1,j)+mom(ke,i,j+1)+mom(ke,i-1,j+1) ) &
2118  + f2h(ke-1,2,i_xvz) &
2119  * 0.25_rp * ( mom(ke-1,i,j)+mom(ke-1,i-1,j)+mom(ke-1,i,j+1)+mom(ke-1,i-1,j+1) ) ) &
2120  / ( f2h(ke-1,1,i_xvz) &
2121  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
2122  + f2h(ke-1,2,i_xvz) &
2123  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
2124  vel = vel * j13g(ke-1,i,j)
2125  flux(ke-1,i,j) = vel / mapf(i,j,+2) &
2126  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
2127 
2128  vel = ( f2h(ks+1,1,i_xvz) &
2129  * 0.25_rp * ( mom(ks+2,i,j)+mom(ks+2,i-1,j)+mom(ks+2,i,j+1)+mom(ks+2,i-1,j+1) ) &
2130  + f2h(ks+1,2,i_xvz) &
2131  * 0.25_rp * ( mom(ks+1,i,j)+mom(ks+1,i-1,j)+mom(ks+1,i,j+1)+mom(ks+1,i-1,j+1) ) ) &
2132  / ( f2h(ks+1,1,i_xvz) &
2133  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) &
2134  + f2h(ks+1,2,i_xvz) &
2135  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) )
2136  vel = vel * j13g(ks+1,i,j)
2137  flux(ks+1,i,j) = vel / mapf(i,j,+2) &
2138  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
2139  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) )
2140 
2141  vel = ( f2h(ke-2,1,i_xvz) &
2142  * 0.25_rp * ( mom(ke-1,i,j)+mom(ke-1,i-1,j)+mom(ke-1,i,j+1)+mom(ke-1,i-1,j+1) ) &
2143  + f2h(ke-2,2,i_xvz) &
2144  * 0.25_rp * ( mom(ke-2,i,j)+mom(ke-2,i-1,j)+mom(ke-2,i,j+1)+mom(ke-2,i-1,j+1) ) ) &
2145  / ( f2h(ke-2,1,i_xvz) &
2146  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) &
2147  + f2h(ke-2,2,i_xvz) &
2148  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) )
2149  vel = vel * j13g(ke-2,i,j)
2150  flux(ke-2,i,j) = vel / mapf(i,j,+2) &
2151  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
2152  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) )
2153 
2154  vel = ( f2h(ks+2,1,i_xvz) &
2155  * 0.25_rp * ( mom(ks+3,i,j)+mom(ks+3,i-1,j)+mom(ks+3,i,j+1)+mom(ks+3,i-1,j+1) ) &
2156  + f2h(ks+2,2,i_xvz) &
2157  * 0.25_rp * ( mom(ks+2,i,j)+mom(ks+2,i-1,j)+mom(ks+2,i,j+1)+mom(ks+2,i-1,j+1) ) ) &
2158  / ( f2h(ks+2,1,i_xvz) &
2159  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i,j+1) ) &
2160  + f2h(ks+2,2,i_xvz) &
2161  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) )
2162  vel = vel * j13g(ks+2,i,j)
2163  flux(ks+2,i,j) = vel / mapf(i,j,+2) &
2164  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
2165  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
2166  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) )
2167 
2168  vel = ( f2h(ke-3,1,i_xvz) &
2169  * 0.25_rp * ( mom(ke-2,i,j)+mom(ke-2,i-1,j)+mom(ke-2,i,j+1)+mom(ke-2,i-1,j+1) ) &
2170  + f2h(ke-3,2,i_xvz) &
2171  * 0.25_rp * ( mom(ke-3,i,j)+mom(ke-3,i-1,j)+mom(ke-3,i,j+1)+mom(ke-3,i-1,j+1) ) ) &
2172  / ( f2h(ke-3,1,i_xvz) &
2173  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) &
2174  + f2h(ke-3,2,i_xvz) &
2175  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i,j+1) ) )
2176  vel = vel * j13g(ke-3,i,j)
2177  flux(ke-3,i,j) = vel / mapf(i,j,+2) &
2178  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
2179  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
2180  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) )
2181 
2182  flux(ke ,i,j) = 0.0_rp
2183  enddo
2184  enddo
2185  !$omp end do nowait
2186 
2187 
2188 
2189  !$omp end parallel
2190  return
2191  end subroutine atmos_dyn_fvm_fluxj13_xvz_cd8
2192 
2193  !-----------------------------------------------------------------------------
2195  subroutine atmos_dyn_fvm_fluxj23_xvz_cd8( &
2196  flux, &
2197  mom, val, DENS, &
2198  GSQRT, J23G, MAPF, &
2199  CDZ, TwoD, &
2200  IIS, IIE, JJS, JJE )
2201  implicit none
2202 
2203  real(rp), intent(inout) :: flux (ka,ia,ja)
2204  real(rp), intent(in) :: mom (ka,ia,ja)
2205  real(rp), intent(in) :: val (ka,ia,ja)
2206  real(rp), intent(in) :: dens (ka,ia,ja)
2207  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2208  real(rp), intent(in) :: j23g (ka,ia,ja)
2209  real(rp), intent(in) :: mapf ( ia,ja,2)
2210  real(rp), intent(in) :: cdz (ka)
2211  logical, intent(in) :: twod
2212  integer, intent(in) :: iis, iie, jjs, jje
2213 
2214  real(rp) :: vel
2215  integer :: k, i, j
2216  !---------------------------------------------------------------------------
2217 
2218  !$omp parallel default(none) private(i,j,k,vel) &
2219  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J23G,MAPF) &
2220  !$omp shared(GSQRT,CDZ,TwoD)
2221 
2222 
2223 
2224  !$omp do OMP_SCHEDULE_ collapse(2)
2225  do j = jjs, jje
2226  do i = iis, iie
2227  do k = ks+3, ke-4
2228  vel = ( f2h(k,1,i_xvz) &
2229  * mom(k+1,i,j) &
2230  + f2h(k,2,i_xvz) &
2231  * mom(k,i,j) ) &
2232  / ( f2h(k,1,i_xvz) &
2233  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
2234  + f2h(k,2,i_xvz) &
2235  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
2236  vel = vel * j23g(k,i,j)
2237  flux(k,i,j) = vel / mapf(i,j,+1) &
2238  * ( f81 * ( val(k+1,i,j)+val(k,i,j) ) &
2239  + f82 * ( val(k+2,i,j)+val(k-1,i,j) ) &
2240  + f83 * ( val(k+3,i,j)+val(k-2,i,j) ) &
2241  + f84 * ( val(k+4,i,j)+val(k-3,i,j) ) )
2242  enddo
2243  enddo
2244  enddo
2245  !$omp end do nowait
2246 
2247  !$omp do OMP_SCHEDULE_ collapse(2)
2248  do j = jjs, jje
2249  do i = iis, iie
2250  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
2251  ! The flux at KS-1 can be non-zero.
2252  ! To reduce calculations, all the fluxes are set to zero.
2253  flux(ks-1,i,j) = 0.0_rp
2254 
2255  vel = ( f2h(ks,1,i_xvz) &
2256  * mom(ks+1,i,j) &
2257  + f2h(ks,2,i_xvz) &
2258  * mom(ks,i,j) ) &
2259  / ( f2h(ks,1,i_xvz) &
2260  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
2261  + f2h(ks,2,i_xvz) &
2262  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
2263  vel = vel * j23g(ks,i,j)
2264  flux(ks,i,j) = vel / mapf(i,j,+1) &
2265  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
2266 
2267  vel = ( f2h(ke-1,1,i_xvz) &
2268  * mom(ke,i,j) &
2269  + f2h(ke-1,2,i_xvz) &
2270  * mom(ke-1,i,j) ) &
2271  / ( f2h(ke-1,1,i_xvz) &
2272  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
2273  + f2h(ke-1,2,i_xvz) &
2274  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
2275  vel = vel * j23g(ke-1,i,j)
2276  flux(ke-1,i,j) = vel / mapf(i,j,+1) &
2277  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
2278 
2279  vel = ( f2h(ks+1,1,i_xvz) &
2280  * mom(ks+2,i,j) &
2281  + f2h(ks+1,2,i_xvz) &
2282  * mom(ks+1,i,j) ) &
2283  / ( f2h(ks+1,1,i_xvz) &
2284  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) &
2285  + f2h(ks+1,2,i_xvz) &
2286  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) )
2287  vel = vel * j23g(ks+1,i,j)
2288  flux(ks+1,i,j) = vel / mapf(i,j,+1) &
2289  * ( f41 * ( val(ks+2,i,j)+val(ks+1,i,j) ) &
2290  + f42 * ( val(ks+3,i,j)+val(ks,i,j) ) )
2291 
2292  vel = ( f2h(ke-2,1,i_xvz) &
2293  * mom(ke-1,i,j) &
2294  + f2h(ke-2,2,i_xvz) &
2295  * mom(ke-2,i,j) ) &
2296  / ( f2h(ke-2,1,i_xvz) &
2297  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) &
2298  + f2h(ke-2,2,i_xvz) &
2299  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) )
2300  vel = vel * j23g(ke-2,i,j)
2301  flux(ke-2,i,j) = vel / mapf(i,j,+1) &
2302  * ( f41 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
2303  + f42 * ( val(ke,i,j)+val(ke-3,i,j) ) )
2304 
2305  vel = ( f2h(ks+2,1,i_xvz) &
2306  * mom(ks+3,i,j) &
2307  + f2h(ks+2,2,i_xvz) &
2308  * mom(ks+2,i,j) ) &
2309  / ( f2h(ks+2,1,i_xvz) &
2310  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i,j+1) ) &
2311  + f2h(ks+2,2,i_xvz) &
2312  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) )
2313  vel = vel * j23g(ks+2,i,j)
2314  flux(ks+2,i,j) = vel / mapf(i,j,+1) &
2315  * ( f61 * ( val(ks+3,i,j)+val(ks+2,i,j) ) &
2316  + f62 * ( val(ks+4,i,j)+val(ks+1,i,j) ) &
2317  + f63 * ( val(ks+5,i,j)+val(ks,i,j) ) )
2318 
2319  vel = ( f2h(ke-3,1,i_xvz) &
2320  * mom(ke-2,i,j) &
2321  + f2h(ke-3,2,i_xvz) &
2322  * mom(ke-3,i,j) ) &
2323  / ( f2h(ke-3,1,i_xvz) &
2324  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) &
2325  + f2h(ke-3,2,i_xvz) &
2326  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i,j+1) ) )
2327  vel = vel * j23g(ke-3,i,j)
2328  flux(ke-3,i,j) = vel / mapf(i,j,+1) &
2329  * ( f61 * ( val(ke-2,i,j)+val(ke-3,i,j) ) &
2330  + f62 * ( val(ke-1,i,j)+val(ke-4,i,j) ) &
2331  + f63 * ( val(ke,i,j)+val(ke-5,i,j) ) )
2332 
2333  flux(ke ,i,j) = 0.0_rp
2334  enddo
2335  enddo
2336  !$omp end do nowait
2337 
2338 
2339 
2340  !$omp end parallel
2341  return
2342  end subroutine atmos_dyn_fvm_fluxj23_xvz_cd8
2343 
2344  !-----------------------------------------------------------------------------
2346  subroutine atmos_dyn_fvm_fluxx_xvz_cd8( &
2347  flux, &
2348  mom, val, DENS, &
2349  GSQRT, MAPF, &
2350  num_diff, &
2351  CDZ, TwoD, &
2352  IIS, IIE, JJS, JJE )
2353  implicit none
2354 
2355  real(rp), intent(inout) :: flux (ka,ia,ja)
2356  real(rp), intent(in) :: mom (ka,ia,ja)
2357  real(rp), intent(in) :: val (ka,ia,ja)
2358  real(rp), intent(in) :: dens (ka,ia,ja)
2359  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2360  real(rp), intent(in) :: mapf ( ia,ja,2)
2361  real(rp), intent(in) :: num_diff(ka,ia,ja)
2362  real(rp), intent(in) :: cdz (ka)
2363  logical, intent(in) :: twod
2364  integer, intent(in) :: iis, iie, jjs, jje
2365 
2366  real(rp) :: vel
2367  integer :: k, i, j
2368  !---------------------------------------------------------------------------
2369 
2370 
2371 
2372  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
2373  !$omp private(vel) &
2374  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff)
2375  do j = jjs, jje
2376  do i = iis-1, iie
2377  do k = ks, ke
2378 #ifdef DEBUG
2379  call check( __line__, mom(k,i ,j) )
2380  call check( __line__, mom(k,i,j-1) )
2381 
2382  call check( __line__, val(k,i,j) )
2383  call check( __line__, val(k,i+1,j) )
2384 
2385  call check( __line__, val(k,i-1,j) )
2386  call check( __line__, val(k,i+2,j) )
2387 
2388  call check( __line__, val(k,i-2,j) )
2389  call check( __line__, val(k,i+3,j) )
2390 
2391  call check( __line__, val(k,i-3,j) )
2392  call check( __line__, val(k,i+4,j) )
2393 
2394 #endif
2395  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j+1) ) ) &
2396  / ( 0.25_rp * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) )
2397  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
2398  * ( f81 * ( val(k,i+1,j)+val(k,i,j) ) &
2399  + f82 * ( val(k,i+2,j)+val(k,i-1,j) ) &
2400  + f83 * ( val(k,i+3,j)+val(k,i-2,j) ) &
2401  + f84 * ( val(k,i+4,j)+val(k,i-3,j) ) ) &
2402  + gsqrt(k,i,j) * num_diff(k,i,j)
2403  enddo
2404  enddo
2405  enddo
2406 #ifdef DEBUG
2407  k = iundef; i = iundef; j = iundef
2408 #endif
2409 
2410 
2411 
2412  return
2413  end subroutine atmos_dyn_fvm_fluxx_xvz_cd8
2414 
2415  !-----------------------------------------------------------------------------
2417  subroutine atmos_dyn_fvm_fluxy_xvz_cd8( &
2418  flux, &
2419  mom, val, DENS, &
2420  GSQRT, MAPF, &
2421  num_diff, &
2422  CDZ, TwoD, &
2423  IIS, IIE, JJS, JJE )
2424  implicit none
2425 
2426  real(rp), intent(inout) :: flux (ka,ia,ja)
2427  real(rp), intent(in) :: mom (ka,ia,ja)
2428  real(rp), intent(in) :: val (ka,ia,ja)
2429  real(rp), intent(in) :: dens (ka,ia,ja)
2430  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2431  real(rp), intent(in) :: mapf ( ia,ja,2)
2432  real(rp), intent(in) :: num_diff(ka,ia,ja)
2433  real(rp), intent(in) :: cdz (ka)
2434  logical, intent(in) :: twod
2435  integer, intent(in) :: iis, iie, jjs, jje
2436 
2437  real(rp) :: vel
2438  integer :: k, i, j
2439  !---------------------------------------------------------------------------
2440 
2441  ! note that y-index is added by -1
2442 
2443 
2444 
2445  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
2446  !$omp private(vel) &
2447  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff)
2448  do j = jjs, jje+1
2449  do i = iis, iie
2450  do k = ks, ke
2451 #ifdef DEBUG
2452  call check( __line__, mom(k,i ,j) )
2453  call check( __line__, mom(k,i,j-1) )
2454 
2455  call check( __line__, val(k,i,j-1) )
2456  call check( __line__, val(k,i,j) )
2457 
2458  call check( __line__, val(k,i,j-2) )
2459  call check( __line__, val(k,i,j+1) )
2460 
2461  call check( __line__, val(k,i,j-3) )
2462  call check( __line__, val(k,i,j+2) )
2463 
2464  call check( __line__, val(k,i,j-4) )
2465  call check( __line__, val(k,i,j+3) )
2466 
2467 #endif
2468  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j-1) ) ) &
2469  / ( dens(k,i,j) )
2470  flux(k,i,j-1) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
2471  * ( f81 * ( val(k,i,j)+val(k,i,j-1) ) &
2472  + f82 * ( val(k,i,j+1)+val(k,i,j-2) ) &
2473  + f83 * ( val(k,i,j+2)+val(k,i,j-3) ) &
2474  + f84 * ( val(k,i,j+3)+val(k,i,j-4) ) ) &
2475  + gsqrt(k,i,j) * num_diff(k,i,j)
2476  enddo
2477  enddo
2478  enddo
2479 #ifdef DEBUG
2480  k = iundef; i = iundef; j = iundef
2481 #endif
2482 
2483 
2484 
2485  return
2486  end subroutine atmos_dyn_fvm_fluxy_xvz_cd8
2487 
2488 
2489 
2490 
2491 
2492 
2493 
2495 
2496 !--
2497 ! vi:set readonly sw=4 ts=8
2498 !
2499 !Local Variables:
2500 !mode: f90
2501 !buffer-read-only: t
2502 !End:
2503 !
2504 !++
scale_atmos_grid_cartesc_index::ke
integer, public ke
end point of inner domain: z, local
Definition: scale_atmos_grid_cartesC_index.F90:52
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxz_uyz_cd8
subroutine, public atmos_dyn_fvm_fluxz_uyz_cd8(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation z-flux at UY
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:978
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxj13_xvz_cd8
subroutine, public atmos_dyn_fvm_fluxj13_xvz_cd8(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J13-flux at XVZ
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:2050
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxx_xyz_cd8
subroutine, public atmos_dyn_fvm_fluxx_xyz_cd8(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation X-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:337
scale_index
module Index
Definition: scale_index.F90:11
scale_const::const_undef2
integer, parameter, public const_undef2
undefined value (INT2)
Definition: scale_const.F90:38
scale_precision
module PRECISION
Definition: scale_precision.F90:14
scale_atmos_grid_cartesc_index::ka
integer, public ka
Definition: scale_atmos_grid_cartesC_index.F90:47
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxj23_xyw_cd8
subroutine, public atmos_dyn_fvm_fluxj23_xyw_cd8(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J23-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:704
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxx_uyz_cd8
subroutine, public atmos_dyn_fvm_fluxx_uyz_cd8(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation X-flux at UY
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:1694
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:33
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxz_xyz_cd8
subroutine, public atmos_dyn_fvm_fluxz_xyz_cd8(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation z-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:200
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxz_xvz_cd8
subroutine, public atmos_dyn_fvm_fluxz_xvz_cd8(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation z-flux at XV
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:1883
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxj13_uyz_cd8
subroutine, public atmos_dyn_fvm_fluxj13_uyz_cd8(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J13-flux at UYZ
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:1274
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxy_xyw_cd8
subroutine, public atmos_dyn_fvm_fluxy_xyw_cd8(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation Y-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:888
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_flux_valuew_z_cd8
subroutine, public atmos_dyn_fvm_flux_valuew_z_cd8(valW, mflx, val, GSQRT, CDZ)
value at XYW
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:111
scale_atmos_grid_cartesc_index::i_xyz
integer, public i_xyz
Definition: scale_atmos_grid_cartesC_index.F90:90
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxx_xyw_cd8
subroutine, public atmos_dyn_fvm_fluxx_xyw_cd8(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation X-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:799
scale_atmos_grid_cartesc_index::i_uyz
integer, public i_uyz
Definition: scale_atmos_grid_cartesC_index.F90:94
scale_io
module STDIO
Definition: scale_io.F90:10
scale_tracer::k
real(rp), public k
Definition: scale_tracer.F90:44
scale_atmos_grid_cartesc_index
module atmosphere / grid / cartesC index
Definition: scale_atmos_grid_cartesC_index.F90:12
scale_const
module CONSTANT
Definition: scale_const.F90:11
scale_atmos_grid_cartesc_index::ia
integer, public ia
Definition: scale_atmos_grid_cartesC_index.F90:48
scale_debug::check
subroutine, public check(current_line, v)
Undefined value checker.
Definition: scale_debug.F90:56
scale_atmos_grid_cartesc_index::i_xvz
integer, public i_xvz
Definition: scale_atmos_grid_cartesC_index.F90:95
scale_prof
module profiler
Definition: scale_prof.F90:11
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxy_xvz_cd8
subroutine, public atmos_dyn_fvm_fluxy_xvz_cd8(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation Y-flux at XV
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:2424
scale_atmos_grid_cartesc_index::ja
integer, public ja
Definition: scale_atmos_grid_cartesC_index.F90:49
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxx_xvz_cd8
subroutine, public atmos_dyn_fvm_fluxx_xvz_cd8(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation X-flux at XV
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:2353
scale_tracer
module TRACER
Definition: scale_tracer.F90:12
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxj23_uyz_cd8
subroutine, public atmos_dyn_fvm_fluxj23_uyz_cd8(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J23-flux at UYZ
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:1425
scale_atmos_grid_cartesc_index::ks
integer, public ks
start point of inner domain: z, local
Definition: scale_atmos_grid_cartesC_index.F90:51
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxz_xyw_cd8
subroutine, public atmos_dyn_fvm_fluxz_xyw_cd8(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, FDZ, dtrk, IIS, IIE, JJS, JJE)
calculation z-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:462
scale_debug
module DEBUG
Definition: scale_debug.F90:11
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxy_uyz_cd8
subroutine, public atmos_dyn_fvm_fluxy_uyz_cd8(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation Y-flux at UY
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:1767
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxy_xyz_cd8
subroutine, public atmos_dyn_fvm_fluxy_xyz_cd8(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation Y-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:398
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:41
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxj13_xyw_cd8
subroutine, public atmos_dyn_fvm_fluxj13_xyw_cd8(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J13-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:611
scale_atmos_dyn_fvm_flux_cd8
module scale_atmos_dyn_fvm_flux_cd8
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:16
scale_atmos_dyn_fvm_flux_cd8::atmos_dyn_fvm_fluxj23_xvz_cd8
subroutine, public atmos_dyn_fvm_fluxj23_xvz_cd8(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J23-flux at XVZ
Definition: scale_atmos_dyn_fvm_flux_cd8.F90:2201