SCALE-RM
scale_atmos_dyn_fvm_flux_ud7.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"
16 module scale_atmos_dyn_fvm_flux_ud7
17  !-----------------------------------------------------------------------------
18  !
19  !++ used modules
20  !
21  use scale_precision
22  use scale_io
23  use scale_prof
24  use scale_atmos_grid_cartesc_index
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  !
42  public :: atmos_dyn_fvm_flux_valuew_z_ud7
43 
44  public :: atmos_dyn_fvm_fluxz_xyz_ud7
45  public :: atmos_dyn_fvm_fluxx_xyz_ud7
46  public :: atmos_dyn_fvm_fluxy_xyz_ud7
47 
48  public :: atmos_dyn_fvm_fluxz_xyw_ud7
49  public :: atmos_dyn_fvm_fluxj13_xyw_ud7
50  public :: atmos_dyn_fvm_fluxj23_xyw_ud7
51  public :: atmos_dyn_fvm_fluxx_xyw_ud7
52  public :: atmos_dyn_fvm_fluxy_xyw_ud7
53 
54  public :: atmos_dyn_fvm_fluxz_uyz_ud7
55  public :: atmos_dyn_fvm_fluxj13_uyz_ud7
56  public :: atmos_dyn_fvm_fluxj23_uyz_ud7
57  public :: atmos_dyn_fvm_fluxx_uyz_ud7
58  public :: atmos_dyn_fvm_fluxy_uyz_ud7
59 
60  public :: atmos_dyn_fvm_fluxz_xvz_ud7
61  public :: atmos_dyn_fvm_fluxj13_xvz_ud7
62  public :: atmos_dyn_fvm_fluxj23_xvz_ud7
63  public :: atmos_dyn_fvm_fluxx_xvz_ud7
64  public :: atmos_dyn_fvm_fluxy_xvz_ud7
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 :: F31 = -1.0_rp/12.0_rp
88  real(RP), parameter :: F32 = 7.0_rp/12.0_rp
89  real(RP), parameter :: F33 = 3.0_rp/12.0_rp
90 
91 
92  real(RP), parameter :: F51 = 1.0_rp/60.0_rp
93  real(RP), parameter :: F52 = -8.0_rp/60.0_rp
94  real(RP), parameter :: F53 = 37.0_rp/60.0_rp
95  real(RP), parameter :: F54 = -5.0_rp/60.0_rp
96  real(RP), parameter :: F55 = 10.0_rp/60.0_rp
97 
98 
99  real(RP), parameter :: F71 = -3.0_rp/840.0_rp
100  real(RP), parameter :: F72 = 29.0_rp/840.0_rp
101  real(RP), parameter :: F73 = -139.0_rp/840.0_rp
102  real(RP), parameter :: F74 = 533.0_rp/840.0_rp
103  real(RP), parameter :: F75 = 21.0_rp/840.0_rp
104  real(RP), parameter :: F76 = -63.0_rp/840.0_rp
105  real(RP), parameter :: F77 = 105.0_rp/840.0_rp
106 
107 
108 
109 contains
110  !-----------------------------------------------------------------------------
112 !OCL SERIAL
113  subroutine atmos_dyn_fvm_flux_valuew_z_ud7( &
114  valW, &
115  mflx, val, GSQRT, &
116  CDZ )
117  !$acc routine vector
118  implicit none
119 
120  real(rp), intent(out) :: valw (ka)
121  real(rp), intent(in) :: mflx (ka)
122  real(rp), intent(in) :: val (ka)
123  real(rp), intent(in) :: gsqrt(ka)
124  real(rp), intent(in) :: cdz (ka)
125 
126  integer :: k
127  !---------------------------------------------------------------------------
128 
129  do k = ks+3, ke-4
130 #ifdef DEBUG
131  call check( __line__, mflx(k) )
132 
133  call check( __line__, val(k) )
134  call check( __line__, val(k+1) )
135 
136  call check( __line__, val(k-1) )
137  call check( __line__, val(k+2) )
138 
139  call check( __line__, val(k-2) )
140  call check( __line__, val(k+3) )
141 
142  call check( __line__, val(k-3) )
143  call check( __line__, val(k+4) )
144 
145 #endif
146  valw(k) = ( f71 * ( val(k+4)+val(k-3) ) &
147  + f72 * ( val(k+3)+val(k-2) ) &
148  + f73 * ( val(k+2)+val(k-1) ) &
149  + f74 * ( val(k+1)+val(k) ) ) &
150  - ( f71 * ( val(k+4)-val(k-3) ) &
151  + f75 * ( val(k+3)-val(k-2) ) &
152  + f76 * ( val(k+2)-val(k-1) ) &
153  + f77 * ( val(k+1)-val(k) ) ) * sign(1.0_rp,mflx(k))
154  enddo
155 #ifdef DEBUG
156  k = iundef
157 #endif
158 
159 #ifdef DEBUG
160 
161  call check( __line__, mflx(ks) )
162  call check( __line__, val(ks ) )
163  call check( __line__, val(ks+1) )
164  call check( __line__, mflx(ke-1) )
165  call check( __line__, val(ke ) )
166  call check( __line__, val(ke-1) )
167 
168  call check( __line__, mflx(ks+1) )
169  call check( __line__, val(ks+2 ) )
170  call check( __line__, val(ks+3) )
171  call check( __line__, mflx(ke-2) )
172  call check( __line__, val(ke-2 ) )
173  call check( __line__, val(ke-3) )
174 
175  call check( __line__, mflx(ks+2) )
176  call check( __line__, val(ks+4 ) )
177  call check( __line__, val(ks+5) )
178  call check( __line__, mflx(ke-3) )
179  call check( __line__, val(ke-4 ) )
180  call check( __line__, val(ke-5) )
181 
182 #endif
183 
184  valw(ks) = f2 * ( val(ks+1)+val(ks) ) &
185  * ( 0.5_rp + sign(0.5_rp,mflx(ks)) ) &
186  + ( 2.0_rp * val(ks) + 5.0_rp * val(ks+1) - val(ks+2) ) / 6.0_rp &
187  * ( 0.5_rp - sign(0.5_rp,mflx(ks)) )
188  valw(ke-1) = ( 2.0_rp * val(ke) + 5.0_rp * val(ke-1) - val(ke-2) ) / 6.0_rp &
189  * ( 0.5_rp + sign(0.5_rp,mflx(ke-1)) ) &
190  + f2 * ( val(ke)+val(ke-1) ) &
191  * ( 0.5_rp - sign(0.5_rp,mflx(ke-1)) )
192 
193  valw(ks+1) = ( 2.0_rp * val(ks+2) + 5.0_rp * val(ks+1) - val(ks) ) / 6.0_rp &
194  * ( 0.5_rp + sign(0.5_rp,mflx(ks+1)) ) &
195  + ( - 3.0_rp * val(ks) &
196  + 27.0_rp * val(ks+1) &
197  + 47.0_rp * val(ks+2) &
198  - 13.0_rp * val(ks+3) &
199  + 2.0_rp * val(ks+4) ) / 60.0_rp &
200  * ( 0.5_rp - sign(0.5_rp,mflx(ks+1)) )
201  valw(ke-2) = ( - 3.0_rp * val(ke) &
202  + 27.0_rp * val(ke-1) &
203  + 47.0_rp * val(ke-2) &
204  - 13.0_rp * val(ke-3) &
205  + 2.0_rp * val(ke-4) ) / 60.0_rp &
206  * ( 0.5_rp + sign(0.5_rp,mflx(ke-2)) ) &
207  + ( 2.0_rp * val(ke-2) + 5.0_rp * val(ke-1) - val(ke) ) / 6.0_rp &
208  * ( 0.5_rp - sign(0.5_rp,mflx(ke-2)) )
209 
210  valw(ks+2) = ( - 3.0_rp * val(ks+4) &
211  + 27.0_rp * val(ks+3) &
212  + 47.0_rp * val(ks+2) &
213  - 13.0_rp * val(ks+1) &
214  + 2.0_rp * val(ks) ) / 60.0_rp &
215  * ( 0.5_rp + sign(0.5_rp,mflx(ks+2)) ) &
216  + ( 4.0_rp * val(ks) &
217  - 38.0_rp * val(ks+1) &
218  + 214.0_rp * val(ks+2) &
219  + 319.0_rp * val(ks+3) &
220  - 101.0_rp * val(ks+4) &
221  + 25.0_rp * val(ks+5) &
222  - 3.0_rp * val(ks+6) ) / 420.0_rp &
223  * ( 0.5_rp - sign(0.5_rp,mflx(ks+2)) )
224  valw(ke-3) = ( 4.0_rp * val(ke) &
225  - 38.0_rp * val(ke-1) &
226  + 214.0_rp * val(ke-2) &
227  + 319.0_rp * val(ke-3) &
228  - 101.0_rp * val(ke-4) &
229  + 25.0_rp * val(ke-5) &
230  - 3.0_rp * val(ke-6) ) / 420.0_rp &
231  * ( 0.5_rp + sign(0.5_rp,mflx(ke-3)) ) &
232  + ( - 3.0_rp * val(ke-4) &
233  + 27.0_rp * val(ke-3) &
234  + 47.0_rp * val(ke-2) &
235  - 13.0_rp * val(ke-1) &
236  + 2.0_rp * val(ke) ) / 60.0_rp &
237  * ( 0.5_rp - sign(0.5_rp,mflx(ke-3)) )
238 
239 
240  return
241  end subroutine atmos_dyn_fvm_flux_valuew_z_ud7
242 
243  !-----------------------------------------------------------------------------
245  subroutine atmos_dyn_fvm_fluxz_xyz_ud7( &
246  flux, &
247  mflx, val, GSQRT, &
248  num_diff, &
249  CDZ, &
250  IIS, IIE, JJS, JJE )
251  use scale_const, only: &
252  eps => const_eps
253  implicit none
254 
255  real(rp), intent(inout) :: flux (ka,ia,ja)
256  real(rp), intent(in) :: mflx (ka,ia,ja)
257  real(rp), intent(in) :: val (ka,ia,ja)
258  real(rp), intent(in) :: gsqrt (ka,ia,ja)
259  real(rp), intent(in) :: num_diff(ka,ia,ja)
260  real(rp), intent(in) :: cdz (ka)
261  integer, intent(in) :: iis, iie, jjs, jje
262 
263  real(rp) :: vel
264  integer :: k, i, j
265  !---------------------------------------------------------------------------
266 
267  !$omp parallel default(none) private(i,j,k, vel) &
268  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mflx,val,flux,GSQRT,num_diff,EPS)
269 
270  !$acc data copy(flux) copyin(mflx, val, GSQRT, num_diff, CDZ)
271 
272  !$omp do OMP_SCHEDULE_ collapse(2)
273  !$acc kernels
274  do j = jjs, jje
275  do i = iis, iie
276  do k = ks+3, ke-4
277 #ifdef DEBUG
278  call check( __line__, mflx(k,i,j) )
279 
280  call check( __line__, val(k,i,j) )
281  call check( __line__, val(k+1,i,j) )
282 
283  call check( __line__, val(k-1,i,j) )
284  call check( __line__, val(k+2,i,j) )
285 
286  call check( __line__, val(k-2,i,j) )
287  call check( __line__, val(k+3,i,j) )
288 
289  call check( __line__, val(k-3,i,j) )
290  call check( __line__, val(k+4,i,j) )
291 
292 #endif
293  vel = mflx(k,i,j)
294  flux(k,i,j) = vel &
295  * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
296  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
297  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
298  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
299  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
300  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
301  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
302  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
303  + gsqrt(k,i,j) * num_diff(k,i,j)
304  enddo
305  enddo
306  enddo
307  !$acc end kernels
308  !$omp end do nowait
309 #ifdef DEBUG
310  k = iundef; i = iundef; j = iundef
311 #endif
312 
313  !$omp do OMP_SCHEDULE_ collapse(2)
314  !$acc kernels
315  do j = jjs, jje
316  do i = iis, iie
317 #ifdef DEBUG
318 
319  call check( __line__, mflx(ks,i,j) )
320  call check( __line__, val(ks ,i,j) )
321  call check( __line__, val(ks+1,i,j) )
322  call check( __line__, mflx(ke-1,i,j) )
323  call check( __line__, val(ke ,i,j) )
324  call check( __line__, val(ke-1,i,j) )
325 
326  call check( __line__, mflx(ks+1,i,j) )
327  call check( __line__, val(ks+2 ,i,j) )
328  call check( __line__, val(ks+3,i,j) )
329  call check( __line__, mflx(ke-2,i,j) )
330  call check( __line__, val(ke-2 ,i,j) )
331  call check( __line__, val(ke-3,i,j) )
332 
333  call check( __line__, mflx(ks+2,i,j) )
334  call check( __line__, val(ks+4 ,i,j) )
335  call check( __line__, val(ks+5,i,j) )
336  call check( __line__, mflx(ke-3,i,j) )
337  call check( __line__, val(ke-4 ,i,j) )
338  call check( __line__, val(ke-5,i,j) )
339 
340 #endif
341  flux(ks-1,i,j) = 0.0_rp
342 
343  vel = mflx(ks,i,j)
344  flux(ks,i,j) = vel &
345  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
346  * ( 0.5_rp + sign(0.5_rp,vel) ) &
347  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
348  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
349  + gsqrt(ks,i,j) * num_diff(ks,i,j)
350  vel = mflx(ke-1,i,j)
351  flux(ke-1,i,j) = vel &
352  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
353  * ( 0.5_rp + sign(0.5_rp,vel) ) &
354  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
355  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
356  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
357 
358  vel = mflx(ks+1,i,j)
359  flux(ks+1,i,j) = vel &
360  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
361  * ( 0.5_rp + sign(0.5_rp,vel) ) &
362  + ( - 3.0_rp * val(ks,i,j) &
363  + 27.0_rp * val(ks+1,i,j) &
364  + 47.0_rp * val(ks+2,i,j) &
365  - 13.0_rp * val(ks+3,i,j) &
366  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
367  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
368  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j)
369  vel = mflx(ke-2,i,j)
370  flux(ke-2,i,j) = vel &
371  * ( ( - 3.0_rp * val(ke,i,j) &
372  + 27.0_rp * val(ke-1,i,j) &
373  + 47.0_rp * val(ke-2,i,j) &
374  - 13.0_rp * val(ke-3,i,j) &
375  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
376  * ( 0.5_rp + sign(0.5_rp,vel) ) &
377  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
378  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
379  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j)
380 
381  vel = mflx(ks+2,i,j)
382  flux(ks+2,i,j) = vel &
383  * ( ( - 3.0_rp * val(ks+4,i,j) &
384  + 27.0_rp * val(ks+3,i,j) &
385  + 47.0_rp * val(ks+2,i,j) &
386  - 13.0_rp * val(ks+1,i,j) &
387  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
388  * ( 0.5_rp + sign(0.5_rp,vel) ) &
389  + ( 4.0_rp * val(ks,i,j) &
390  - 38.0_rp * val(ks+1,i,j) &
391  + 214.0_rp * val(ks+2,i,j) &
392  + 319.0_rp * val(ks+3,i,j) &
393  - 101.0_rp * val(ks+4,i,j) &
394  + 25.0_rp * val(ks+5,i,j) &
395  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
396  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
397  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j)
398  vel = mflx(ke-3,i,j)
399  flux(ke-3,i,j) = vel &
400  * ( ( 4.0_rp * val(ke,i,j) &
401  - 38.0_rp * val(ke-1,i,j) &
402  + 214.0_rp * val(ke-2,i,j) &
403  + 319.0_rp * val(ke-3,i,j) &
404  - 101.0_rp * val(ke-4,i,j) &
405  + 25.0_rp * val(ke-5,i,j) &
406  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
407  * ( 0.5_rp + sign(0.5_rp,vel) ) &
408  + ( - 3.0_rp * val(ke-4,i,j) &
409  + 27.0_rp * val(ke-3,i,j) &
410  + 47.0_rp * val(ke-2,i,j) &
411  - 13.0_rp * val(ke-1,i,j) &
412  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
413  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
414  + gsqrt(ke-3,i,j) * num_diff(ke-3,i,j)
415 
416  flux(ke ,i,j) = 0.0_rp
417  enddo
418  enddo
419  !$acc end kernels
420  !$omp end do nowait
421 
422  !$acc end data
423 
424  !$omp end parallel
425 #ifdef DEBUG
426  k = iundef; i = iundef; j = iundef
427 #endif
428 
429  return
430  end subroutine atmos_dyn_fvm_fluxz_xyz_ud7
431 
432  !-----------------------------------------------------------------------------
434  subroutine atmos_dyn_fvm_fluxx_xyz_ud7( &
435  flux, &
436  mflx, val, GSQRT, &
437  num_diff, &
438  CDZ, &
439  IIS, IIE, JJS, JJE )
440  implicit none
441 
442  real(rp), intent(inout) :: flux (ka,ia,ja)
443  real(rp), intent(in) :: mflx (ka,ia,ja)
444  real(rp), intent(in) :: val (ka,ia,ja)
445  real(rp), intent(in) :: gsqrt (ka,ia,ja)
446  real(rp), intent(in) :: num_diff(ka,ia,ja)
447  real(rp), intent(in) :: cdz(ka)
448  integer, intent(in) :: iis, iie, jjs, jje
449 
450  real(rp) :: vel
451  integer :: k, i, j
452  !---------------------------------------------------------------------------
453 
454  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
455  !$omp private(vel) &
456  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mflx,val,flux,GSQRT,num_diff)
457  !$acc kernels
458  do j = jjs, jje
459  do i = iis-1, iie
460  do k = ks, ke
461 #ifdef DEBUG
462  call check( __line__, mflx(k,i,j) )
463 
464  call check( __line__, val(k,i,j) )
465  call check( __line__, val(k,i+1,j) )
466 
467  call check( __line__, val(k,i-1,j) )
468  call check( __line__, val(k,i+2,j) )
469 
470  call check( __line__, val(k,i-2,j) )
471  call check( __line__, val(k,i+3,j) )
472 
473  call check( __line__, val(k,i-3,j) )
474  call check( __line__, val(k,i+4,j) )
475 
476 #endif
477  vel = mflx(k,i,j)
478  flux(k,i,j) = vel &
479  * ( ( f71 * ( val(k,i+4,j)+val(k,i-3,j) ) &
480  + f72 * ( val(k,i+3,j)+val(k,i-2,j) ) &
481  + f73 * ( val(k,i+2,j)+val(k,i-1,j) ) &
482  + f74 * ( val(k,i+1,j)+val(k,i,j) ) ) &
483  - ( f71 * ( val(k,i+4,j)-val(k,i-3,j) ) &
484  + f75 * ( val(k,i+3,j)-val(k,i-2,j) ) &
485  + f76 * ( val(k,i+2,j)-val(k,i-1,j) ) &
486  + f77 * ( val(k,i+1,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
487  + gsqrt(k,i,j) * num_diff(k,i,j)
488  enddo
489  enddo
490  enddo
491  !$acc end kernels
492 #ifdef DEBUG
493  k = iundef; i = iundef; j = iundef
494 #endif
495 
496  return
497  end subroutine atmos_dyn_fvm_fluxx_xyz_ud7
498 
499  !-----------------------------------------------------------------------------
501  subroutine atmos_dyn_fvm_fluxy_xyz_ud7( &
502  flux, &
503  mflx, val, GSQRT, &
504  num_diff, &
505  CDZ, &
506  IIS, IIE, JJS, JJE )
507  implicit none
508 
509  real(rp), intent(inout) :: flux (ka,ia,ja)
510  real(rp), intent(in) :: mflx (ka,ia,ja)
511  real(rp), intent(in) :: val (ka,ia,ja)
512  real(rp), intent(in) :: gsqrt (ka,ia,ja)
513  real(rp), intent(in) :: num_diff(ka,ia,ja)
514  real(rp), intent(in) :: cdz(ka)
515  integer, intent(in) :: iis, iie, jjs, jje
516 
517  real(rp) :: vel
518  integer :: k, i, j
519  !---------------------------------------------------------------------------
520 
521  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
522  !$omp private(vel) &
523  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mflx,val,flux,GSQRT,num_diff)
524  !$acc kernels
525  do j = jjs-1, jje
526  do i = iis, iie
527  do k = ks, ke
528 #ifdef DEBUG
529  call check( __line__, mflx(k,i,j) )
530 
531  call check( __line__, val(k,i,j) )
532  call check( __line__, val(k,i,j+1) )
533 
534  call check( __line__, val(k,i,j-1) )
535  call check( __line__, val(k,i,j+2) )
536 
537  call check( __line__, val(k,i,j-2) )
538  call check( __line__, val(k,i,j+3) )
539 
540  call check( __line__, val(k,i,j-3) )
541  call check( __line__, val(k,i,j+4) )
542 
543 #endif
544  vel = mflx(k,i,j)
545  flux(k,i,j) = vel &
546  * ( ( f71 * ( val(k,i,j+4)+val(k,i,j-3) ) &
547  + f72 * ( val(k,i,j+3)+val(k,i,j-2) ) &
548  + f73 * ( val(k,i,j+2)+val(k,i,j-1) ) &
549  + f74 * ( val(k,i,j+1)+val(k,i,j) ) ) &
550  - ( f71 * ( val(k,i,j+4)-val(k,i,j-3) ) &
551  + f75 * ( val(k,i,j+3)-val(k,i,j-2) ) &
552  + f76 * ( val(k,i,j+2)-val(k,i,j-1) ) &
553  + f77 * ( val(k,i,j+1)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
554  + gsqrt(k,i,j) * num_diff(k,i,j)
555  enddo
556  enddo
557  enddo
558  !$acc end kernels
559 #ifdef DEBUG
560  k = iundef; i = iundef; j = iundef
561 #endif
562 
563  return
564  end subroutine atmos_dyn_fvm_fluxy_xyz_ud7
565 
566 
567  !-----------------------------------------------------------------------------
569  subroutine atmos_dyn_fvm_fluxz_xyw_ud7( &
570  flux, &
571  mom, val, DENS, &
572  GSQRT, J33G, &
573  num_diff, &
574  CDZ, FDZ, &
575  dtrk, &
576  IIS, IIE, JJS, JJE )
577  implicit none
578 
579  real(rp), intent(inout) :: flux (ka,ia,ja)
580  real(rp), intent(in) :: mom (ka,ia,ja)
581  real(rp), intent(in) :: val (ka,ia,ja)
582  real(rp), intent(in) :: dens (ka,ia,ja)
583  real(rp), intent(in) :: gsqrt (ka,ia,ja)
584  real(rp), intent(in) :: j33g
585  real(rp), intent(in) :: num_diff(ka,ia,ja)
586  real(rp), intent(in) :: cdz (ka)
587  real(rp), intent(in) :: fdz (ka-1)
588  real(rp), intent(in) :: dtrk
589  integer, intent(in) :: iis, iie, jjs, jje
590 
591  real(rp) :: vel
592  integer :: k, i, j
593  !---------------------------------------------------------------------------
594 
595  ! note than z-index is added by -1
596 
597  !$omp parallel default(none) private(i,j,k,vel) &
598  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,flux,J33G,GSQRT,num_diff,DENS,FDZ,dtrk)
599 
600  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, num_diff, CDZ, FDZ)
601 
602  !$omp do OMP_SCHEDULE_ collapse(2)
603  !$acc kernels
604  do j = jjs, jje
605  do i = iis, iie
606  do k = ks+4, ke-3
607 #ifdef DEBUG
608  call check( __line__, mom(k-1,i,j) )
609  call check( __line__, mom(k ,i,j) )
610 
611  call check( __line__, val(k-1,i,j) )
612  call check( __line__, val(k,i,j) )
613 
614  call check( __line__, val(k-2,i,j) )
615  call check( __line__, val(k+1,i,j) )
616 
617  call check( __line__, val(k-3,i,j) )
618  call check( __line__, val(k+2,i,j) )
619 
620  call check( __line__, val(k-4,i,j) )
621  call check( __line__, val(k+3,i,j) )
622 
623 #endif
624  vel = ( 0.5_rp * ( mom(k-1,i,j) &
625  + mom(k,i,j) ) ) &
626  / dens(k,i,j)
627  flux(k-1,i,j) = j33g * vel &
628  * ( ( f71 * ( val(k+3,i,j)+val(k-4,i,j) ) &
629  + f72 * ( val(k+2,i,j)+val(k-3,i,j) ) &
630  + f73 * ( val(k+1,i,j)+val(k-2,i,j) ) &
631  + f74 * ( val(k,i,j)+val(k-1,i,j) ) ) &
632  - ( f71 * ( val(k+3,i,j)-val(k-4,i,j) ) &
633  + f75 * ( val(k+2,i,j)-val(k-3,i,j) ) &
634  + f76 * ( val(k+1,i,j)-val(k-2,i,j) ) &
635  + f77 * ( val(k,i,j)-val(k-1,i,j) ) ) * sign(1.0_rp,vel) ) &
636  + gsqrt(k,i,j) * num_diff(k,i,j)
637  enddo
638  enddo
639  enddo
640  !$acc end kernels
641  !$omp end do nowait
642 #ifdef DEBUG
643  k = iundef; i = iundef; j = iundef
644 #endif
645 
646  !$omp do OMP_SCHEDULE_ collapse(2)
647  !$acc kernels
648  do j = jjs, jje
649  do i = iis, iie
650 #ifdef DEBUG
651 
652  call check( __line__, val(ks,i,j) )
653  call check( __line__, val(ks+1,i,j) )
654  call check( __line__, val(ks+2,i,j) )
655  call check( __line__, val(ks+3,i,j) )
656  call check( __line__, val(ks+4,i,j) )
657  call check( __line__, val(ks+5,i,j) )
658  call check( __line__, val(ks+6,i,j) )
659 
660 
661  call check( __line__, val(ke-6,i,j) )
662  call check( __line__, val(ke-5,i,j) )
663  call check( __line__, val(ke-4,i,j) )
664  call check( __line__, val(ke-3,i,j) )
665  call check( __line__, val(ke-2,i,j) )
666  call check( __line__, val(ke-1,i,j) )
667 
668 #endif
669  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS.
670  ! The flux at KS can be non-zero.
671  ! To reduce calculations, all the fluxes are set to zero.
672  flux(ks-1,i,j) = 0.0_rp ! k = KS
673 
674  vel = ( 0.5_rp * ( mom(ks,i,j) &
675  + mom(ks+1,i,j) ) ) &
676  / dens(ks+1,i,j)
677  flux(ks,i,j) = j33g * vel &
678  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
679  * ( 0.5_rp + sign(0.5_rp,vel) ) &
680  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
681  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
682  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j) ! k = KS+1
683 
684  vel = ( 0.5_rp * ( mom(ks+1,i,j) &
685  + mom(ks+2,i,j) ) ) &
686  / dens(ks+2,i,j)
687  flux(ks+1,i,j) = j33g * vel &
688  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
689  * ( 0.5_rp + sign(0.5_rp,vel) ) &
690  + ( - 3.0_rp * val(ks,i,j) &
691  + 27.0_rp * val(ks+1,i,j) &
692  + 47.0_rp * val(ks+2,i,j) &
693  - 13.0_rp * val(ks+3,i,j) &
694  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
695  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
696  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j) ! k = KS+2
697 
698  vel = ( 0.5_rp * ( mom(ks+2,i,j) &
699  + mom(ks+3,i,j) ) ) &
700  / dens(ks+3,i,j)
701  flux(ks+2,i,j) = j33g * vel &
702  * ( ( - 3.0_rp * val(ks+4,i,j) &
703  + 27.0_rp * val(ks+3,i,j) &
704  + 47.0_rp * val(ks+2,i,j) &
705  - 13.0_rp * val(ks+1,i,j) &
706  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
707  * ( 0.5_rp + sign(0.5_rp,vel) ) &
708  + ( 4.0_rp * val(ks,i,j) &
709  - 38.0_rp * val(ks+1,i,j) &
710  + 214.0_rp * val(ks+2,i,j) &
711  + 319.0_rp * val(ks+3,i,j) &
712  - 101.0_rp * val(ks+4,i,j) &
713  + 25.0_rp * val(ks+5,i,j) &
714  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
715  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
716  + gsqrt(ks+3,i,j) * num_diff(ks+3,i,j) ! k = KS+3
717 
718 
719 
720  vel = ( 0.5_rp * ( mom(ke-3,i,j) &
721  + mom(ke-2,i,j) ) ) &
722  / dens(ke-1,i,j)
723  flux(ke-3,i,j) = j33g * vel &
724  * ( ( 4.0_rp * val(ke,i,j) &
725  - 38.0_rp * val(ke-1,i,j) &
726  + 214.0_rp * val(ke-2,i,j) &
727  + 319.0_rp * val(ke-3,i,j) &
728  - 101.0_rp * val(ke-4,i,j) &
729  + 25.0_rp * val(ke-5,i,j) &
730  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
731  * ( 0.5_rp + sign(0.5_rp,vel) ) &
732  + ( - 3.0_rp * val(ke-4,i,j) &
733  + 27.0_rp * val(ke-3,i,j) &
734  + 47.0_rp * val(ke-2,i,j) &
735  - 13.0_rp * val(ke-1,i,j) &
736  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
737  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
738  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j) ! k = KE-2
739 
740  vel = ( 0.5_rp * ( mom(ke-2,i,j) &
741  + mom(ke-1,i,j) ) ) &
742  / dens(ke-1,i,j)
743  flux(ke-2,i,j) = j33g * vel &
744  * ( ( - 3.0_rp * val(ke,i,j) &
745  + 27.0_rp * val(ke-1,i,j) &
746  + 47.0_rp * val(ke-2,i,j) &
747  - 13.0_rp * val(ke-3,i,j) &
748  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
749  * ( 0.5_rp + sign(0.5_rp,vel) ) &
750  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
751  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
752  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j) ! k = KE-1
753 
754  flux(ke-1,i,j) = 0.0_rp ! k = KE
755  flux(ke ,i,j) = 0.0_rp ! k = KE+1
756  enddo
757  enddo
758  !$acc end kernels
759  !$omp end do nowait
760 
761  !$acc end data
762 
763  !$omp end parallel
764 
765  return
766  end subroutine atmos_dyn_fvm_fluxz_xyw_ud7
767 
768 
769  !-----------------------------------------------------------------------------
771  subroutine atmos_dyn_fvm_fluxj13_xyw_ud7( &
772  flux, &
773  mom, val, DENS, &
774  GSQRT, J13G, MAPF, &
775  CDZ, TwoD, &
776  IIS, IIE, JJS, JJE )
777  implicit none
778 
779  real(rp), intent(inout) :: flux (ka,ia,ja)
780  real(rp), intent(in) :: mom (ka,ia,ja)
781  real(rp), intent(in) :: val (ka,ia,ja)
782  real(rp), intent(in) :: dens (ka,ia,ja)
783  real(rp), intent(in) :: gsqrt (ka,ia,ja)
784  real(rp), intent(in) :: j13g (ka,ia,ja)
785  real(rp), intent(in) :: mapf ( ia,ja,2)
786  real(rp), intent(in) :: cdz (ka)
787  logical, intent(in) :: twod
788  integer, intent(in) :: iis, iie, jjs, jje
789 
790  real(rp) :: vel
791  integer :: k, i, j
792  !---------------------------------------------------------------------------
793 
794  !$omp parallel default(none) private(i,j,k,vel) &
795  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J13G,MAPF)
796 
797  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, J13G, MAPF, CDZ)
798 
799  !$omp do OMP_SCHEDULE_ collapse(2)
800  !$acc kernels
801  do j = jjs, jje
802  do i = iis, iie
803  do k = ks+3, ke-3
804  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i-1,j) ) ) &
805  / dens(k,i,j)
806  vel = vel * j13g(k,i,j)
807  flux(k-1,i,j) = vel / mapf(i,j,+2) &
808  * ( ( f71 * ( val(k+3,i,j)+val(k-4,i,j) ) &
809  + f72 * ( val(k+2,i,j)+val(k-3,i,j) ) &
810  + f73 * ( val(k+1,i,j)+val(k-2,i,j) ) &
811  + f74 * ( val(k,i,j)+val(k-1,i,j) ) ) &
812  - ( f71 * ( val(k+3,i,j)-val(k-4,i,j) ) &
813  + f75 * ( val(k+2,i,j)-val(k-3,i,j) ) &
814  + f76 * ( val(k+1,i,j)-val(k-2,i,j) ) &
815  + f77 * ( val(k,i,j)-val(k-1,i,j) ) ) * sign(1.0_rp,vel) )
816  enddo
817  enddo
818  enddo
819  !$acc end kernels
820  !$omp end do nowait
821 
822  !$omp do OMP_SCHEDULE_ collapse(2)
823  !$acc kernels
824  do j = jjs, jje
825  do i = iis, iie
826  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS.
827  ! The flux at KS can be non-zero.
828  ! To reduce calculations, all the fluxes are set to zero.
829  flux(ks-1,i,j) = 0.0_rp ! k = KS
830 
831  ! physically incorrect but for numerical stability
832  vel = ( ( 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i-1,j) ) ) / dens(ks+1,i,j) &
833  + ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i-1,j) ) ) / dens(ks ,i,j) ) * 0.5_rp
834 ! vel = ( 0.5_RP * ( mom(KS+1,i,j)+mom(KS+1,i-1,j) ) ) &
835 ! / DENS(KS+1,i,j)
836  vel = vel * j13g(ks+1,i,j)
837  flux(ks,i,j) = vel / mapf(i,j,+2) &
838  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
839  * ( 0.5_rp + sign(0.5_rp,vel) ) &
840  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
841  * ( 0.5_rp - sign(0.5_rp,vel) ) ) ! k = KS+1
842 
843  vel = ( 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i-1,j) ) ) &
844  / dens(ks+2,i,j)
845  vel = vel * j13g(ks,i,j)
846  flux(ks+1,i,j) = vel / mapf(i,j,+2) &
847  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
848  * ( 0.5_rp + sign(0.5_rp,vel) ) &
849  + ( - 3.0_rp * val(ks,i,j) &
850  + 27.0_rp * val(ks+1,i,j) &
851  + 47.0_rp * val(ks+2,i,j) &
852  - 13.0_rp * val(ks+3,i,j) &
853  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
854  * ( 0.5_rp - sign(0.5_rp,vel) ) ) ! k = KS+2
855 
856 
857  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i-1,j) ) ) &
858  / dens(ke-1,i,j)
859  vel = vel * j13g(ke-1,i,j)
860  flux(ke-2,i,j) = vel / mapf(i,j,+2) &
861  * ( ( 2.0_rp * val(ke-1,i,j) + 5.0_rp * val(ke-2,i,j) - val(ke-3,i,j) ) / 6.0_rp &
862  * ( 0.5_rp + sign(0.5_rp,vel) ) &
863  + f2 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
864  * ( 0.5_rp - sign(0.5_rp,vel) ) ) ! k = KE-3
865 
866  vel = ( 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i-1,j) ) ) &
867  / dens(ke-2,i,j)
868  vel = vel * j13g(ke-2,i,j)
869  flux(ke-3,i,j) = vel / mapf(i,j,+2) &
870  * ( ( - 3.0_rp * val(ke-1,i,j) &
871  + 27.0_rp * val(ke-2,i,j) &
872  + 47.0_rp * val(ke-3,i,j) &
873  - 13.0_rp * val(ke-4,i,j) &
874  + 2.0_rp * val(ke-5,i,j) ) / 60.0_rp &
875  * ( 0.5_rp + sign(0.5_rp,vel) ) &
876  + ( 2.0_rp * val(ke-3,i,j) + 5.0_rp * val(ke-2,i,j) - val(ke-1,i,j) ) / 6.0_rp &
877  * ( 0.5_rp - sign(0.5_rp,vel) ) ) ! k = KE-4
878 
879  flux(ke-1,i,j) = 0.0_rp
880  enddo
881  enddo
882  !$acc end kernels
883  !$omp end do nowait
884 
885  !$acc end data
886 
887  !$omp end parallel
888 
889  return
890  end subroutine atmos_dyn_fvm_fluxj13_xyw_ud7
891 
892  !-----------------------------------------------------------------------------
894  subroutine atmos_dyn_fvm_fluxj23_xyw_ud7( &
895  flux, &
896  mom, val, DENS, &
897  GSQRT, J23G, MAPF, &
898  CDZ, TwoD, &
899  IIS, IIE, JJS, JJE )
900  implicit none
901 
902  real(rp), intent(inout) :: flux (ka,ia,ja)
903  real(rp), intent(in) :: mom (ka,ia,ja)
904  real(rp), intent(in) :: val (ka,ia,ja)
905  real(rp), intent(in) :: dens (ka,ia,ja)
906  real(rp), intent(in) :: gsqrt (ka,ia,ja)
907  real(rp), intent(in) :: j23g (ka,ia,ja)
908  real(rp), intent(in) :: mapf ( ia,ja,2)
909  real(rp), intent(in) :: cdz (ka)
910  logical, intent(in) :: twod
911  integer, intent(in) :: iis, iie, jjs, jje
912 
913  real(rp) :: vel
914  integer :: k, i, j
915  !---------------------------------------------------------------------------
916 
917  !$omp parallel default(none) private(i,j,k,vel) &
918  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J23G,MAPF)
919 
920  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, J23G, MAPF, CDZ)
921 
922  !$omp do OMP_SCHEDULE_ collapse(2)
923  !$acc kernels
924  do j = jjs, jje
925  do i = iis, iie
926  do k = ks+3, ke-3
927  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j-1) ) ) &
928  / dens(k,i,j)
929  vel = vel * j23g(k,i,j)
930  flux(k-1,i,j) = vel / mapf(i,j,+1) &
931  * ( ( f71 * ( val(k+3,i,j)+val(k-4,i,j) ) &
932  + f72 * ( val(k+2,i,j)+val(k-3,i,j) ) &
933  + f73 * ( val(k+1,i,j)+val(k-2,i,j) ) &
934  + f74 * ( val(k,i,j)+val(k-1,i,j) ) ) &
935  - ( f71 * ( val(k+3,i,j)-val(k-4,i,j) ) &
936  + f75 * ( val(k+2,i,j)-val(k-3,i,j) ) &
937  + f76 * ( val(k+1,i,j)-val(k-2,i,j) ) &
938  + f77 * ( val(k,i,j)-val(k-1,i,j) ) ) * sign(1.0_rp,vel) )
939  enddo
940  enddo
941  enddo
942  !$acc end kernels
943  !$omp end do nowait
944 
945  !$omp do OMP_SCHEDULE_ collapse(2)
946  !$acc kernels
947  do j = jjs, jje
948  do i = iis, iie
949  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS.
950  ! The flux at KS can be non-zero.
951  ! To reduce calculations, all the fluxes are set to zero.
952  flux(ks-1,i,j) = 0.0_rp ! k = KS
953 
954  ! physically incorrect but for numerical stability
955  vel = ( ( 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i,j-1) ) ) / dens(ks+1,i,j) &
956  + ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i,j-1) ) ) / dens(ks ,i,j) ) * 0.5_rp
957 ! vel = ( 0.5_RP * ( mom(KS+1,i,j)+mom(KS+1,i,j-1) ) ) &
958 ! / DENS(KS+1,i,j)
959  vel = vel * j23g(ks+1,i,j)
960  flux(ks,i,j) = vel / mapf(i,j,+1) &
961  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
962  * ( 0.5_rp + sign(0.5_rp,vel) ) &
963  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
964  * ( 0.5_rp - sign(0.5_rp,vel) ) ) ! k = KS+1
965 
966  vel = ( 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i,j-1) ) ) &
967  / dens(ks+2,i,j)
968  vel = vel * j23g(ks,i,j)
969  flux(ks+1,i,j) = vel / mapf(i,j,+1) &
970  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
971  * ( 0.5_rp + sign(0.5_rp,vel) ) &
972  + ( - 3.0_rp * val(ks,i,j) &
973  + 27.0_rp * val(ks+1,i,j) &
974  + 47.0_rp * val(ks+2,i,j) &
975  - 13.0_rp * val(ks+3,i,j) &
976  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
977  * ( 0.5_rp - sign(0.5_rp,vel) ) ) ! k = KS+2
978 
979 
980  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j-1) ) ) &
981  / dens(ke-1,i,j)
982  vel = vel * j23g(ke-1,i,j)
983  flux(ke-2,i,j) = vel / mapf(i,j,+1) &
984  * ( ( 2.0_rp * val(ke-1,i,j) + 5.0_rp * val(ke-2,i,j) - val(ke-3,i,j) ) / 6.0_rp &
985  * ( 0.5_rp + sign(0.5_rp,vel) ) &
986  + f2 * ( val(ke-1,i,j)+val(ke-2,i,j) ) &
987  * ( 0.5_rp - sign(0.5_rp,vel) ) ) ! k = KE-3
988 
989  vel = ( 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i,j-1) ) ) &
990  / dens(ke-2,i,j)
991  vel = vel * j23g(ke-2,i,j)
992  flux(ke-3,i,j) = vel / mapf(i,j,+1) &
993  * ( ( - 3.0_rp * val(ke-1,i,j) &
994  + 27.0_rp * val(ke-2,i,j) &
995  + 47.0_rp * val(ke-3,i,j) &
996  - 13.0_rp * val(ke-4,i,j) &
997  + 2.0_rp * val(ke-5,i,j) ) / 60.0_rp &
998  * ( 0.5_rp + sign(0.5_rp,vel) ) &
999  + ( 2.0_rp * val(ke-3,i,j) + 5.0_rp * val(ke-2,i,j) - val(ke-1,i,j) ) / 6.0_rp &
1000  * ( 0.5_rp - sign(0.5_rp,vel) ) ) ! k = KE-4
1001 
1002  flux(ke-1,i,j) = 0.0_rp
1003  enddo
1004  enddo
1005  !$acc end kernels
1006  !$omp end do nowait
1007 
1008  !$acc end data
1009 
1010  !$omp end parallel
1011 
1012  return
1013  end subroutine atmos_dyn_fvm_fluxj23_xyw_ud7
1014 
1015 
1016  !-----------------------------------------------------------------------------
1018  subroutine atmos_dyn_fvm_fluxx_xyw_ud7( &
1019  flux, &
1020  mom, val, DENS, &
1021  GSQRT, MAPF, &
1022  num_diff, &
1023  CDZ, TwoD, &
1024  IIS, IIE, JJS, JJE )
1025  implicit none
1026 
1027  real(rp), intent(inout) :: flux (ka,ia,ja)
1028  real(rp), intent(in) :: mom (ka,ia,ja)
1029  real(rp), intent(in) :: val (ka,ia,ja)
1030  real(rp), intent(in) :: dens (ka,ia,ja)
1031  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1032  real(rp), intent(in) :: mapf ( ia,ja,2)
1033  real(rp), intent(in) :: num_diff(ka,ia,ja)
1034  real(rp), intent(in) :: cdz (ka)
1035  logical, intent(in) :: twod
1036  integer, intent(in) :: iis, iie, jjs, jje
1037 
1038  real(rp) :: vel
1039  integer :: k, i, j
1040  !---------------------------------------------------------------------------
1041 
1042  !$omp parallel default(none) private(i,j,k,vel) &
1043  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff) &
1044  !$omp shared(CDZ)
1045 
1046  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, MAPF, num_diff, CDZ)
1047 
1048  !$omp do OMP_SCHEDULE_ collapse(2)
1049  !$acc kernels
1050  do j = jjs, jje
1051  do i = iis-1, iie
1052  do k = ks, ke-1
1053 #ifdef DEBUG
1054  call check( __line__, mom(k ,i,j) )
1055  call check( __line__, mom(k+1,i,j) )
1056 
1057  call check( __line__, val(k,i,j) )
1058  call check( __line__, val(k,i+1,j) )
1059 
1060  call check( __line__, val(k,i-1,j) )
1061  call check( __line__, val(k,i+2,j) )
1062 
1063  call check( __line__, val(k,i-2,j) )
1064  call check( __line__, val(k,i+3,j) )
1065 
1066  call check( __line__, val(k,i-3,j) )
1067  call check( __line__, val(k,i+4,j) )
1068 
1069 #endif
1070  vel = ( f2h(k,1,i_uyz) &
1071  * mom(k+1,i,j) &
1072  + f2h(k,2,i_uyz) &
1073  * mom(k,i,j) ) &
1074  / ( f2h(k,1,i_uyz) &
1075  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
1076  + f2h(k,2,i_uyz) &
1077  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
1078  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
1079  * ( ( f71 * ( val(k,i+4,j)+val(k,i-3,j) ) &
1080  + f72 * ( val(k,i+3,j)+val(k,i-2,j) ) &
1081  + f73 * ( val(k,i+2,j)+val(k,i-1,j) ) &
1082  + f74 * ( val(k,i+1,j)+val(k,i,j) ) ) &
1083  - ( f71 * ( val(k,i+4,j)-val(k,i-3,j) ) &
1084  + f75 * ( val(k,i+3,j)-val(k,i-2,j) ) &
1085  + f76 * ( val(k,i+2,j)-val(k,i-1,j) ) &
1086  + f77 * ( val(k,i+1,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
1087  + gsqrt(k,i,j) * num_diff(k,i,j)
1088  enddo
1089  enddo
1090  enddo
1091  !$acc end kernels
1092  !$omp end do nowait
1093 #ifdef DEBUG
1094  k = iundef; i = iundef; j = iundef
1095 #endif
1096 
1097  !$omp do OMP_SCHEDULE_ collapse(2)
1098  !$acc kernels
1099  do j = jjs, jje
1100  do i = iis-1, iie
1101  flux(ke,i,j) = 0.0_rp
1102  enddo
1103  enddo
1104  !$acc end kernels
1105  !$omp end do nowait
1106 
1107  !$acc end data
1108 
1109  !$omp end parallel
1110 #ifdef DEBUG
1111  k = iundef; i = iundef; j = iundef
1112 #endif
1113 
1114  return
1115  end subroutine atmos_dyn_fvm_fluxx_xyw_ud7
1116 
1117  !-----------------------------------------------------------------------------
1119  subroutine atmos_dyn_fvm_fluxy_xyw_ud7( &
1120  flux, &
1121  mom, val, DENS, &
1122  GSQRT, MAPF, &
1123  num_diff, &
1124  CDZ, TwoD, &
1125  IIS, IIE, JJS, JJE )
1126  implicit none
1127 
1128  real(rp), intent(inout) :: flux (ka,ia,ja)
1129  real(rp), intent(in) :: mom (ka,ia,ja)
1130  real(rp), intent(in) :: val (ka,ia,ja)
1131  real(rp), intent(in) :: dens (ka,ia,ja)
1132  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1133  real(rp), intent(in) :: mapf ( ia,ja,2)
1134  real(rp), intent(in) :: num_diff(ka,ia,ja)
1135  real(rp), intent(in) :: cdz (ka)
1136  logical, intent(in) :: twod
1137  integer, intent(in) :: iis, iie, jjs, jje
1138 
1139  real(rp) :: vel
1140  integer :: k, i, j
1141  !---------------------------------------------------------------------------
1142 
1143  !$omp parallel default(none) private(i,j,k,vel) &
1144  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff) &
1145  !$omp shared(CDZ)
1146 
1147  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, MAPF, num_diff, CDZ)
1148 
1149  !$omp do OMP_SCHEDULE_ collapse(2)
1150  !$acc kernels
1151  do j = jjs-1, jje
1152  do i = iis, iie
1153  do k = ks, ke-1
1154 #ifdef DEBUG
1155  call check( __line__, mom(k ,i,j) )
1156  call check( __line__, mom(k+1,i,j) )
1157 
1158  call check( __line__, val(k,i,j) )
1159  call check( __line__, val(k,i,j+1) )
1160 
1161  call check( __line__, val(k,i,j-1) )
1162  call check( __line__, val(k,i,j+2) )
1163 
1164  call check( __line__, val(k,i,j-2) )
1165  call check( __line__, val(k,i,j+3) )
1166 
1167  call check( __line__, val(k,i,j-3) )
1168  call check( __line__, val(k,i,j+4) )
1169 
1170 #endif
1171  vel = ( f2h(k,1,i_xvz) &
1172  * mom(k+1,i,j) &
1173  + f2h(k,2,i_xvz) &
1174  * mom(k,i,j) ) &
1175  / ( f2h(k,1,i_xvz) &
1176  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
1177  + f2h(k,2,i_xvz) &
1178  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
1179  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
1180  * ( ( f71 * ( val(k,i,j+4)+val(k,i,j-3) ) &
1181  + f72 * ( val(k,i,j+3)+val(k,i,j-2) ) &
1182  + f73 * ( val(k,i,j+2)+val(k,i,j-1) ) &
1183  + f74 * ( val(k,i,j+1)+val(k,i,j) ) ) &
1184  - ( f71 * ( val(k,i,j+4)-val(k,i,j-3) ) &
1185  + f75 * ( val(k,i,j+3)-val(k,i,j-2) ) &
1186  + f76 * ( val(k,i,j+2)-val(k,i,j-1) ) &
1187  + f77 * ( val(k,i,j+1)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
1188  + gsqrt(k,i,j) * num_diff(k,i,j)
1189  enddo
1190  enddo
1191  enddo
1192  !$acc end kernels
1193  !$omp end do nowait
1194 #ifdef DEBUG
1195  k = iundef; i = iundef; j = iundef
1196 #endif
1197 
1198  !$omp do OMP_SCHEDULE_ collapse(2)
1199  !$acc kernels
1200  do j = jjs-1, jje
1201  do i = iis, iie
1202  flux(ke,i,j) = 0.0_rp
1203  enddo
1204  enddo
1205  !$acc end kernels
1206  !$omp end do nowait
1207 
1208  !$acc end data
1209 
1210  !$omp end parallel
1211 #ifdef DEBUG
1212  k = iundef; i = iundef; j = iundef
1213 #endif
1214 
1215  return
1216  end subroutine atmos_dyn_fvm_fluxy_xyw_ud7
1217 
1218 
1219  !-----------------------------------------------------------------------------
1221  subroutine atmos_dyn_fvm_fluxz_uyz_ud7( &
1222  flux, &
1223  mom, val, DENS, &
1224  GSQRT, J33G, &
1225  num_diff, &
1226  CDZ, TwoD, &
1227  IIS, IIE, JJS, JJE )
1228  implicit none
1229 
1230  real(rp), intent(inout) :: flux (ka,ia,ja)
1231  real(rp), intent(in) :: mom (ka,ia,ja)
1232  real(rp), intent(in) :: val (ka,ia,ja)
1233  real(rp), intent(in) :: dens (ka,ia,ja)
1234  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1235  real(rp), intent(in) :: j33g
1236  real(rp), intent(in) :: num_diff(ka,ia,ja)
1237  real(rp), intent(in) :: cdz (ka)
1238  logical, intent(in) :: twod
1239  integer, intent(in) :: iis, iie, jjs, jje
1240 
1241  real(rp) :: vel
1242  integer :: k, i, j
1243  !---------------------------------------------------------------------------
1244 
1245  !$omp parallel default(none) private(i,j,k,vel) &
1246  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J33G,GSQRT,num_diff) &
1247  !$omp shared(CDZ,TwoD)
1248 
1249  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, num_diff, CDZ)
1250 
1251 
1252  if ( twod ) then
1253 
1254  !$omp do OMP_SCHEDULE_
1255  !$acc kernels
1256  do j = jjs, jje
1257  do k = ks+3, ke-4
1258  i = iis
1259 #ifdef DEBUG
1260  call check( __line__, mom(k,i,j) )
1261 
1262  call check( __line__, val(k,i,j) )
1263  call check( __line__, val(k+1,i,j) )
1264 
1265  call check( __line__, val(k-1,i,j) )
1266  call check( __line__, val(k+2,i,j) )
1267 
1268  call check( __line__, val(k-2,i,j) )
1269  call check( __line__, val(k+3,i,j) )
1270 
1271  call check( __line__, val(k-3,i,j) )
1272  call check( __line__, val(k+4,i,j) )
1273 
1274 #endif
1275  vel = ( mom(k,i,j) ) &
1276  / ( f2h(k,1,i_xyz) &
1277  * dens(k+1,i,j) &
1278  + f2h(k,2,i_xyz) &
1279  * dens(k,i,j) )
1280  flux(k,i,j) = j33g * vel &
1281  * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
1282  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1283  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1284  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
1285  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
1286  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
1287  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
1288  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
1289  + gsqrt(k,i,j) * num_diff(k,i,j)
1290  enddo
1291  enddo
1292  !$acc end kernels
1293  !$omp end do nowait
1294 #ifdef DEBUG
1295  k = iundef; i = iundef; j = iundef
1296 #endif
1297 
1298  !$omp do OMP_SCHEDULE_
1299  !$acc kernels
1300  do j = jjs, jje
1301  i = iis
1302 #ifdef DEBUG
1303 
1304  call check( __line__, mom(ks,i ,j) )
1305  call check( __line__, val(ks+1,i,j) )
1306  call check( __line__, val(ks,i,j) )
1307 
1308  call check( __line__, mom(ks+1,i ,j) )
1309  call check( __line__, val(ks+3,i,j) )
1310  call check( __line__, val(ks+2,i,j) )
1311 
1312  call check( __line__, mom(ks+2,i ,j) )
1313  call check( __line__, val(ks+5,i,j) )
1314  call check( __line__, val(ks+4,i,j) )
1315 
1316 #endif
1317  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1318  ! The flux at KS-1 can be non-zero.
1319  ! To reduce calculations, all the fluxes are set to zero.
1320  flux(ks-1,i,j) = 0.0_rp
1321 
1322  vel = ( mom(ks,i,j) ) &
1323  / ( f2h(ks,1,i_xyz) &
1324  * dens(ks+1,i,j) &
1325  + f2h(ks,2,i_xyz) &
1326  * dens(ks,i,j) )
1327  flux(ks,i,j) = j33g * vel &
1328  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
1329  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1330  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
1331  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1332  + gsqrt(ks,i,j) * num_diff(ks,i,j)
1333  vel = ( mom(ke-1,i,j) ) &
1334  / ( f2h(ke-1,1,i_xyz) &
1335  * dens(ke,i,j) &
1336  + f2h(ke-1,2,i_xyz) &
1337  * dens(ke-1,i,j) )
1338  flux(ke-1,i,j) = j33g * vel &
1339  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
1340  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1341  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
1342  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1343  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
1344 
1345  vel = ( mom(ks+1,i,j) ) &
1346  / ( f2h(ks+1,1,i_xyz) &
1347  * dens(ks+2,i,j) &
1348  + f2h(ks+1,2,i_xyz) &
1349  * dens(ks+1,i,j) )
1350  flux(ks+1,i,j) = j33g * vel &
1351  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
1352  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1353  + ( - 3.0_rp * val(ks,i,j) &
1354  + 27.0_rp * val(ks+1,i,j) &
1355  + 47.0_rp * val(ks+2,i,j) &
1356  - 13.0_rp * val(ks+3,i,j) &
1357  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
1358  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1359  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j)
1360  vel = ( mom(ke-2,i,j) ) &
1361  / ( f2h(ke-2,1,i_xyz) &
1362  * dens(ke-1,i,j) &
1363  + f2h(ke-2,2,i_xyz) &
1364  * dens(ke-2,i,j) )
1365  flux(ke-2,i,j) = j33g * vel &
1366  * ( ( - 3.0_rp * val(ke,i,j) &
1367  + 27.0_rp * val(ke-1,i,j) &
1368  + 47.0_rp * val(ke-2,i,j) &
1369  - 13.0_rp * val(ke-3,i,j) &
1370  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
1371  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1372  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
1373  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1374  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j)
1375 
1376  vel = ( mom(ks+2,i,j) ) &
1377  / ( f2h(ks+2,1,i_xyz) &
1378  * dens(ks+3,i,j) &
1379  + f2h(ks+2,2,i_xyz) &
1380  * dens(ks+2,i,j) )
1381  flux(ks+2,i,j) = j33g * vel &
1382  * ( ( - 3.0_rp * val(ks+4,i,j) &
1383  + 27.0_rp * val(ks+3,i,j) &
1384  + 47.0_rp * val(ks+2,i,j) &
1385  - 13.0_rp * val(ks+1,i,j) &
1386  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
1387  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1388  + ( 4.0_rp * val(ks,i,j) &
1389  - 38.0_rp * val(ks+1,i,j) &
1390  + 214.0_rp * val(ks+2,i,j) &
1391  + 319.0_rp * val(ks+3,i,j) &
1392  - 101.0_rp * val(ks+4,i,j) &
1393  + 25.0_rp * val(ks+5,i,j) &
1394  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
1395  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1396  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j)
1397  vel = ( mom(ke-3,i,j) ) &
1398  / ( f2h(ke-3,1,i_xyz) &
1399  * dens(ke-2,i,j) &
1400  + f2h(ke-3,2,i_xyz) &
1401  * dens(ke-3,i,j) )
1402  flux(ke-3,i,j) = j33g * vel &
1403  * ( ( 4.0_rp * val(ke,i,j) &
1404  - 38.0_rp * val(ke-1,i,j) &
1405  + 214.0_rp * val(ke-2,i,j) &
1406  + 319.0_rp * val(ke-3,i,j) &
1407  - 101.0_rp * val(ke-4,i,j) &
1408  + 25.0_rp * val(ke-5,i,j) &
1409  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
1410  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1411  + ( - 3.0_rp * val(ke-4,i,j) &
1412  + 27.0_rp * val(ke-3,i,j) &
1413  + 47.0_rp * val(ke-2,i,j) &
1414  - 13.0_rp * val(ke-1,i,j) &
1415  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
1416  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1417  + gsqrt(ke-3,i,j) * num_diff(ke-3,i,j)
1418 
1419  flux(ke,i,j) = 0.0_rp
1420  enddo
1421  !$acc end kernels
1422  !$omp end do nowait
1423 
1424  else
1425 
1426 
1427  !$omp do OMP_SCHEDULE_ collapse(2)
1428  !$acc kernels
1429  do j = jjs, jje
1430  do i = iis, iie
1431  do k = ks+3, ke-4
1432 #ifdef DEBUG
1433  call check( __line__, mom(k,i,j) )
1434  call check( __line__, mom(k,i+1,j) )
1435 
1436  call check( __line__, val(k,i,j) )
1437  call check( __line__, val(k+1,i,j) )
1438 
1439  call check( __line__, val(k-1,i,j) )
1440  call check( __line__, val(k+2,i,j) )
1441 
1442  call check( __line__, val(k-2,i,j) )
1443  call check( __line__, val(k+3,i,j) )
1444 
1445  call check( __line__, val(k-3,i,j) )
1446  call check( __line__, val(k+4,i,j) )
1447 
1448 #endif
1449  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i+1,j) ) ) &
1450  / ( f2h(k,1,i_uyz) &
1451  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
1452  + f2h(k,2,i_uyz) &
1453  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
1454  flux(k,i,j) = j33g * vel &
1455  * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
1456  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1457  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1458  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
1459  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
1460  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
1461  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
1462  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
1463  + gsqrt(k,i,j) * num_diff(k,i,j)
1464  enddo
1465  enddo
1466  enddo
1467  !$acc end kernels
1468  !$omp end do nowait
1469 #ifdef DEBUG
1470  k = iundef; i = iundef; j = iundef
1471 #endif
1472 
1473  !$omp do OMP_SCHEDULE_ collapse(2)
1474  !$acc kernels
1475  do j = jjs, jje
1476  do i = iis, iie
1477 #ifdef DEBUG
1478 
1479  call check( __line__, mom(ks,i ,j) )
1480  call check( __line__, mom(ks,i+1,j) )
1481  call check( __line__, val(ks+1,i,j) )
1482  call check( __line__, val(ks,i,j) )
1483 
1484  call check( __line__, mom(ks+1,i ,j) )
1485  call check( __line__, mom(ks+1,i+1,j) )
1486  call check( __line__, val(ks+3,i,j) )
1487  call check( __line__, val(ks+2,i,j) )
1488 
1489  call check( __line__, mom(ks+2,i ,j) )
1490  call check( __line__, mom(ks+2,i+1,j) )
1491  call check( __line__, val(ks+5,i,j) )
1492  call check( __line__, val(ks+4,i,j) )
1493 
1494 #endif
1495  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1496  ! The flux at KS-1 can be non-zero.
1497  ! To reduce calculations, all the fluxes are set to zero.
1498  flux(ks-1,i,j) = 0.0_rp
1499 
1500  vel = ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i+1,j) ) ) &
1501  / ( f2h(ks,1,i_uyz) &
1502  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
1503  + f2h(ks,2,i_uyz) &
1504  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
1505  flux(ks,i,j) = j33g * vel &
1506  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
1507  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1508  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
1509  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1510  + gsqrt(ks,i,j) * num_diff(ks,i,j)
1511  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i+1,j) ) ) &
1512  / ( f2h(ke-1,1,i_uyz) &
1513  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
1514  + f2h(ke-1,2,i_uyz) &
1515  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
1516  flux(ke-1,i,j) = j33g * vel &
1517  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
1518  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1519  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
1520  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1521  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
1522 
1523  vel = ( 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i+1,j) ) ) &
1524  / ( f2h(ks+1,1,i_uyz) &
1525  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) &
1526  + f2h(ks+1,2,i_uyz) &
1527  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) )
1528  flux(ks+1,i,j) = j33g * vel &
1529  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
1530  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1531  + ( - 3.0_rp * val(ks,i,j) &
1532  + 27.0_rp * val(ks+1,i,j) &
1533  + 47.0_rp * val(ks+2,i,j) &
1534  - 13.0_rp * val(ks+3,i,j) &
1535  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
1536  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1537  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j)
1538  vel = ( 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i+1,j) ) ) &
1539  / ( f2h(ke-2,1,i_uyz) &
1540  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) &
1541  + f2h(ke-2,2,i_uyz) &
1542  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) )
1543  flux(ke-2,i,j) = j33g * vel &
1544  * ( ( - 3.0_rp * val(ke,i,j) &
1545  + 27.0_rp * val(ke-1,i,j) &
1546  + 47.0_rp * val(ke-2,i,j) &
1547  - 13.0_rp * val(ke-3,i,j) &
1548  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
1549  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1550  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
1551  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1552  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j)
1553 
1554  vel = ( 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i+1,j) ) ) &
1555  / ( f2h(ks+2,1,i_uyz) &
1556  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i+1,j) ) &
1557  + f2h(ks+2,2,i_uyz) &
1558  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) )
1559  flux(ks+2,i,j) = j33g * vel &
1560  * ( ( - 3.0_rp * val(ks+4,i,j) &
1561  + 27.0_rp * val(ks+3,i,j) &
1562  + 47.0_rp * val(ks+2,i,j) &
1563  - 13.0_rp * val(ks+1,i,j) &
1564  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
1565  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1566  + ( 4.0_rp * val(ks,i,j) &
1567  - 38.0_rp * val(ks+1,i,j) &
1568  + 214.0_rp * val(ks+2,i,j) &
1569  + 319.0_rp * val(ks+3,i,j) &
1570  - 101.0_rp * val(ks+4,i,j) &
1571  + 25.0_rp * val(ks+5,i,j) &
1572  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
1573  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1574  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j)
1575  vel = ( 0.5_rp * ( mom(ke-3,i,j)+mom(ke-3,i+1,j) ) ) &
1576  / ( f2h(ke-3,1,i_uyz) &
1577  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) &
1578  + f2h(ke-3,2,i_uyz) &
1579  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i+1,j) ) )
1580  flux(ke-3,i,j) = j33g * vel &
1581  * ( ( 4.0_rp * val(ke,i,j) &
1582  - 38.0_rp * val(ke-1,i,j) &
1583  + 214.0_rp * val(ke-2,i,j) &
1584  + 319.0_rp * val(ke-3,i,j) &
1585  - 101.0_rp * val(ke-4,i,j) &
1586  + 25.0_rp * val(ke-5,i,j) &
1587  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
1588  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1589  + ( - 3.0_rp * val(ke-4,i,j) &
1590  + 27.0_rp * val(ke-3,i,j) &
1591  + 47.0_rp * val(ke-2,i,j) &
1592  - 13.0_rp * val(ke-1,i,j) &
1593  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
1594  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1595  + gsqrt(ke-3,i,j) * num_diff(ke-3,i,j)
1596 
1597  flux(ke,i,j) = 0.0_rp
1598  enddo
1599  enddo
1600  !$acc end kernels
1601  !$omp end do nowait
1602 
1603  end if
1604 
1605 
1606  !$acc end data
1607 
1608  !$omp end parallel
1609 #ifdef DEBUG
1610  k = iundef; i = iundef; j = iundef
1611 #endif
1612 
1613  return
1614  end subroutine atmos_dyn_fvm_fluxz_uyz_ud7
1615 
1616  !-----------------------------------------------------------------------------
1618  subroutine atmos_dyn_fvm_fluxj13_uyz_ud7( &
1619  flux, &
1620  mom, val, DENS, &
1621  GSQRT, J13G, MAPF, &
1622  CDZ, TwoD, &
1623  IIS, IIE, JJS, JJE )
1624  implicit none
1625 
1626  real(rp), intent(inout) :: flux (ka,ia,ja)
1627  real(rp), intent(in) :: mom (ka,ia,ja)
1628  real(rp), intent(in) :: val (ka,ia,ja)
1629  real(rp), intent(in) :: dens (ka,ia,ja)
1630  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1631  real(rp), intent(in) :: j13g (ka,ia,ja)
1632  real(rp), intent(in) :: mapf ( ia,ja,2)
1633  real(rp), intent(in) :: cdz (ka)
1634  logical, intent(in) :: twod
1635  integer, intent(in) :: iis, iie, jjs, jje
1636 
1637  real(rp) :: vel
1638  integer :: k, i, j
1639  !---------------------------------------------------------------------------
1640 
1641  !$omp parallel default(none) private(i,j,k,vel) &
1642  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J13G,MAPF) &
1643  !$omp shared(GSQRT,CDZ,TwoD)
1644 
1645  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, J13G, MAPF, CDZ)
1646 
1647 
1648 
1649  !$omp do OMP_SCHEDULE_ collapse(2)
1650  !$acc kernels
1651  do j = jjs, jje
1652  do i = iis, iie
1653  do k = ks+3, ke-4
1654  vel = ( f2h(k,1,i_uyz) &
1655  * mom(k+1,i,j) &
1656  + f2h(k,2,i_uyz) &
1657  * mom(k,i,j) ) &
1658  / ( f2h(k,1,i_uyz) &
1659  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
1660  + f2h(k,2,i_uyz) &
1661  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
1662  vel = vel * j13g(k,i,j)
1663  flux(k,i,j) = vel / mapf(i,j,+2) &
1664  * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
1665  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1666  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1667  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
1668  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
1669  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
1670  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
1671  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) )
1672  enddo
1673  enddo
1674  enddo
1675  !$acc end kernels
1676  !$omp end do nowait
1677 
1678  !$omp do OMP_SCHEDULE_ collapse(2)
1679  !$acc kernels
1680  do j = jjs, jje
1681  do i = iis, iie
1682  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1683  ! The flux at KS-1 can be non-zero.
1684  ! To reduce calculations, all the fluxes are set to zero.
1685  flux(ks-1,i,j) = 0.0_rp
1686 
1687  vel = ( f2h(ks,1,i_uyz) &
1688  * mom(ks+1,i,j) &
1689  + f2h(ks,2,i_uyz) &
1690  * mom(ks,i,j) ) &
1691  / ( f2h(ks,1,i_uyz) &
1692  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
1693  + f2h(ks,2,i_uyz) &
1694  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
1695  vel = vel * j13g(ks,i,j)
1696  flux(ks,i,j) = vel / mapf(i,j,+2) &
1697  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
1698  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1699  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
1700  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1701 
1702  vel = ( f2h(ke-1,1,i_uyz) &
1703  * mom(ke,i,j) &
1704  + f2h(ke-1,2,i_uyz) &
1705  * mom(ke-1,i,j) ) &
1706  / ( f2h(ke-1,1,i_uyz) &
1707  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
1708  + f2h(ke-1,2,i_uyz) &
1709  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
1710  vel = vel * j13g(ke-1,i,j)
1711  flux(ke-1,i,j) = vel / mapf(i,j,+2) &
1712  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
1713  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1714  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
1715  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1716 
1717  vel = ( f2h(ks+1,1,i_uyz) &
1718  * mom(ks+2,i,j) &
1719  + f2h(ks+1,2,i_uyz) &
1720  * mom(ks+1,i,j) ) &
1721  / ( f2h(ks+1,1,i_uyz) &
1722  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) &
1723  + f2h(ks+1,2,i_uyz) &
1724  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) )
1725  vel = vel * j13g(ks+1,i,j)
1726  flux(ks+1,i,j) = vel / mapf(i,j,+2) &
1727  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
1728  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1729  + ( - 3.0_rp * val(ks,i,j) &
1730  + 27.0_rp * val(ks+1,i,j) &
1731  + 47.0_rp * val(ks+2,i,j) &
1732  - 13.0_rp * val(ks+3,i,j) &
1733  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
1734  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1735 
1736  vel = ( f2h(ke-2,1,i_uyz) &
1737  * mom(ke-1,i,j) &
1738  + f2h(ke-2,2,i_uyz) &
1739  * mom(ke-2,i,j) ) &
1740  / ( f2h(ke-2,1,i_uyz) &
1741  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) &
1742  + f2h(ke-2,2,i_uyz) &
1743  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) )
1744  vel = vel * j13g(ke-2,i,j)
1745  flux(ke-2,i,j) = vel / mapf(i,j,+2) &
1746  * ( ( - 3.0_rp * val(ke,i,j) &
1747  + 27.0_rp * val(ke-1,i,j) &
1748  + 47.0_rp * val(ke-2,i,j) &
1749  - 13.0_rp * val(ke-3,i,j) &
1750  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
1751  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1752  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
1753  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1754 
1755  vel = ( f2h(ks+2,1,i_uyz) &
1756  * mom(ks+3,i,j) &
1757  + f2h(ks+2,2,i_uyz) &
1758  * mom(ks+2,i,j) ) &
1759  / ( f2h(ks+2,1,i_uyz) &
1760  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i+1,j) ) &
1761  + f2h(ks+2,2,i_uyz) &
1762  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) )
1763  vel = vel * j13g(ks+2,i,j)
1764  flux(ks+2,i,j) = vel / mapf(i,j,+2) &
1765  * ( ( - 3.0_rp * val(ks+4,i,j) &
1766  + 27.0_rp * val(ks+3,i,j) &
1767  + 47.0_rp * val(ks+2,i,j) &
1768  - 13.0_rp * val(ks+1,i,j) &
1769  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
1770  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1771  + ( 4.0_rp * val(ks,i,j) &
1772  - 38.0_rp * val(ks+1,i,j) &
1773  + 214.0_rp * val(ks+2,i,j) &
1774  + 319.0_rp * val(ks+3,i,j) &
1775  - 101.0_rp * val(ks+4,i,j) &
1776  + 25.0_rp * val(ks+5,i,j) &
1777  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
1778  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1779 
1780  vel = ( f2h(ke-3,1,i_uyz) &
1781  * mom(ke-2,i,j) &
1782  + f2h(ke-3,2,i_uyz) &
1783  * mom(ke-3,i,j) ) &
1784  / ( f2h(ke-3,1,i_uyz) &
1785  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) &
1786  + f2h(ke-3,2,i_uyz) &
1787  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i+1,j) ) )
1788  vel = vel * j13g(ke-3,i,j)
1789  flux(ke-3,i,j) = vel / mapf(i,j,+2) &
1790  * ( ( 4.0_rp * val(ke,i,j) &
1791  - 38.0_rp * val(ke-1,i,j) &
1792  + 214.0_rp * val(ke-2,i,j) &
1793  + 319.0_rp * val(ke-3,i,j) &
1794  - 101.0_rp * val(ke-4,i,j) &
1795  + 25.0_rp * val(ke-5,i,j) &
1796  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
1797  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1798  + ( - 3.0_rp * val(ke-4,i,j) &
1799  + 27.0_rp * val(ke-3,i,j) &
1800  + 47.0_rp * val(ke-2,i,j) &
1801  - 13.0_rp * val(ke-1,i,j) &
1802  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
1803  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1804 
1805  flux(ke ,i,j) = 0.0_rp
1806  enddo
1807  enddo
1808  !$acc end kernels
1809  !$omp end do nowait
1810 
1811 
1812 
1813  !$acc end data
1814 
1815  !$omp end parallel
1816  return
1817  end subroutine atmos_dyn_fvm_fluxj13_uyz_ud7
1818 
1819  !-----------------------------------------------------------------------------
1821  subroutine atmos_dyn_fvm_fluxj23_uyz_ud7( &
1822  flux, &
1823  mom, val, DENS, &
1824  GSQRT, J23G, MAPF, &
1825  CDZ, TwoD, &
1826  IIS, IIE, JJS, JJE )
1827  implicit none
1828 
1829  real(rp), intent(inout) :: flux (ka,ia,ja)
1830  real(rp), intent(in) :: mom (ka,ia,ja)
1831  real(rp), intent(in) :: val (ka,ia,ja)
1832  real(rp), intent(in) :: dens (ka,ia,ja)
1833  real(rp), intent(in) :: gsqrt (ka,ia,ja)
1834  real(rp), intent(in) :: j23g (ka,ia,ja)
1835  real(rp), intent(in) :: mapf ( ia,ja,2)
1836  real(rp), intent(in) :: cdz (ka)
1837  logical, intent(in) :: twod
1838  integer, intent(in) :: iis, iie, jjs, jje
1839 
1840  real(rp) :: vel
1841  integer :: k, i, j
1842  !---------------------------------------------------------------------------
1843 
1844  !$omp parallel default(none) private(i,j,k,vel) &
1845  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J23G,MAPF) &
1846  !$omp shared(GSQRT,CDZ,TwoD)
1847 
1848  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, J23G, MAPF, CDZ)
1849 
1850 
1851  if ( twod ) then
1852 
1853  !$omp do OMP_SCHEDULE_
1854  !$acc kernels
1855  do j = jjs, jje
1856  do k = ks+3, ke-4
1857  i = iis
1858  vel = ( f2h(k,1,i_xyz) &
1859  * 0.5_rp * ( mom(k+1,i,j)+mom(k+1,i,j-1) ) &
1860  + f2h(k,2,i_xyz) &
1861  * 0.5_rp * ( mom(k,i,j)+mom(k,i,j-1) ) ) &
1862  / ( f2h(k,1,i_xyz) &
1863  * dens(k+1,i,j) &
1864  + f2h(k,2,i_xyz) &
1865  * dens(k,i,j) )
1866  vel = vel * j23g(k,i,j)
1867  flux(k,i,j) = vel * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
1868  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
1869  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
1870  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
1871  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
1872  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
1873  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
1874  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) )
1875  enddo
1876  enddo
1877  !$acc end kernels
1878  !$omp end do nowait
1879 
1880  !$omp do OMP_SCHEDULE_
1881  !$acc kernels
1882  do j = jjs, jje
1883  i = iis
1884  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
1885  ! The flux at KS-1 can be non-zero.
1886  ! To reduce calculations, all the fluxes are set to zero.
1887  flux(ks-1,i,j) = 0.0_rp
1888 
1889  vel = ( f2h(ks,1,i_xyz) &
1890  * 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i,j-1) ) &
1891  + f2h(ks,2,i_xyz) &
1892  * 0.5_rp * ( mom(ks,i,j)+mom(ks,i,j-1) ) ) &
1893  / ( f2h(ks,1,i_xyz) &
1894  * dens(ks+1,i,j) &
1895  + f2h(ks,2,i_xyz) &
1896  * dens(ks,i,j) )
1897  vel = vel * j23g(ks,i,j)
1898  flux(ks,i,j) = vel / mapf(i,j,+1) &
1899  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
1900  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1901  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
1902  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1903 
1904  vel = ( f2h(ke-1,1,i_xyz) &
1905  * 0.5_rp * ( mom(ke,i,j)+mom(ke,i,j-1) ) &
1906  + f2h(ke-1,2,i_xyz) &
1907  * 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j-1) ) ) &
1908  / ( f2h(ke-1,1,i_xyz) &
1909  * dens(ke,i,j) &
1910  + f2h(ke-1,2,i_xyz) &
1911  * dens(ke-1,i,j) )
1912  vel = vel * j23g(ke-1,i,j)
1913  flux(ke-1,i,j) = vel / mapf(i,j,+1) &
1914  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
1915  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1916  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
1917  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1918 
1919  vel = ( f2h(ks+1,1,i_xyz) &
1920  * 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i,j-1) ) &
1921  + f2h(ks+1,2,i_xyz) &
1922  * 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i,j-1) ) ) &
1923  / ( f2h(ks+1,1,i_xyz) &
1924  * dens(ks+2,i,j) &
1925  + f2h(ks+1,2,i_xyz) &
1926  * dens(ks+1,i,j) )
1927  vel = vel * j23g(ks+1,i,j)
1928  flux(ks+1,i,j) = vel / mapf(i,j,+1) &
1929  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
1930  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1931  + ( - 3.0_rp * val(ks,i,j) &
1932  + 27.0_rp * val(ks+1,i,j) &
1933  + 47.0_rp * val(ks+2,i,j) &
1934  - 13.0_rp * val(ks+3,i,j) &
1935  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
1936  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1937 
1938  vel = ( f2h(ke-2,1,i_xyz) &
1939  * 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j-1) ) &
1940  + f2h(ke-2,2,i_xyz) &
1941  * 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i,j-1) ) ) &
1942  / ( f2h(ke-2,1,i_xyz) &
1943  * dens(ke-1,i,j) &
1944  + f2h(ke-2,2,i_xyz) &
1945  * dens(ke-2,i,j) )
1946  vel = vel * j23g(ke-2,i,j)
1947  flux(ke-2,i,j) = vel / mapf(i,j,+1) &
1948  * ( ( - 3.0_rp * val(ke,i,j) &
1949  + 27.0_rp * val(ke-1,i,j) &
1950  + 47.0_rp * val(ke-2,i,j) &
1951  - 13.0_rp * val(ke-3,i,j) &
1952  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
1953  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1954  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
1955  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1956 
1957  vel = ( f2h(ks+2,1,i_xyz) &
1958  * 0.5_rp * ( mom(ks+3,i,j)+mom(ks+3,i,j-1) ) &
1959  + f2h(ks+2,2,i_xyz) &
1960  * 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i,j-1) ) ) &
1961  / ( f2h(ks+2,1,i_xyz) &
1962  * dens(ks+3,i,j) &
1963  + f2h(ks+2,2,i_xyz) &
1964  * dens(ks+2,i,j) )
1965  vel = vel * j23g(ks+2,i,j)
1966  flux(ks+2,i,j) = vel / mapf(i,j,+1) &
1967  * ( ( - 3.0_rp * val(ks+4,i,j) &
1968  + 27.0_rp * val(ks+3,i,j) &
1969  + 47.0_rp * val(ks+2,i,j) &
1970  - 13.0_rp * val(ks+1,i,j) &
1971  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
1972  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1973  + ( 4.0_rp * val(ks,i,j) &
1974  - 38.0_rp * val(ks+1,i,j) &
1975  + 214.0_rp * val(ks+2,i,j) &
1976  + 319.0_rp * val(ks+3,i,j) &
1977  - 101.0_rp * val(ks+4,i,j) &
1978  + 25.0_rp * val(ks+5,i,j) &
1979  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
1980  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1981 
1982  vel = ( f2h(ke-3,1,i_xyz) &
1983  * 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i,j-1) ) &
1984  + f2h(ke-3,2,i_xyz) &
1985  * 0.5_rp * ( mom(ke-3,i,j)+mom(ke-3,i,j-1) ) ) &
1986  / ( f2h(ke-3,1,i_xyz) &
1987  * dens(ke-2,i,j) &
1988  + f2h(ke-3,2,i_xyz) &
1989  * dens(ke-3,i,j) )
1990  vel = vel * j23g(ke-3,i,j)
1991  flux(ke-3,i,j) = vel / mapf(i,j,+1) &
1992  * ( ( 4.0_rp * val(ke,i,j) &
1993  - 38.0_rp * val(ke-1,i,j) &
1994  + 214.0_rp * val(ke-2,i,j) &
1995  + 319.0_rp * val(ke-3,i,j) &
1996  - 101.0_rp * val(ke-4,i,j) &
1997  + 25.0_rp * val(ke-5,i,j) &
1998  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
1999  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2000  + ( - 3.0_rp * val(ke-4,i,j) &
2001  + 27.0_rp * val(ke-3,i,j) &
2002  + 47.0_rp * val(ke-2,i,j) &
2003  - 13.0_rp * val(ke-1,i,j) &
2004  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
2005  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2006 
2007  flux(ke ,i,j) = 0.0_rp
2008  enddo
2009  !$acc end kernels
2010  !$omp end do nowait
2011 
2012  else
2013 
2014 
2015  !$omp do OMP_SCHEDULE_ collapse(2)
2016  !$acc kernels
2017  do j = jjs, jje
2018  do i = iis, iie
2019  do k = ks+3, ke-4
2020  vel = ( f2h(k,1,i_uyz) &
2021  * 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) ) &
2022  + f2h(k,2,i_uyz) &
2023  * 0.25_rp * ( mom(k,i,j)+mom(k,i+1,j)+mom(k,i,j-1)+mom(k,i+1,j-1) ) ) &
2024  / ( f2h(k,1,i_uyz) &
2025  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
2026  + f2h(k,2,i_uyz) &
2027  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
2028  vel = vel * j23g(k,i,j)
2029  flux(k,i,j) = vel / mapf(i,j,+1) &
2030  * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
2031  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
2032  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
2033  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
2034  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
2035  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
2036  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
2037  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) )
2038  enddo
2039  enddo
2040  enddo
2041  !$acc end kernels
2042  !$omp end do nowait
2043 
2044  !$omp do OMP_SCHEDULE_ collapse(2)
2045  !$acc kernels
2046  do j = jjs, jje
2047  do i = iis, iie
2048  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
2049  ! The flux at KS-1 can be non-zero.
2050  ! To reduce calculations, all the fluxes are set to zero.
2051  flux(ks-1,i,j) = 0.0_rp
2052 
2053  vel = ( f2h(ks,1,i_uyz) &
2054  * 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) ) &
2055  + f2h(ks,2,i_uyz) &
2056  * 0.25_rp * ( mom(ks,i,j)+mom(ks,i+1,j)+mom(ks,i,j-1)+mom(ks,i+1,j-1) ) ) &
2057  / ( f2h(ks,1,i_uyz) &
2058  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
2059  + f2h(ks,2,i_uyz) &
2060  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
2061  vel = vel * j23g(ks,i,j)
2062  flux(ks,i,j) = vel / mapf(i,j,+1) &
2063  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
2064  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2065  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
2066  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2067 
2068  vel = ( f2h(ke-1,1,i_uyz) &
2069  * 0.25_rp * ( mom(ke,i,j)+mom(ke,i+1,j)+mom(ke,i,j-1)+mom(ke,i+1,j-1) ) &
2070  + f2h(ke-1,2,i_uyz) &
2071  * 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) ) ) &
2072  / ( f2h(ke-1,1,i_uyz) &
2073  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
2074  + f2h(ke-1,2,i_uyz) &
2075  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
2076  vel = vel * j23g(ke-1,i,j)
2077  flux(ke-1,i,j) = vel / mapf(i,j,+1) &
2078  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
2079  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2080  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
2081  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2082 
2083  vel = ( f2h(ks+1,1,i_uyz) &
2084  * 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) ) &
2085  + f2h(ks+1,2,i_uyz) &
2086  * 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) ) ) &
2087  / ( f2h(ks+1,1,i_uyz) &
2088  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) &
2089  + f2h(ks+1,2,i_uyz) &
2090  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) )
2091  vel = vel * j23g(ks+1,i,j)
2092  flux(ks+1,i,j) = vel / mapf(i,j,+1) &
2093  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
2094  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2095  + ( - 3.0_rp * val(ks,i,j) &
2096  + 27.0_rp * val(ks+1,i,j) &
2097  + 47.0_rp * val(ks+2,i,j) &
2098  - 13.0_rp * val(ks+3,i,j) &
2099  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
2100  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2101 
2102  vel = ( f2h(ke-2,1,i_uyz) &
2103  * 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) ) &
2104  + f2h(ke-2,2,i_uyz) &
2105  * 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) ) ) &
2106  / ( f2h(ke-2,1,i_uyz) &
2107  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) &
2108  + f2h(ke-2,2,i_uyz) &
2109  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) )
2110  vel = vel * j23g(ke-2,i,j)
2111  flux(ke-2,i,j) = vel / mapf(i,j,+1) &
2112  * ( ( - 3.0_rp * val(ke,i,j) &
2113  + 27.0_rp * val(ke-1,i,j) &
2114  + 47.0_rp * val(ke-2,i,j) &
2115  - 13.0_rp * val(ke-3,i,j) &
2116  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
2117  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2118  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
2119  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2120 
2121  vel = ( f2h(ks+2,1,i_uyz) &
2122  * 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) ) &
2123  + f2h(ks+2,2,i_uyz) &
2124  * 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) ) ) &
2125  / ( f2h(ks+2,1,i_uyz) &
2126  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i+1,j) ) &
2127  + f2h(ks+2,2,i_uyz) &
2128  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i+1,j) ) )
2129  vel = vel * j23g(ks+2,i,j)
2130  flux(ks+2,i,j) = vel / mapf(i,j,+1) &
2131  * ( ( - 3.0_rp * val(ks+4,i,j) &
2132  + 27.0_rp * val(ks+3,i,j) &
2133  + 47.0_rp * val(ks+2,i,j) &
2134  - 13.0_rp * val(ks+1,i,j) &
2135  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
2136  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2137  + ( 4.0_rp * val(ks,i,j) &
2138  - 38.0_rp * val(ks+1,i,j) &
2139  + 214.0_rp * val(ks+2,i,j) &
2140  + 319.0_rp * val(ks+3,i,j) &
2141  - 101.0_rp * val(ks+4,i,j) &
2142  + 25.0_rp * val(ks+5,i,j) &
2143  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
2144  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2145 
2146  vel = ( f2h(ke-3,1,i_uyz) &
2147  * 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) ) &
2148  + f2h(ke-3,2,i_uyz) &
2149  * 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) ) ) &
2150  / ( f2h(ke-3,1,i_uyz) &
2151  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i+1,j) ) &
2152  + f2h(ke-3,2,i_uyz) &
2153  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i+1,j) ) )
2154  vel = vel * j23g(ke-3,i,j)
2155  flux(ke-3,i,j) = vel / mapf(i,j,+1) &
2156  * ( ( 4.0_rp * val(ke,i,j) &
2157  - 38.0_rp * val(ke-1,i,j) &
2158  + 214.0_rp * val(ke-2,i,j) &
2159  + 319.0_rp * val(ke-3,i,j) &
2160  - 101.0_rp * val(ke-4,i,j) &
2161  + 25.0_rp * val(ke-5,i,j) &
2162  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
2163  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2164  + ( - 3.0_rp * val(ke-4,i,j) &
2165  + 27.0_rp * val(ke-3,i,j) &
2166  + 47.0_rp * val(ke-2,i,j) &
2167  - 13.0_rp * val(ke-1,i,j) &
2168  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
2169  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2170 
2171  flux(ke ,i,j) = 0.0_rp
2172  enddo
2173  enddo
2174  !$acc end kernels
2175  !$omp end do nowait
2176 
2177 
2178  end if
2179 
2180 
2181  !$acc end data
2182 
2183  !$omp end parallel
2184  return
2185  end subroutine atmos_dyn_fvm_fluxj23_uyz_ud7
2186 
2187  !-----------------------------------------------------------------------------
2189  subroutine atmos_dyn_fvm_fluxx_uyz_ud7( &
2190  flux, &
2191  mom, val, DENS, &
2192  GSQRT, MAPF, &
2193  num_diff, &
2194  CDZ, TwoD, &
2195  IIS, IIE, JJS, JJE )
2196  implicit none
2197 
2198  real(rp), intent(inout) :: flux (ka,ia,ja)
2199  real(rp), intent(in) :: mom (ka,ia,ja)
2200  real(rp), intent(in) :: val (ka,ia,ja)
2201  real(rp), intent(in) :: dens (ka,ia,ja)
2202  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2203  real(rp), intent(in) :: mapf ( ia,ja,2)
2204  real(rp), intent(in) :: num_diff(ka,ia,ja)
2205  real(rp), intent(in) :: cdz (ka)
2206  logical, intent(in) :: twod
2207  integer, intent(in) :: iis, iie, jjs, jje
2208 
2209  real(rp) :: vel
2210  integer :: k, i, j
2211  !---------------------------------------------------------------------------
2212 
2213  ! note that x-index is added by -1
2214 
2215 
2216 
2217  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
2218  !$omp private(vel) &
2219  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff)
2220  !$acc kernels
2221  do j = jjs, jje
2222  do i = iis, iie+1
2223  do k = ks, ke
2224 #ifdef DEBUG
2225  call check( __line__, mom(k,i ,j) )
2226  call check( __line__, mom(k,i-1,j) )
2227 
2228  call check( __line__, val(k,i-1,j) )
2229  call check( __line__, val(k,i,j) )
2230 
2231  call check( __line__, val(k,i-2,j) )
2232  call check( __line__, val(k,i+1,j) )
2233 
2234  call check( __line__, val(k,i-3,j) )
2235  call check( __line__, val(k,i+2,j) )
2236 
2237  call check( __line__, val(k,i-4,j) )
2238  call check( __line__, val(k,i+3,j) )
2239 
2240 #endif
2241  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i-1,j) ) ) &
2242  / ( dens(k,i,j) )
2243  flux(k,i-1,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
2244  * ( ( f71 * ( val(k,i+3,j)+val(k,i-4,j) ) &
2245  + f72 * ( val(k,i+2,j)+val(k,i-3,j) ) &
2246  + f73 * ( val(k,i+1,j)+val(k,i-2,j) ) &
2247  + f74 * ( val(k,i,j)+val(k,i-1,j) ) ) &
2248  - ( f71 * ( val(k,i+3,j)-val(k,i-4,j) ) &
2249  + f75 * ( val(k,i+2,j)-val(k,i-3,j) ) &
2250  + f76 * ( val(k,i+1,j)-val(k,i-2,j) ) &
2251  + f77 * ( val(k,i,j)-val(k,i-1,j) ) ) * sign(1.0_rp,vel) ) &
2252  + gsqrt(k,i,j) * num_diff(k,i,j)
2253  enddo
2254  enddo
2255  enddo
2256  !$acc end kernels
2257 #ifdef DEBUG
2258  k = iundef; i = iundef; j = iundef
2259 #endif
2260 
2261 
2262 
2263  return
2264  end subroutine atmos_dyn_fvm_fluxx_uyz_ud7
2265 
2266  !-----------------------------------------------------------------------------
2268  subroutine atmos_dyn_fvm_fluxy_uyz_ud7( &
2269  flux, &
2270  mom, val, DENS, &
2271  GSQRT, MAPF, &
2272  num_diff, &
2273  CDZ, TwoD, &
2274  IIS, IIE, JJS, JJE )
2275  implicit none
2276 
2277  real(rp), intent(inout) :: flux (ka,ia,ja)
2278  real(rp), intent(in) :: mom (ka,ia,ja)
2279  real(rp), intent(in) :: val (ka,ia,ja)
2280  real(rp), intent(in) :: dens (ka,ia,ja)
2281  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2282  real(rp), intent(in) :: mapf ( ia,ja,2)
2283  real(rp), intent(in) :: num_diff(ka,ia,ja)
2284  real(rp), intent(in) :: cdz (ka)
2285  logical, intent(in) :: twod
2286  integer, intent(in) :: iis, iie, jjs, jje
2287 
2288  real(rp) :: vel
2289  integer :: k, i, j
2290  !---------------------------------------------------------------------------
2291 
2292 
2293 
2294  if ( twod ) then
2295 
2296  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ &
2297  !$omp private(vel) &
2298  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff,TwoD)
2299  !$acc kernels
2300  do j = jjs-1, jje
2301  do k = ks, ke
2302  i = iis
2303 #ifdef DEBUG
2304  call check( __line__, mom(k,i ,j) )
2305 
2306  call check( __line__, val(k,i,j) )
2307  call check( __line__, val(k,i,j+1) )
2308 
2309  call check( __line__, val(k,i,j-1) )
2310  call check( __line__, val(k,i,j+2) )
2311 
2312  call check( __line__, val(k,i,j-2) )
2313  call check( __line__, val(k,i,j+3) )
2314 
2315  call check( __line__, val(k,i,j-3) )
2316  call check( __line__, val(k,i,j+4) )
2317 
2318 #endif
2319  vel = ( mom(k,i,j) ) &
2320  / ( 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
2321  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
2322  * ( ( f71 * ( val(k,i,j+4)+val(k,i,j-3) ) &
2323  + f72 * ( val(k,i,j+3)+val(k,i,j-2) ) &
2324  + f73 * ( val(k,i,j+2)+val(k,i,j-1) ) &
2325  + f74 * ( val(k,i,j+1)+val(k,i,j) ) ) &
2326  - ( f71 * ( val(k,i,j+4)-val(k,i,j-3) ) &
2327  + f75 * ( val(k,i,j+3)-val(k,i,j-2) ) &
2328  + f76 * ( val(k,i,j+2)-val(k,i,j-1) ) &
2329  + f77 * ( val(k,i,j+1)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
2330  + gsqrt(k,i,j) * num_diff(k,i,j)
2331  enddo
2332  enddo
2333  !$acc end kernels
2334 #ifdef DEBUG
2335  k = iundef; i = iundef; j = iundef
2336 #endif
2337 
2338  else
2339 
2340 
2341  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
2342  !$omp private(vel) &
2343  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff)
2344  !$acc kernels
2345  do j = jjs-1, jje
2346  do i = iis, iie
2347  do k = ks, ke
2348 #ifdef DEBUG
2349  call check( __line__, mom(k,i ,j) )
2350  call check( __line__, mom(k,i-1,j) )
2351 
2352  call check( __line__, val(k,i,j) )
2353  call check( __line__, val(k,i,j+1) )
2354 
2355  call check( __line__, val(k,i,j-1) )
2356  call check( __line__, val(k,i,j+2) )
2357 
2358  call check( __line__, val(k,i,j-2) )
2359  call check( __line__, val(k,i,j+3) )
2360 
2361  call check( __line__, val(k,i,j-3) )
2362  call check( __line__, val(k,i,j+4) )
2363 
2364 #endif
2365  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i+1,j) ) ) &
2366  / ( 0.25_rp * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) )
2367  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
2368  * ( ( f71 * ( val(k,i,j+4)+val(k,i,j-3) ) &
2369  + f72 * ( val(k,i,j+3)+val(k,i,j-2) ) &
2370  + f73 * ( val(k,i,j+2)+val(k,i,j-1) ) &
2371  + f74 * ( val(k,i,j+1)+val(k,i,j) ) ) &
2372  - ( f71 * ( val(k,i,j+4)-val(k,i,j-3) ) &
2373  + f75 * ( val(k,i,j+3)-val(k,i,j-2) ) &
2374  + f76 * ( val(k,i,j+2)-val(k,i,j-1) ) &
2375  + f77 * ( val(k,i,j+1)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
2376  + gsqrt(k,i,j) * num_diff(k,i,j)
2377  enddo
2378  enddo
2379  enddo
2380  !$acc end kernels
2381 #ifdef DEBUG
2382  k = iundef; i = iundef; j = iundef
2383 #endif
2384 
2385 
2386  end if
2387 
2388 
2389  return
2390  end subroutine atmos_dyn_fvm_fluxy_uyz_ud7
2391 
2392 
2393 
2394  !-----------------------------------------------------------------------------
2396  subroutine atmos_dyn_fvm_fluxz_xvz_ud7( &
2397  flux, &
2398  mom, val, DENS, &
2399  GSQRT, J33G, &
2400  num_diff, &
2401  CDZ, TwoD, &
2402  IIS, IIE, JJS, JJE )
2403  implicit none
2404 
2405  real(rp), intent(inout) :: flux (ka,ia,ja)
2406  real(rp), intent(in) :: mom (ka,ia,ja)
2407  real(rp), intent(in) :: val (ka,ia,ja)
2408  real(rp), intent(in) :: dens (ka,ia,ja)
2409  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2410  real(rp), intent(in) :: j33g
2411  real(rp), intent(in) :: num_diff(ka,ia,ja)
2412  real(rp), intent(in) :: cdz (ka)
2413  logical, intent(in) :: twod
2414  integer, intent(in) :: iis, iie, jjs, jje
2415 
2416  real(rp) :: vel
2417  integer :: k, i, j
2418  !---------------------------------------------------------------------------
2419 
2420  !$omp parallel default(none) private(i,j,k,vel) &
2421  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J33G,GSQRT,num_diff) &
2422  !$omp shared(CDZ,TwoD)
2423 
2424  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, num_diff, CDZ)
2425 
2426 
2427  !$omp do OMP_SCHEDULE_ collapse(2)
2428  !$acc kernels
2429  do j = jjs, jje
2430  do i = iis, iie
2431  do k = ks+3, ke-4
2432 #ifdef DEBUG
2433  call check( __line__, mom(k,i,j) )
2434  call check( __line__, mom(k,i,j+1) )
2435 
2436  call check( __line__, val(k,i,j) )
2437  call check( __line__, val(k+1,i,j) )
2438 
2439  call check( __line__, val(k-1,i,j) )
2440  call check( __line__, val(k+2,i,j) )
2441 
2442  call check( __line__, val(k-2,i,j) )
2443  call check( __line__, val(k+3,i,j) )
2444 
2445  call check( __line__, val(k-3,i,j) )
2446  call check( __line__, val(k+4,i,j) )
2447 
2448 #endif
2449  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j+1) ) ) &
2450  / ( f2h(k,1,i_xvz) &
2451  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
2452  + f2h(k,2,i_xvz) &
2453  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
2454  flux(k,i,j) = j33g * vel &
2455  * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
2456  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
2457  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
2458  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
2459  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
2460  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
2461  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
2462  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
2463  + gsqrt(k,i,j) * num_diff(k,i,j)
2464  enddo
2465  enddo
2466  enddo
2467  !$acc end kernels
2468  !$omp end do nowait
2469 #ifdef DEBUG
2470  k = iundef; i = iundef; j = iundef
2471 #endif
2472 
2473  !$omp do OMP_SCHEDULE_ collapse(2)
2474  !$acc kernels
2475  do j = jjs, jje
2476  do i = iis, iie
2477 #ifdef DEBUG
2478 
2479  call check( __line__, mom(ks,i ,j) )
2480  call check( __line__, mom(ks,i,j+1) )
2481  call check( __line__, val(ks+1,i,j) )
2482  call check( __line__, val(ks,i,j) )
2483 
2484  call check( __line__, mom(ks+1,i ,j) )
2485  call check( __line__, mom(ks+1,i,j+1) )
2486  call check( __line__, val(ks+3,i,j) )
2487  call check( __line__, val(ks+2,i,j) )
2488 
2489  call check( __line__, mom(ks+2,i ,j) )
2490  call check( __line__, mom(ks+2,i,j+1) )
2491  call check( __line__, val(ks+5,i,j) )
2492  call check( __line__, val(ks+4,i,j) )
2493 
2494 #endif
2495  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
2496  ! The flux at KS-1 can be non-zero.
2497  ! To reduce calculations, all the fluxes are set to zero.
2498  flux(ks-1,i,j) = 0.0_rp
2499 
2500  vel = ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i,j+1) ) ) &
2501  / ( f2h(ks,1,i_xvz) &
2502  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
2503  + f2h(ks,2,i_xvz) &
2504  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
2505  flux(ks,i,j) = j33g * vel &
2506  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
2507  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2508  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
2509  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
2510  + gsqrt(ks,i,j) * num_diff(ks,i,j)
2511  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j+1) ) ) &
2512  / ( f2h(ke-1,1,i_xvz) &
2513  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
2514  + f2h(ke-1,2,i_xvz) &
2515  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
2516  flux(ke-1,i,j) = j33g * vel &
2517  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
2518  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2519  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
2520  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
2521  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
2522 
2523  vel = ( 0.5_rp * ( mom(ks+1,i,j)+mom(ks+1,i,j+1) ) ) &
2524  / ( f2h(ks+1,1,i_xvz) &
2525  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) &
2526  + f2h(ks+1,2,i_xvz) &
2527  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) )
2528  flux(ks+1,i,j) = j33g * vel &
2529  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
2530  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2531  + ( - 3.0_rp * val(ks,i,j) &
2532  + 27.0_rp * val(ks+1,i,j) &
2533  + 47.0_rp * val(ks+2,i,j) &
2534  - 13.0_rp * val(ks+3,i,j) &
2535  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
2536  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
2537  + gsqrt(ks+1,i,j) * num_diff(ks+1,i,j)
2538  vel = ( 0.5_rp * ( mom(ke-2,i,j)+mom(ke-2,i,j+1) ) ) &
2539  / ( f2h(ke-2,1,i_xvz) &
2540  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) &
2541  + f2h(ke-2,2,i_xvz) &
2542  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) )
2543  flux(ke-2,i,j) = j33g * vel &
2544  * ( ( - 3.0_rp * val(ke,i,j) &
2545  + 27.0_rp * val(ke-1,i,j) &
2546  + 47.0_rp * val(ke-2,i,j) &
2547  - 13.0_rp * val(ke-3,i,j) &
2548  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
2549  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2550  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
2551  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
2552  + gsqrt(ke-2,i,j) * num_diff(ke-2,i,j)
2553 
2554  vel = ( 0.5_rp * ( mom(ks+2,i,j)+mom(ks+2,i,j+1) ) ) &
2555  / ( f2h(ks+2,1,i_xvz) &
2556  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i,j+1) ) &
2557  + f2h(ks+2,2,i_xvz) &
2558  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) )
2559  flux(ks+2,i,j) = j33g * vel &
2560  * ( ( - 3.0_rp * val(ks+4,i,j) &
2561  + 27.0_rp * val(ks+3,i,j) &
2562  + 47.0_rp * val(ks+2,i,j) &
2563  - 13.0_rp * val(ks+1,i,j) &
2564  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
2565  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2566  + ( 4.0_rp * val(ks,i,j) &
2567  - 38.0_rp * val(ks+1,i,j) &
2568  + 214.0_rp * val(ks+2,i,j) &
2569  + 319.0_rp * val(ks+3,i,j) &
2570  - 101.0_rp * val(ks+4,i,j) &
2571  + 25.0_rp * val(ks+5,i,j) &
2572  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
2573  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
2574  + gsqrt(ks+2,i,j) * num_diff(ks+2,i,j)
2575  vel = ( 0.5_rp * ( mom(ke-3,i,j)+mom(ke-3,i,j+1) ) ) &
2576  / ( f2h(ke-3,1,i_xvz) &
2577  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) &
2578  + f2h(ke-3,2,i_xvz) &
2579  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i,j+1) ) )
2580  flux(ke-3,i,j) = j33g * vel &
2581  * ( ( 4.0_rp * val(ke,i,j) &
2582  - 38.0_rp * val(ke-1,i,j) &
2583  + 214.0_rp * val(ke-2,i,j) &
2584  + 319.0_rp * val(ke-3,i,j) &
2585  - 101.0_rp * val(ke-4,i,j) &
2586  + 25.0_rp * val(ke-5,i,j) &
2587  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
2588  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2589  + ( - 3.0_rp * val(ke-4,i,j) &
2590  + 27.0_rp * val(ke-3,i,j) &
2591  + 47.0_rp * val(ke-2,i,j) &
2592  - 13.0_rp * val(ke-1,i,j) &
2593  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
2594  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
2595  + gsqrt(ke-3,i,j) * num_diff(ke-3,i,j)
2596 
2597  flux(ke,i,j) = 0.0_rp
2598  enddo
2599  enddo
2600  !$acc end kernels
2601  !$omp end do nowait
2602 
2603 
2604  !$acc end data
2605 
2606  !$omp end parallel
2607 #ifdef DEBUG
2608  k = iundef; i = iundef; j = iundef
2609 #endif
2610 
2611  return
2612  end subroutine atmos_dyn_fvm_fluxz_xvz_ud7
2613 
2614  !-----------------------------------------------------------------------------
2616  subroutine atmos_dyn_fvm_fluxj13_xvz_ud7( &
2617  flux, &
2618  mom, val, DENS, &
2619  GSQRT, J13G, MAPF, &
2620  CDZ, TwoD, &
2621  IIS, IIE, JJS, JJE )
2622  implicit none
2623 
2624  real(rp), intent(inout) :: flux (ka,ia,ja)
2625  real(rp), intent(in) :: mom (ka,ia,ja)
2626  real(rp), intent(in) :: val (ka,ia,ja)
2627  real(rp), intent(in) :: dens (ka,ia,ja)
2628  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2629  real(rp), intent(in) :: j13g (ka,ia,ja)
2630  real(rp), intent(in) :: mapf ( ia,ja,2)
2631  real(rp), intent(in) :: cdz (ka)
2632  logical, intent(in) :: twod
2633  integer, intent(in) :: iis, iie, jjs, jje
2634 
2635  real(rp) :: vel
2636  integer :: k, i, j
2637  !---------------------------------------------------------------------------
2638 
2639  !$omp parallel default(none) private(i,j,k,vel) &
2640  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J13G,MAPF) &
2641  !$omp shared(GSQRT,CDZ,TwoD)
2642 
2643  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, J13G, MAPF, CDZ)
2644 
2645 
2646 
2647  !$omp do OMP_SCHEDULE_ collapse(2)
2648  !$acc kernels
2649  do j = jjs, jje
2650  do i = iis, iie
2651  do k = ks+3, ke-4
2652  vel = ( f2h(k,1,i_xvz) &
2653  * 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) ) &
2654  + f2h(k,2,i_xvz) &
2655  * 0.25_rp * ( mom(k,i,j)+mom(k,i-1,j)+mom(k,i,j+1)+mom(k,i-1,j+1) ) ) &
2656  / ( f2h(k,1,i_xvz) &
2657  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
2658  + f2h(k,2,i_xvz) &
2659  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
2660  vel = vel * j13g(k,i,j)
2661  flux(k,i,j) = vel / mapf(i,j,+2) &
2662  * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
2663  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
2664  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
2665  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
2666  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
2667  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
2668  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
2669  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) )
2670  enddo
2671  enddo
2672  enddo
2673  !$acc end kernels
2674  !$omp end do nowait
2675 
2676  !$omp do OMP_SCHEDULE_ collapse(2)
2677  !$acc kernels
2678  do j = jjs, jje
2679  do i = iis, iie
2680  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
2681  ! The flux at KS-1 can be non-zero.
2682  ! To reduce calculations, all the fluxes are set to zero.
2683  flux(ks-1,i,j) = 0.0_rp
2684 
2685  vel = ( f2h(ks,1,i_xvz) &
2686  * 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) ) &
2687  + f2h(ks,2,i_xvz) &
2688  * 0.25_rp * ( mom(ks,i,j)+mom(ks,i-1,j)+mom(ks,i,j+1)+mom(ks,i-1,j+1) ) ) &
2689  / ( f2h(ks,1,i_xvz) &
2690  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
2691  + f2h(ks,2,i_xvz) &
2692  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
2693  vel = vel * j13g(ks,i,j)
2694  flux(ks,i,j) = vel / mapf(i,j,+2) &
2695  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
2696  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2697  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
2698  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2699 
2700  vel = ( f2h(ke-1,1,i_xvz) &
2701  * 0.25_rp * ( mom(ke,i,j)+mom(ke,i-1,j)+mom(ke,i,j+1)+mom(ke,i-1,j+1) ) &
2702  + f2h(ke-1,2,i_xvz) &
2703  * 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) ) ) &
2704  / ( f2h(ke-1,1,i_xvz) &
2705  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
2706  + f2h(ke-1,2,i_xvz) &
2707  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
2708  vel = vel * j13g(ke-1,i,j)
2709  flux(ke-1,i,j) = vel / mapf(i,j,+2) &
2710  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
2711  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2712  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
2713  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2714 
2715  vel = ( f2h(ks+1,1,i_xvz) &
2716  * 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) ) &
2717  + f2h(ks+1,2,i_xvz) &
2718  * 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) ) ) &
2719  / ( f2h(ks+1,1,i_xvz) &
2720  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) &
2721  + f2h(ks+1,2,i_xvz) &
2722  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) )
2723  vel = vel * j13g(ks+1,i,j)
2724  flux(ks+1,i,j) = vel / mapf(i,j,+2) &
2725  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
2726  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2727  + ( - 3.0_rp * val(ks,i,j) &
2728  + 27.0_rp * val(ks+1,i,j) &
2729  + 47.0_rp * val(ks+2,i,j) &
2730  - 13.0_rp * val(ks+3,i,j) &
2731  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
2732  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2733 
2734  vel = ( f2h(ke-2,1,i_xvz) &
2735  * 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) ) &
2736  + f2h(ke-2,2,i_xvz) &
2737  * 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) ) ) &
2738  / ( f2h(ke-2,1,i_xvz) &
2739  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) &
2740  + f2h(ke-2,2,i_xvz) &
2741  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) )
2742  vel = vel * j13g(ke-2,i,j)
2743  flux(ke-2,i,j) = vel / mapf(i,j,+2) &
2744  * ( ( - 3.0_rp * val(ke,i,j) &
2745  + 27.0_rp * val(ke-1,i,j) &
2746  + 47.0_rp * val(ke-2,i,j) &
2747  - 13.0_rp * val(ke-3,i,j) &
2748  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
2749  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2750  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
2751  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2752 
2753  vel = ( f2h(ks+2,1,i_xvz) &
2754  * 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) ) &
2755  + f2h(ks+2,2,i_xvz) &
2756  * 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) ) ) &
2757  / ( f2h(ks+2,1,i_xvz) &
2758  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i,j+1) ) &
2759  + f2h(ks+2,2,i_xvz) &
2760  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) )
2761  vel = vel * j13g(ks+2,i,j)
2762  flux(ks+2,i,j) = vel / mapf(i,j,+2) &
2763  * ( ( - 3.0_rp * val(ks+4,i,j) &
2764  + 27.0_rp * val(ks+3,i,j) &
2765  + 47.0_rp * val(ks+2,i,j) &
2766  - 13.0_rp * val(ks+1,i,j) &
2767  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
2768  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2769  + ( 4.0_rp * val(ks,i,j) &
2770  - 38.0_rp * val(ks+1,i,j) &
2771  + 214.0_rp * val(ks+2,i,j) &
2772  + 319.0_rp * val(ks+3,i,j) &
2773  - 101.0_rp * val(ks+4,i,j) &
2774  + 25.0_rp * val(ks+5,i,j) &
2775  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
2776  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2777 
2778  vel = ( f2h(ke-3,1,i_xvz) &
2779  * 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) ) &
2780  + f2h(ke-3,2,i_xvz) &
2781  * 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) ) ) &
2782  / ( f2h(ke-3,1,i_xvz) &
2783  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) &
2784  + f2h(ke-3,2,i_xvz) &
2785  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i,j+1) ) )
2786  vel = vel * j13g(ke-3,i,j)
2787  flux(ke-3,i,j) = vel / mapf(i,j,+2) &
2788  * ( ( 4.0_rp * val(ke,i,j) &
2789  - 38.0_rp * val(ke-1,i,j) &
2790  + 214.0_rp * val(ke-2,i,j) &
2791  + 319.0_rp * val(ke-3,i,j) &
2792  - 101.0_rp * val(ke-4,i,j) &
2793  + 25.0_rp * val(ke-5,i,j) &
2794  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
2795  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2796  + ( - 3.0_rp * val(ke-4,i,j) &
2797  + 27.0_rp * val(ke-3,i,j) &
2798  + 47.0_rp * val(ke-2,i,j) &
2799  - 13.0_rp * val(ke-1,i,j) &
2800  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
2801  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2802 
2803  flux(ke ,i,j) = 0.0_rp
2804  enddo
2805  enddo
2806  !$acc end kernels
2807  !$omp end do nowait
2808 
2809 
2810 
2811  !$acc end data
2812 
2813  !$omp end parallel
2814  return
2815  end subroutine atmos_dyn_fvm_fluxj13_xvz_ud7
2816 
2817  !-----------------------------------------------------------------------------
2819  subroutine atmos_dyn_fvm_fluxj23_xvz_ud7( &
2820  flux, &
2821  mom, val, DENS, &
2822  GSQRT, J23G, MAPF, &
2823  CDZ, TwoD, &
2824  IIS, IIE, JJS, JJE )
2825  implicit none
2826 
2827  real(rp), intent(inout) :: flux (ka,ia,ja)
2828  real(rp), intent(in) :: mom (ka,ia,ja)
2829  real(rp), intent(in) :: val (ka,ia,ja)
2830  real(rp), intent(in) :: dens (ka,ia,ja)
2831  real(rp), intent(in) :: gsqrt (ka,ia,ja)
2832  real(rp), intent(in) :: j23g (ka,ia,ja)
2833  real(rp), intent(in) :: mapf ( ia,ja,2)
2834  real(rp), intent(in) :: cdz (ka)
2835  logical, intent(in) :: twod
2836  integer, intent(in) :: iis, iie, jjs, jje
2837 
2838  real(rp) :: vel
2839  integer :: k, i, j
2840  !---------------------------------------------------------------------------
2841 
2842  !$omp parallel default(none) private(i,j,k,vel) &
2843  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,J23G,MAPF) &
2844  !$omp shared(GSQRT,CDZ,TwoD)
2845 
2846  !$acc data copy(flux) copyin(mom, val, DENS, GSQRT, J23G, MAPF, CDZ)
2847 
2848 
2849 
2850  !$omp do OMP_SCHEDULE_ collapse(2)
2851  !$acc kernels
2852  do j = jjs, jje
2853  do i = iis, iie
2854  do k = ks+3, ke-4
2855  vel = ( f2h(k,1,i_xvz) &
2856  * mom(k+1,i,j) &
2857  + f2h(k,2,i_xvz) &
2858  * mom(k,i,j) ) &
2859  / ( f2h(k,1,i_xvz) &
2860  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
2861  + f2h(k,2,i_xvz) &
2862  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
2863  vel = vel * j23g(k,i,j)
2864  flux(k,i,j) = vel / mapf(i,j,+1) &
2865  * ( ( f71 * ( val(k+4,i,j)+val(k-3,i,j) ) &
2866  + f72 * ( val(k+3,i,j)+val(k-2,i,j) ) &
2867  + f73 * ( val(k+2,i,j)+val(k-1,i,j) ) &
2868  + f74 * ( val(k+1,i,j)+val(k,i,j) ) ) &
2869  - ( f71 * ( val(k+4,i,j)-val(k-3,i,j) ) &
2870  + f75 * ( val(k+3,i,j)-val(k-2,i,j) ) &
2871  + f76 * ( val(k+2,i,j)-val(k-1,i,j) ) &
2872  + f77 * ( val(k+1,i,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) )
2873  enddo
2874  enddo
2875  enddo
2876  !$acc end kernels
2877  !$omp end do nowait
2878 
2879  !$omp do OMP_SCHEDULE_ collapse(2)
2880  !$acc kernels
2881  do j = jjs, jje
2882  do i = iis, iie
2883  ! The boundary condition is qflx_hi + qflxJ13 + qfluxJ23 = 0 at KS-1.
2884  ! The flux at KS-1 can be non-zero.
2885  ! To reduce calculations, all the fluxes are set to zero.
2886  flux(ks-1,i,j) = 0.0_rp
2887 
2888  vel = ( f2h(ks,1,i_xvz) &
2889  * mom(ks+1,i,j) &
2890  + f2h(ks,2,i_xvz) &
2891  * mom(ks,i,j) ) &
2892  / ( f2h(ks,1,i_xvz) &
2893  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
2894  + f2h(ks,2,i_xvz) &
2895  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
2896  vel = vel * j23g(ks,i,j)
2897  flux(ks,i,j) = vel / mapf(i,j,+1) &
2898  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) &
2899  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2900  + ( 2.0_rp * val(ks,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks+2,i,j) ) / 6.0_rp &
2901  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2902 
2903  vel = ( f2h(ke-1,1,i_xvz) &
2904  * mom(ke,i,j) &
2905  + f2h(ke-1,2,i_xvz) &
2906  * mom(ke-1,i,j) ) &
2907  / ( f2h(ke-1,1,i_xvz) &
2908  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
2909  + f2h(ke-1,2,i_xvz) &
2910  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
2911  vel = vel * j23g(ke-1,i,j)
2912  flux(ke-1,i,j) = vel / mapf(i,j,+1) &
2913  * ( ( 2.0_rp * val(ke,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke-2,i,j) ) / 6.0_rp &
2914  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2915  + f2 * ( val(ke,i,j)+val(ke-1,i,j) ) &
2916  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2917 
2918  vel = ( f2h(ks+1,1,i_xvz) &
2919  * mom(ks+2,i,j) &
2920  + f2h(ks+1,2,i_xvz) &
2921  * mom(ks+1,i,j) ) &
2922  / ( f2h(ks+1,1,i_xvz) &
2923  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) &
2924  + f2h(ks+1,2,i_xvz) &
2925  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) )
2926  vel = vel * j23g(ks+1,i,j)
2927  flux(ks+1,i,j) = vel / mapf(i,j,+1) &
2928  * ( ( 2.0_rp * val(ks+2,i,j) + 5.0_rp * val(ks+1,i,j) - val(ks,i,j) ) / 6.0_rp &
2929  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2930  + ( - 3.0_rp * val(ks,i,j) &
2931  + 27.0_rp * val(ks+1,i,j) &
2932  + 47.0_rp * val(ks+2,i,j) &
2933  - 13.0_rp * val(ks+3,i,j) &
2934  + 2.0_rp * val(ks+4,i,j) ) / 60.0_rp &
2935  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2936 
2937  vel = ( f2h(ke-2,1,i_xvz) &
2938  * mom(ke-1,i,j) &
2939  + f2h(ke-2,2,i_xvz) &
2940  * mom(ke-2,i,j) ) &
2941  / ( f2h(ke-2,1,i_xvz) &
2942  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) &
2943  + f2h(ke-2,2,i_xvz) &
2944  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) )
2945  vel = vel * j23g(ke-2,i,j)
2946  flux(ke-2,i,j) = vel / mapf(i,j,+1) &
2947  * ( ( - 3.0_rp * val(ke,i,j) &
2948  + 27.0_rp * val(ke-1,i,j) &
2949  + 47.0_rp * val(ke-2,i,j) &
2950  - 13.0_rp * val(ke-3,i,j) &
2951  + 2.0_rp * val(ke-4,i,j) ) / 60.0_rp &
2952  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2953  + ( 2.0_rp * val(ke-2,i,j) + 5.0_rp * val(ke-1,i,j) - val(ke,i,j) ) / 6.0_rp &
2954  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2955 
2956  vel = ( f2h(ks+2,1,i_xvz) &
2957  * mom(ks+3,i,j) &
2958  + f2h(ks+2,2,i_xvz) &
2959  * mom(ks+2,i,j) ) &
2960  / ( f2h(ks+2,1,i_xvz) &
2961  * 0.5_rp * ( dens(ks+3,i,j)+dens(ks+3,i,j+1) ) &
2962  + f2h(ks+2,2,i_xvz) &
2963  * 0.5_rp * ( dens(ks+2,i,j)+dens(ks+2,i,j+1) ) )
2964  vel = vel * j23g(ks+2,i,j)
2965  flux(ks+2,i,j) = vel / mapf(i,j,+1) &
2966  * ( ( - 3.0_rp * val(ks+4,i,j) &
2967  + 27.0_rp * val(ks+3,i,j) &
2968  + 47.0_rp * val(ks+2,i,j) &
2969  - 13.0_rp * val(ks+1,i,j) &
2970  + 2.0_rp * val(ks,i,j) ) / 60.0_rp &
2971  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2972  + ( 4.0_rp * val(ks,i,j) &
2973  - 38.0_rp * val(ks+1,i,j) &
2974  + 214.0_rp * val(ks+2,i,j) &
2975  + 319.0_rp * val(ks+3,i,j) &
2976  - 101.0_rp * val(ks+4,i,j) &
2977  + 25.0_rp * val(ks+5,i,j) &
2978  - 3.0_rp * val(ks+6,i,j) ) / 420.0_rp &
2979  * ( 0.5_rp - sign(0.5_rp,vel) ) )
2980 
2981  vel = ( f2h(ke-3,1,i_xvz) &
2982  * mom(ke-2,i,j) &
2983  + f2h(ke-3,2,i_xvz) &
2984  * mom(ke-3,i,j) ) &
2985  / ( f2h(ke-3,1,i_xvz) &
2986  * 0.5_rp * ( dens(ke-2,i,j)+dens(ke-2,i,j+1) ) &
2987  + f2h(ke-3,2,i_xvz) &
2988  * 0.5_rp * ( dens(ke-3,i,j)+dens(ke-3,i,j+1) ) )
2989  vel = vel * j23g(ke-3,i,j)
2990  flux(ke-3,i,j) = vel / mapf(i,j,+1) &
2991  * ( ( 4.0_rp * val(ke,i,j) &
2992  - 38.0_rp * val(ke-1,i,j) &
2993  + 214.0_rp * val(ke-2,i,j) &
2994  + 319.0_rp * val(ke-3,i,j) &
2995  - 101.0_rp * val(ke-4,i,j) &
2996  + 25.0_rp * val(ke-5,i,j) &
2997  - 3.0_rp * val(ke-6,i,j) ) / 420.0_rp &
2998  * ( 0.5_rp + sign(0.5_rp,vel) ) &
2999  + ( - 3.0_rp * val(ke-4,i,j) &
3000  + 27.0_rp * val(ke-3,i,j) &
3001  + 47.0_rp * val(ke-2,i,j) &
3002  - 13.0_rp * val(ke-1,i,j) &
3003  + 2.0_rp * val(ke,i,j) ) / 60.0_rp &
3004  * ( 0.5_rp - sign(0.5_rp,vel) ) )
3005 
3006  flux(ke ,i,j) = 0.0_rp
3007  enddo
3008  enddo
3009  !$acc end kernels
3010  !$omp end do nowait
3011 
3012 
3013 
3014  !$acc end data
3015 
3016  !$omp end parallel
3017  return
3018  end subroutine atmos_dyn_fvm_fluxj23_xvz_ud7
3019 
3020  !-----------------------------------------------------------------------------
3022  subroutine atmos_dyn_fvm_fluxx_xvz_ud7( &
3023  flux, &
3024  mom, val, DENS, &
3025  GSQRT, MAPF, &
3026  num_diff, &
3027  CDZ, TwoD, &
3028  IIS, IIE, JJS, JJE )
3029  implicit none
3030 
3031  real(rp), intent(inout) :: flux (ka,ia,ja)
3032  real(rp), intent(in) :: mom (ka,ia,ja)
3033  real(rp), intent(in) :: val (ka,ia,ja)
3034  real(rp), intent(in) :: dens (ka,ia,ja)
3035  real(rp), intent(in) :: gsqrt (ka,ia,ja)
3036  real(rp), intent(in) :: mapf ( ia,ja,2)
3037  real(rp), intent(in) :: num_diff(ka,ia,ja)
3038  real(rp), intent(in) :: cdz (ka)
3039  logical, intent(in) :: twod
3040  integer, intent(in) :: iis, iie, jjs, jje
3041 
3042  real(rp) :: vel
3043  integer :: k, i, j
3044  !---------------------------------------------------------------------------
3045 
3046 
3047 
3048  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
3049  !$omp private(vel) &
3050  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff)
3051  !$acc kernels
3052  do j = jjs, jje
3053  do i = iis-1, iie
3054  do k = ks, ke
3055 #ifdef DEBUG
3056  call check( __line__, mom(k,i ,j) )
3057  call check( __line__, mom(k,i,j-1) )
3058 
3059  call check( __line__, val(k,i,j) )
3060  call check( __line__, val(k,i+1,j) )
3061 
3062  call check( __line__, val(k,i-1,j) )
3063  call check( __line__, val(k,i+2,j) )
3064 
3065  call check( __line__, val(k,i-2,j) )
3066  call check( __line__, val(k,i+3,j) )
3067 
3068  call check( __line__, val(k,i-3,j) )
3069  call check( __line__, val(k,i+4,j) )
3070 
3071 #endif
3072  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j+1) ) ) &
3073  / ( 0.25_rp * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) )
3074  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
3075  * ( ( f71 * ( val(k,i+4,j)+val(k,i-3,j) ) &
3076  + f72 * ( val(k,i+3,j)+val(k,i-2,j) ) &
3077  + f73 * ( val(k,i+2,j)+val(k,i-1,j) ) &
3078  + f74 * ( val(k,i+1,j)+val(k,i,j) ) ) &
3079  - ( f71 * ( val(k,i+4,j)-val(k,i-3,j) ) &
3080  + f75 * ( val(k,i+3,j)-val(k,i-2,j) ) &
3081  + f76 * ( val(k,i+2,j)-val(k,i-1,j) ) &
3082  + f77 * ( val(k,i+1,j)-val(k,i,j) ) ) * sign(1.0_rp,vel) ) &
3083  + gsqrt(k,i,j) * num_diff(k,i,j)
3084  enddo
3085  enddo
3086  enddo
3087  !$acc end kernels
3088 #ifdef DEBUG
3089  k = iundef; i = iundef; j = iundef
3090 #endif
3091 
3092 
3093 
3094  return
3095  end subroutine atmos_dyn_fvm_fluxx_xvz_ud7
3096 
3097  !-----------------------------------------------------------------------------
3099  subroutine atmos_dyn_fvm_fluxy_xvz_ud7( &
3100  flux, &
3101  mom, val, DENS, &
3102  GSQRT, MAPF, &
3103  num_diff, &
3104  CDZ, TwoD, &
3105  IIS, IIE, JJS, JJE )
3106  implicit none
3107 
3108  real(rp), intent(inout) :: flux (ka,ia,ja)
3109  real(rp), intent(in) :: mom (ka,ia,ja)
3110  real(rp), intent(in) :: val (ka,ia,ja)
3111  real(rp), intent(in) :: dens (ka,ia,ja)
3112  real(rp), intent(in) :: gsqrt (ka,ia,ja)
3113  real(rp), intent(in) :: mapf ( ia,ja,2)
3114  real(rp), intent(in) :: num_diff(ka,ia,ja)
3115  real(rp), intent(in) :: cdz (ka)
3116  logical, intent(in) :: twod
3117  integer, intent(in) :: iis, iie, jjs, jje
3118 
3119  real(rp) :: vel
3120  integer :: k, i, j
3121  !---------------------------------------------------------------------------
3122 
3123  ! note that y-index is added by -1
3124 
3125 
3126 
3127  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
3128  !$omp private(vel) &
3129  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,mom,val,DENS,flux,GSQRT,MAPF,num_diff)
3130  !$acc kernels
3131  do j = jjs, jje+1
3132  do i = iis, iie
3133  do k = ks, ke
3134 #ifdef DEBUG
3135  call check( __line__, mom(k,i ,j) )
3136  call check( __line__, mom(k,i,j-1) )
3137 
3138  call check( __line__, val(k,i,j-1) )
3139  call check( __line__, val(k,i,j) )
3140 
3141  call check( __line__, val(k,i,j-2) )
3142  call check( __line__, val(k,i,j+1) )
3143 
3144  call check( __line__, val(k,i,j-3) )
3145  call check( __line__, val(k,i,j+2) )
3146 
3147  call check( __line__, val(k,i,j-4) )
3148  call check( __line__, val(k,i,j+3) )
3149 
3150 #endif
3151  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j-1) ) ) &
3152  / ( dens(k,i,j) )
3153  flux(k,i,j-1) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
3154  * ( ( f71 * ( val(k,i,j+3)+val(k,i,j-4) ) &
3155  + f72 * ( val(k,i,j+2)+val(k,i,j-3) ) &
3156  + f73 * ( val(k,i,j+1)+val(k,i,j-2) ) &
3157  + f74 * ( val(k,i,j)+val(k,i,j-1) ) ) &
3158  - ( f71 * ( val(k,i,j+3)-val(k,i,j-4) ) &
3159  + f75 * ( val(k,i,j+2)-val(k,i,j-3) ) &
3160  + f76 * ( val(k,i,j+1)-val(k,i,j-2) ) &
3161  + f77 * ( val(k,i,j)-val(k,i,j-1) ) ) * sign(1.0_rp,vel) ) &
3162  + gsqrt(k,i,j) * num_diff(k,i,j)
3163  enddo
3164  enddo
3165  enddo
3166  !$acc end kernels
3167 #ifdef DEBUG
3168  k = iundef; i = iundef; j = iundef
3169 #endif
3170 
3171 
3172 
3173  return
3174  end subroutine atmos_dyn_fvm_fluxy_xvz_ud7
3175 
3176 
3177 
3178 
3179 
3180 
3181 
3182 end module scale_atmos_dyn_fvm_flux_ud7
3183 
3184 !--
3185 ! vi:set readonly sw=4 ts=8
3186 !
3187 !Local Variables:
3188 !mode: f90
3189 !buffer-read-only: t
3190 !End:
3191 !
3192 !++
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_ud7::atmos_dyn_fvm_fluxy_xvz_ud7
subroutine, public atmos_dyn_fvm_fluxy_xvz_ud7(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_ud7.F90:3106
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxz_xyz_ud7
subroutine, public atmos_dyn_fvm_fluxz_xyz_ud7(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation z-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:251
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:40
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxx_uyz_ud7
subroutine, public atmos_dyn_fvm_fluxx_uyz_ud7(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_ud7.F90:2196
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxj23_uyz_ud7
subroutine, public atmos_dyn_fvm_fluxj23_uyz_ud7(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J23-flux at UYZ
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:1827
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_ud7
module scale_atmos_dyn_fvm_flux_ud7
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:16
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:35
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxz_uyz_ud7
subroutine, public atmos_dyn_fvm_fluxz_uyz_ud7(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_ud7.F90:1228
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxy_xyw_ud7
subroutine, public atmos_dyn_fvm_fluxy_xyw_ud7(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_ud7.F90:1126
scale_atmos_grid_cartesc_index::i_xyz
integer, public i_xyz
Definition: scale_atmos_grid_cartesC_index.F90:91
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxy_uyz_ud7
subroutine, public atmos_dyn_fvm_fluxy_uyz_ud7(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_ud7.F90:2275
scale_atmos_grid_cartesc_index::i_uyz
integer, public i_uyz
Definition: scale_atmos_grid_cartesC_index.F90:95
scale_io
module STDIO
Definition: scale_io.F90:10
scale_tracer::k
real(rp), public k
Definition: scale_tracer.F90:45
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_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxj23_xyw_ud7
subroutine, public atmos_dyn_fvm_fluxj23_xyw_ud7(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J23-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:900
scale_debug::check
subroutine, public check(current_line, v)
Undefined value checker.
Definition: scale_debug.F90:59
scale_atmos_grid_cartesc_index::i_xvz
integer, public i_xvz
Definition: scale_atmos_grid_cartesC_index.F90:96
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_flux_valuew_z_ud7
subroutine, public atmos_dyn_fvm_flux_valuew_z_ud7(valW, mflx, val, GSQRT, CDZ)
value at XYW
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:117
scale_prof
module profiler
Definition: scale_prof.F90:11
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxj13_uyz_ud7
subroutine, public atmos_dyn_fvm_fluxj13_uyz_ud7(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J13-flux at UYZ
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:1624
scale_atmos_grid_cartesc_index::ja
integer, public ja
Definition: scale_atmos_grid_cartesC_index.F90:49
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxz_xvz_ud7
subroutine, public atmos_dyn_fvm_fluxz_xvz_ud7(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_ud7.F90:2403
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxj13_xvz_ud7
subroutine, public atmos_dyn_fvm_fluxj13_xvz_ud7(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J13-flux at XVZ
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:2622
scale_tracer
module TRACER
Definition: scale_tracer.F90:12
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_debug
module DEBUG
Definition: scale_debug.F90:11
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxz_xyw_ud7
subroutine, public atmos_dyn_fvm_fluxz_xyw_ud7(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_ud7.F90:577
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxx_xyw_ud7
subroutine, public atmos_dyn_fvm_fluxx_xyw_ud7(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_ud7.F90:1025
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxy_xyz_ud7
subroutine, public atmos_dyn_fvm_fluxy_xyz_ud7(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation Y-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:507
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxj13_xyw_ud7
subroutine, public atmos_dyn_fvm_fluxj13_xyw_ud7(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J13-flux at XYW
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:777
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:43
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxj23_xvz_ud7
subroutine, public atmos_dyn_fvm_fluxj23_xvz_ud7(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, TwoD, IIS, IIE, JJS, JJE)
calculation J23-flux at XVZ
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:2825
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxx_xvz_ud7
subroutine, public atmos_dyn_fvm_fluxx_xvz_ud7(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_ud7.F90:3029
scale_atmos_dyn_fvm_flux_ud7::atmos_dyn_fvm_fluxx_xyz_ud7
subroutine, public atmos_dyn_fvm_fluxx_xyz_ud7(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation X-flux at XYZ
Definition: scale_atmos_dyn_fvm_flux_ud7.F90:440
Generated by   doxygen 1.8.17