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