SCALE-RM
scale_atmos_dyn_fvm_flux_ud3Koren1993.F90
Go to the documentation of this file.
1 
2 !-------------------------------------------------------------------------------
14 !-------------------------------------------------------------------------------
15 ! Warning: This file was generated from atmos-rm/dynamics/scale_atmos_dyn_fvm_flux_udcd.F90.erb.
16 ! Do not edit this file.
17 !-------------------------------------------------------------------------------
19  !-----------------------------------------------------------------------------
20  !
21  !++ used modules
22  !
23  use scale_precision
24  use scale_stdio
25  use scale_prof
27  use scale_index
28  use scale_tracer
29  use scale_process
30 #ifdef DEBUG
31  use scale_debug, only: &
32  check
33  use scale_const, only: &
34  undef => const_undef, &
35  iundef => const_undef2
36 #endif
37  !-----------------------------------------------------------------------------
38  implicit none
39  private
40  !-----------------------------------------------------------------------------
41  !
42  !++ Public procedure
43  !
45 
49 
55 
61 
67 
68  !-----------------------------------------------------------------------------
69  !
70  !++ Public parameters & variables
71  !
72  !-----------------------------------------------------------------------------
73  !
74  !++ Private procedure
75  !
76 #if 1
77 #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)))
78 #else
79 #define F2H(k,p,q) 0.5_RP
80 #endif
81  !-----------------------------------------------------------------------------
82  !
83  !++ Private parameters & variables
84  !
85 
86  real(RP), parameter :: f2 = 0.5_rp
87 
88 
89 
90  real(RP), parameter :: f31 = -1.0_rp/12.0_rp
91  real(RP), parameter :: f32 = 7.0_rp/12.0_rp
92  real(RP), parameter :: f33 = 3.0_rp/12.0_rp
93 
94 
95 
96 
97 contains
98  !-----------------------------------------------------------------------------
100 !OCL SERIAL
102  valW, &
103  mflx, val, GSQRT, &
104  CDZ )
105  implicit none
106 
107  real(RP), intent(out) :: valW (ka)
108  real(RP), intent(in) :: mflx (ka)
109  real(RP), intent(in) :: val (ka)
110  real(RP), intent(in) :: GSQRT(ka)
111  real(RP), intent(in) :: CDZ (ka)
112 
113  integer :: k
114  !---------------------------------------------------------------------------
115 
116  do k = ks+1, ke-2
117 #ifdef DEBUG
118  call check( __line__, mflx(k) )
119 
120  call check( __line__, val(k) )
121  call check( __line__, val(k+1) )
122 
123  call check( __line__, val(k-1) )
124  call check( __line__, val(k+2) )
125 
126 #endif
127  valw(k) = ( val(k) &
128  + 0.5_rp * phi(val(k+1),val(k),val(k-1)) * ( val(k)-val(k-1) ) ) &
129  * ( 0.5_rp + sign(0.5_rp,mflx(k)) ) &
130  + ( val(k+1) &
131  + 0.5_rp * phi(val(k),val(k+1),val(k+2)) * ( val(k+1)-val(k+2) ) ) &
132  * ( 0.5_rp - sign(0.5_rp,mflx(k)) )
133  enddo
134 #ifdef DEBUG
135  k = iundef
136 #endif
137 
138 #ifdef DEBUG
139 
140  call check( __line__, mflx(ks) )
141  call check( __line__, val(ks ) )
142  call check( __line__, val(ks+1) )
143  call check( __line__, mflx(ke-1) )
144  call check( __line__, val(ke ) )
145  call check( __line__, val(ke-1) )
146 
147 #endif
148 
149  valw(ks) = f2 * ( val(ks+1)+val(ks) )
150  valw(ke-1) = f2 * ( val(ke)+val(ke-1) )
151 
152  return
154 
155  !-----------------------------------------------------------------------------
158  flux, &
159  mflx, val, GSQRT, &
160 
161  num_diff, &
162 
163  CDZ, &
164  IIS, IIE, JJS, JJE )
165  implicit none
166 
167  real(RP), intent(out) :: flux (ka,ia,ja)
168  real(RP), intent(in) :: mflx (ka,ia,ja)
169  real(RP), intent(in) :: val (ka,ia,ja)
170  real(RP), intent(in) :: GSQRT (ka,ia,ja)
171 
172  real(RP), intent(in) :: num_diff(ka,ia,ja)
173 
174  real(RP), intent(in) :: CDZ (ka)
175  integer, intent(in) :: IIS, IIE, JJS, JJE
176 
177  real(RP) :: vel
178  integer :: k, i, j
179  !---------------------------------------------------------------------------
180 
181  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
182  do j = jjs, jje
183  do i = iis, iie
184  do k = ks+1, ke-2
185 #ifdef DEBUG
186  call check( __line__, mflx(k,i,j) )
187 
188  call check( __line__, val(k,i,j) )
189  call check( __line__, val(k+1,i,j) )
190 
191  call check( __line__, val(k-1,i,j) )
192  call check( __line__, val(k+2,i,j) )
193 
194 #endif
195  vel = mflx(k,i,j)
196  flux(k,i,j) = vel &
197  * ( ( val(k,i,j) &
198  + 0.5_rp * phi(val(k+1,i,j),val(k,i,j),val(k-1,i,j)) * ( val(k,i,j)-val(k-1,i,j) ) ) &
199  * ( 0.5_rp + sign(0.5_rp,vel) ) &
200  + ( val(k+1,i,j) &
201  + 0.5_rp * phi(val(k,i,j),val(k+1,i,j),val(k+2,i,j)) * ( val(k+1,i,j)-val(k+2,i,j) ) ) &
202  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
203  + gsqrt(k,i,j) * num_diff(k,i,j)
204  enddo
205  enddo
206  enddo
207 #ifdef DEBUG
208  k = iundef; i = iundef; j = iundef
209 #endif
210 
211  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
212  do j = jjs, jje
213  do i = iis, iie
214 #ifdef DEBUG
215 
216  call check( __line__, mflx(ks,i,j) )
217  call check( __line__, val(ks ,i,j) )
218  call check( __line__, val(ks+1,i,j) )
219  call check( __line__, mflx(ke-1,i,j) )
220  call check( __line__, val(ke ,i,j) )
221  call check( __line__, val(ke-1,i,j) )
222 
223 #endif
224  flux(ks-1,i,j) = 0.0_rp
225 
226  vel = mflx(ks,i,j)
227  flux(ks,i,j) = vel &
228  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) &
229  + gsqrt(ks,i,j) * num_diff(ks,i,j)
230  vel = mflx(ke-1,i,j)
231  flux(ke-1,i,j) = vel &
232  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) ) &
233  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
234 
235  flux(ke ,i,j) = 0.0_rp
236  enddo
237  enddo
238 #ifdef DEBUG
239  k = iundef; i = iundef; j = iundef
240 #endif
241 
242  return
244 
245  !-----------------------------------------------------------------------------
248  flux, &
249  mflx, val, GSQRT, &
250 
251  num_diff, &
252 
253  CDZ, &
254  IIS, IIE, JJS, JJE )
255  implicit none
256 
257  real(RP), intent(out) :: flux (ka,ia,ja)
258  real(RP), intent(in) :: mflx (ka,ia,ja)
259  real(RP), intent(in) :: val (ka,ia,ja)
260  real(RP), intent(in) :: GSQRT (ka,ia,ja)
261 
262  real(RP), intent(in) :: num_diff(ka,ia,ja)
263 
264  real(RP), intent(in) :: CDZ(ka)
265  integer, intent(in) :: IIS, IIE, JJS, JJE
266 
267  real(RP) :: vel
268  integer :: k, i, j
269  !---------------------------------------------------------------------------
270 
271  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
272  do j = jjs, jje
273  do i = iis-1, iie
274  do k = ks, ke
275 #ifdef DEBUG
276  call check( __line__, mflx(k,i,j) )
277 
278  call check( __line__, val(k,i,j) )
279  call check( __line__, val(k,i+1,j) )
280 
281  call check( __line__, val(k,i-1,j) )
282  call check( __line__, val(k,i+2,j) )
283 
284 #endif
285  vel = mflx(k,i,j)
286  flux(k,i,j) = vel &
287  * ( ( val(k,i,j) &
288  + 0.5_rp * phi(val(k,i+1,j),val(k,i,j),val(k,i-1,j)) * ( val(k,i,j)-val(k,i-1,j) ) ) &
289  * ( 0.5_rp + sign(0.5_rp,vel) ) &
290  + ( val(k,i+1,j) &
291  + 0.5_rp * phi(val(k,i,j),val(k,i+1,j),val(k,i+2,j)) * ( val(k,i+1,j)-val(k,i+2,j) ) ) &
292  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
293  + gsqrt(k,i,j) * num_diff(k,i,j)
294  enddo
295  enddo
296  enddo
297 #ifdef DEBUG
298  k = iundef; i = iundef; j = iundef
299 #endif
300 
301  return
303 
304  !-----------------------------------------------------------------------------
307  flux, &
308  mflx, val, GSQRT, &
309 
310  num_diff, &
311 
312  CDZ, &
313  IIS, IIE, JJS, JJE )
314  implicit none
315 
316  real(RP), intent(out) :: flux (ka,ia,ja)
317  real(RP), intent(in) :: mflx (ka,ia,ja)
318  real(RP), intent(in) :: val (ka,ia,ja)
319  real(RP), intent(in) :: GSQRT (ka,ia,ja)
320 
321  real(RP), intent(in) :: num_diff(ka,ia,ja)
322 
323  real(RP), intent(in) :: CDZ(ka)
324  integer, intent(in) :: IIS, IIE, JJS, JJE
325 
326  real(RP) :: vel
327  integer :: k, i, j
328  !---------------------------------------------------------------------------
329 
330  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
331  do j = jjs-1, jje
332  do i = iis, iie
333  do k = ks, ke
334 #ifdef DEBUG
335  call check( __line__, mflx(k,i,j) )
336 
337  call check( __line__, val(k,i,j) )
338  call check( __line__, val(k,i,j+1) )
339 
340  call check( __line__, val(k,i,j-1) )
341  call check( __line__, val(k,i,j+2) )
342 
343 #endif
344  vel = mflx(k,i,j)
345  flux(k,i,j) = vel &
346  * ( ( val(k,i,j) &
347  + 0.5_rp * phi(val(k,i,j+1),val(k,i,j),val(k,i,j-1)) * ( val(k,i,j)-val(k,i,j-1) ) ) &
348  * ( 0.5_rp + sign(0.5_rp,vel) ) &
349  + ( val(k,i,j+1) &
350  + 0.5_rp * phi(val(k,i,j),val(k,i,j+1),val(k,i,j+2)) * ( val(k,i,j+1)-val(k,i,j+2) ) ) &
351  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
352  + gsqrt(k,i,j) * num_diff(k,i,j)
353  enddo
354  enddo
355  enddo
356 #ifdef DEBUG
357  k = iundef; i = iundef; j = iundef
358 #endif
359 
360  return
362 
363 
364  !-----------------------------------------------------------------------------
367  flux, &
368  mom, val, DENS, &
369  GSQRT, J33G, &
370 
371  num_diff, &
372 
373  CDZ, FDZ, &
374  dtrk, &
375  IIS, IIE, JJS, JJE )
376  implicit none
377 
378  real(RP), intent(out) :: flux (ka,ia,ja)
379  real(RP), intent(in) :: mom (ka,ia,ja)
380  real(RP), intent(in) :: val (ka,ia,ja)
381  real(RP), intent(in) :: DENS (ka,ia,ja)
382  real(RP), intent(in) :: GSQRT (ka,ia,ja)
383  real(RP), intent(in) :: J33G
384 
385  real(RP), intent(in) :: num_diff(ka,ia,ja)
386 
387  real(RP), intent(in) :: CDZ (ka)
388  real(RP), intent(in) :: FDZ (ka-1)
389  real(RP), intent(in) :: dtrk
390  integer, intent(in) :: IIS, IIE, JJS, JJE
391 
392  real(RP) :: vel
393  real(RP) :: sw
394  integer :: k, i, j
395  !---------------------------------------------------------------------------
396 
397  ! note than z-index is added by -1
398 
399  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
400  do j = jjs, jje
401  do i = iis, iie
402  do k = ks+1, ke-1
403 #ifdef DEBUG
404  call check( __line__, mom(k-1,i,j) )
405  call check( __line__, mom(k ,i,j) )
406 
407  call check( __line__, val(k-1,i,j) )
408  call check( __line__, val(k,i,j) )
409 
410  call check( __line__, val(k-2,i,j) )
411  call check( __line__, val(k+1,i,j) )
412 
413 #endif
414  vel = ( 0.5_rp * ( mom(k-1,i,j) &
415  + mom(k,i,j) ) ) &
416  / dens(k,i,j)
417  flux(k-1,i,j) = j33g * vel &
418  * ( ( val(k-1,i,j) &
419  + 0.5_rp * phi(val(k,i,j),val(k-1,i,j),val(k-2,i,j)) * ( val(k-1,i,j)-val(k-2,i,j) ) ) &
420  * ( 0.5_rp + sign(0.5_rp,vel) ) &
421  + ( val(k,i,j) &
422  + 0.5_rp * phi(val(k-1,i,j),val(k,i,j),val(k+1,i,j)) * ( val(k,i,j)-val(k+1,i,j) ) ) &
423  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
424  + gsqrt(k,i,j) * num_diff(k,i,j)
425  enddo
426  enddo
427  enddo
428 #ifdef DEBUG
429  k = iundef; i = iundef; j = iundef
430 #endif
431 
432  !$omp parallel do private(i,j,sw) OMP_SCHEDULE_ collapse(2)
433  do j = jjs, jje
434  do i = iis, iie
435 #ifdef DEBUG
436 
437 #endif
438  flux(ks-1,i,j) = 0.0_rp ! k = KS
439 
440  ! if w>0; min(f,w*dz/dt)
441  ! else ; max(f,w*dz/dt) = -min(-f,-w*dz/dt)
442  sw = sign( 1.0_rp, mom(ks,i,j) )
443  flux(ks ,i,j) = sw * min( sw*flux(ks,i,j), sw*val(ks,i,j)*gsqrt(ks,i,j)*fdz(ks)/dtrk )
444 
445 
446  flux(ke-1,i,j) = 0.0_rp ! k = KE
447  flux(ke ,i,j) = 0.0_rp ! k = KE+1
448  enddo
449  enddo
450 
451  return
453 
454 
455  !-----------------------------------------------------------------------------
458  flux, &
459  mom, val, DENS, &
460  GSQRT, J13G, MAPF, &
461  CDZ, &
462  IIS, IIE, JJS, JJE )
463  implicit none
464 
465  real(RP), intent(out) :: flux (ka,ia,ja)
466  real(RP), intent(in) :: mom (ka,ia,ja)
467  real(RP), intent(in) :: val (ka,ia,ja)
468  real(RP), intent(in) :: DENS (ka,ia,ja)
469  real(RP), intent(in) :: GSQRT (ka,ia,ja)
470  real(RP), intent(in) :: J13G (ka,ia,ja)
471  real(RP), intent(in) :: MAPF ( ia,ja,2)
472  real(RP), intent(in) :: CDZ (ka)
473  integer, intent(in) :: IIS, IIE, JJS, JJE
474 
475  real(RP) :: vel
476  integer :: k, i, j
477  !---------------------------------------------------------------------------
478 
479  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
480  do j = jjs, jje
481  do i = iis, iie
482  do k = ks+1, ke-1
483  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i-1,j) ) ) &
484  / dens(k,i,j)
485  flux(k-1,i,j) = j13g(k,i,j) / mapf(i,j,+2) * vel &
486  * ( ( val(k-1,i,j) &
487  + 0.5_rp * phi(val(k,i,j),val(k-1,i,j),val(k-2,i,j)) * ( val(k-1,i,j)-val(k-2,i,j) ) ) &
488  * ( 0.5_rp + sign(0.5_rp,vel) ) &
489  + ( val(k,i,j) &
490  + 0.5_rp * phi(val(k-1,i,j),val(k,i,j),val(k+1,i,j)) * ( val(k,i,j)-val(k+1,i,j) ) ) &
491  * ( 0.5_rp - sign(0.5_rp,vel) ) )
492  enddo
493  enddo
494  enddo
495 
496  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
497  do j = jjs, jje
498  do i = iis, iie
499  flux(ks-1,i,j) = 0.0_rp
500 
501  flux(ke-1,i,j) = 0.0_rp
502  enddo
503  enddo
504 
505  return
507 
508  !-----------------------------------------------------------------------------
511  flux, &
512  mom, val, DENS, &
513  GSQRT, J23G, MAPF, &
514  CDZ, &
515  IIS, IIE, JJS, JJE )
516  implicit none
517 
518  real(RP), intent(out) :: flux (ka,ia,ja)
519  real(RP), intent(in) :: mom (ka,ia,ja)
520  real(RP), intent(in) :: val (ka,ia,ja)
521  real(RP), intent(in) :: DENS (ka,ia,ja)
522  real(RP), intent(in) :: GSQRT (ka,ia,ja)
523  real(RP), intent(in) :: J23G (ka,ia,ja)
524  real(RP), intent(in) :: MAPF ( ia,ja,2)
525  real(RP), intent(in) :: CDZ (ka)
526  integer, intent(in) :: IIS, IIE, JJS, JJE
527 
528  real(RP) :: vel
529  integer :: k, i, j
530  !---------------------------------------------------------------------------
531 
532  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
533  do j = jjs, jje
534  do i = iis, iie
535  do k = ks+1, ke-1
536  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j-1) ) ) &
537  / dens(k,i,j)
538  flux(k-1,i,j) = j23g(k,i,j) / mapf(i,j,+1) * vel &
539  * ( ( val(k-1,i,j) &
540  + 0.5_rp * phi(val(k,i,j),val(k-1,i,j),val(k-2,i,j)) * ( val(k-1,i,j)-val(k-2,i,j) ) ) &
541  * ( 0.5_rp + sign(0.5_rp,vel) ) &
542  + ( val(k,i,j) &
543  + 0.5_rp * phi(val(k-1,i,j),val(k,i,j),val(k+1,i,j)) * ( val(k,i,j)-val(k+1,i,j) ) ) &
544  * ( 0.5_rp - sign(0.5_rp,vel) ) )
545  enddo
546  enddo
547  enddo
548 
549  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
550  do j = jjs, jje
551  do i = iis, iie
552  flux(ks-1,i,j) = 0.0_rp
553 
554  flux(ke-1,i,j) = 0.0_rp
555  enddo
556  enddo
557 
558  return
560 
561 
562  !-----------------------------------------------------------------------------
565  flux, &
566  mom, val, DENS, &
567  GSQRT, MAPF, &
568 
569  num_diff, &
570 
571  CDZ, &
572  IIS, IIE, JJS, JJE )
573  implicit none
574 
575  real(RP), intent(out) :: flux (ka,ia,ja)
576  real(RP), intent(in) :: mom (ka,ia,ja)
577  real(RP), intent(in) :: val (ka,ia,ja)
578  real(RP), intent(in) :: DENS (ka,ia,ja)
579  real(RP), intent(in) :: GSQRT (ka,ia,ja)
580  real(RP), intent(in) :: MAPF ( ia,ja,2)
581 
582  real(RP), intent(in) :: num_diff(ka,ia,ja)
583 
584  real(RP), intent(in) :: CDZ (ka)
585  integer, intent(in) :: IIS, IIE, JJS, JJE
586 
587  real(RP) :: vel
588  integer :: k, i, j
589  !---------------------------------------------------------------------------
590 
591  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
592  do j = jjs, jje
593  do i = iis-1, iie
594  do k = ks, ke-1
595 #ifdef DEBUG
596  call check( __line__, mom(k ,i,j) )
597  call check( __line__, mom(k+1,i,j) )
598 
599  call check( __line__, val(k,i,j) )
600  call check( __line__, val(k,i+1,j) )
601 
602  call check( __line__, val(k,i-1,j) )
603  call check( __line__, val(k,i+2,j) )
604 
605 #endif
606  vel = ( f2h(k,1,i_uyz) &
607  * mom(k+1,i,j) &
608  + f2h(k,2,i_uyz) &
609  * mom(k,i,j) ) &
610  / ( f2h(k,1,i_xyz) &
611  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
612  + f2h(k,2,i_xyz) &
613  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
614  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
615  * ( ( val(k,i,j) &
616  + 0.5_rp * phi(val(k,i+1,j),val(k,i,j),val(k,i-1,j)) * ( val(k,i,j)-val(k,i-1,j) ) ) &
617  * ( 0.5_rp + sign(0.5_rp,vel) ) &
618  + ( val(k,i+1,j) &
619  + 0.5_rp * phi(val(k,i,j),val(k,i+1,j),val(k,i+2,j)) * ( val(k,i+1,j)-val(k,i+2,j) ) ) &
620  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
621  + gsqrt(k,i,j) * num_diff(k,i,j)
622  enddo
623  enddo
624  enddo
625 #ifdef DEBUG
626  k = iundef; i = iundef; j = iundef
627 #endif
628 
629  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
630  do j = jjs, jje
631  do i = iis-1, iie
632  flux(ke,i,j) = 0.0_rp
633  enddo
634  enddo
635 #ifdef DEBUG
636  k = iundef; i = iundef; j = iundef
637 #endif
638 
639  return
641 
642  !-----------------------------------------------------------------------------
645  flux, &
646  mom, val, DENS, &
647  GSQRT, MAPF, &
648 
649  num_diff, &
650 
651  CDZ, &
652  IIS, IIE, JJS, JJE )
653  implicit none
654 
655  real(RP), intent(out) :: flux (ka,ia,ja)
656  real(RP), intent(in) :: mom (ka,ia,ja)
657  real(RP), intent(in) :: val (ka,ia,ja)
658  real(RP), intent(in) :: DENS (ka,ia,ja)
659  real(RP), intent(in) :: GSQRT (ka,ia,ja)
660  real(RP), intent(in) :: MAPF ( ia,ja,2)
661 
662  real(RP), intent(in) :: num_diff(ka,ia,ja)
663 
664  real(RP), intent(in) :: CDZ (ka)
665  integer, intent(in) :: IIS, IIE, JJS, JJE
666 
667  real(RP) :: vel
668  integer :: k, i, j
669  !---------------------------------------------------------------------------
670 
671  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
672  do j = jjs-1, jje
673  do i = iis, iie
674  do k = ks, ke-1
675 #ifdef DEBUG
676  call check( __line__, mom(k ,i,j) )
677  call check( __line__, mom(k+1,i,j) )
678 
679  call check( __line__, val(k,i,j) )
680  call check( __line__, val(k,i,j+1) )
681 
682  call check( __line__, val(k,i,j-1) )
683  call check( __line__, val(k,i,j+2) )
684 
685 #endif
686  vel = ( f2h(k,1,i_xvz) &
687  * mom(k+1,i,j) &
688  + f2h(k,2,i_xvz) &
689  * mom(k,i,j) ) &
690  / ( f2h(k,1,i_xyz) &
691  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
692  + f2h(k,2,i_xyz) &
693  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
694  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
695  * ( ( val(k,i,j) &
696  + 0.5_rp * phi(val(k,i,j+1),val(k,i,j),val(k,i,j-1)) * ( val(k,i,j)-val(k,i,j-1) ) ) &
697  * ( 0.5_rp + sign(0.5_rp,vel) ) &
698  + ( val(k,i,j+1) &
699  + 0.5_rp * phi(val(k,i,j),val(k,i,j+1),val(k,i,j+2)) * ( val(k,i,j+1)-val(k,i,j+2) ) ) &
700  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
701  + gsqrt(k,i,j) * num_diff(k,i,j)
702  enddo
703  enddo
704  enddo
705 #ifdef DEBUG
706  k = iundef; i = iundef; j = iundef
707 #endif
708 
709  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
710  do j = jjs-1, jje
711  do i = iis, iie
712  flux(ke,i,j) = 0.0_rp
713  enddo
714  enddo
715 #ifdef DEBUG
716  k = iundef; i = iundef; j = iundef
717 #endif
718 
719  return
721 
722 
723  !-----------------------------------------------------------------------------
726  flux, &
727  mom, val, DENS, &
728  GSQRT, J33G, &
729 
730  num_diff, &
731 
732  CDZ, &
733  IIS, IIE, JJS, JJE )
734  implicit none
735 
736  real(RP), intent(out) :: flux (ka,ia,ja)
737  real(RP), intent(in) :: mom (ka,ia,ja)
738  real(RP), intent(in) :: val (ka,ia,ja)
739  real(RP), intent(in) :: DENS (ka,ia,ja)
740  real(RP), intent(in) :: GSQRT (ka,ia,ja)
741  real(RP), intent(in) :: J33G
742 
743  real(RP), intent(in) :: num_diff(ka,ia,ja)
744 
745  real(RP), intent(in) :: CDZ (ka)
746  integer, intent(in) :: IIS, IIE, JJS, JJE
747 
748  real(RP) :: vel
749  integer :: k, i, j
750  !---------------------------------------------------------------------------
751 
752  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
753  do j = jjs, jje
754  do i = iis, iie
755  do k = ks+1, ke-2
756 #ifdef DEBUG
757  call check( __line__, mom(k,i,j) )
758  call check( __line__, mom(k,i+1,j) )
759 
760  call check( __line__, val(k,i,j) )
761  call check( __line__, val(k+1,i,j) )
762 
763  call check( __line__, val(k-1,i,j) )
764  call check( __line__, val(k+2,i,j) )
765 
766 #endif
767  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i+1,j) ) ) &
768  / ( f2h(k,1,i_xyz) &
769  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
770  + f2h(k,2,i_xyz) &
771  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
772  flux(k,i,j) = j33g * vel &
773  * ( ( val(k,i,j) &
774  + 0.5_rp * phi(val(k+1,i,j),val(k,i,j),val(k-1,i,j)) * ( val(k,i,j)-val(k-1,i,j) ) ) &
775  * ( 0.5_rp + sign(0.5_rp,vel) ) &
776  + ( val(k+1,i,j) &
777  + 0.5_rp * phi(val(k,i,j),val(k+1,i,j),val(k+2,i,j)) * ( val(k+1,i,j)-val(k+2,i,j) ) ) &
778  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
779  + gsqrt(k,i,j) * num_diff(k,i,j)
780  enddo
781  enddo
782  enddo
783 #ifdef DEBUG
784  k = iundef; i = iundef; j = iundef
785 #endif
786 
787  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
788  do j = jjs, jje
789  do i = iis, iie
790 #ifdef DEBUG
791 
792  call check( __line__, mom(ks,i ,j) )
793  call check( __line__, mom(ks,i+1,j) )
794  call check( __line__, val(ks+1,i,j) )
795  call check( __line__, val(ks,i,j) )
796 
797 #endif
798  flux(ks-1,i,j) = 0.0_rp
799 
800  vel = ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i+1,j) ) ) &
801  / ( f2h(ks,1,i_xyz) &
802  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
803  + f2h(ks,2,i_xyz) &
804  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
805  flux(ks,i,j) = j33g * vel &
806  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) &
807  + gsqrt(ks,i,j) * num_diff(ks,i,j)
808  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i+1,j) ) ) &
809  / ( f2h(ke-1,1,i_xyz) &
810  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
811  + f2h(ke-1,2,i_xyz) &
812  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
813  flux(ke-1,i,j) = j33g * vel &
814  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) ) &
815  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
816 
817  flux(ke,i,j) = 0.0_rp
818  enddo
819  enddo
820 #ifdef DEBUG
821  k = iundef; i = iundef; j = iundef
822 #endif
823 
824  return
826 
827  !-----------------------------------------------------------------------------
830  flux, &
831  mom, val, DENS, &
832  GSQRT, J13G, MAPF, &
833  CDZ, &
834  IIS, IIE, JJS, JJE )
835  implicit none
836 
837  real(RP), intent(out) :: flux (ka,ia,ja)
838  real(RP), intent(in) :: mom (ka,ia,ja)
839  real(RP), intent(in) :: val (ka,ia,ja)
840  real(RP), intent(in) :: DENS (ka,ia,ja)
841  real(RP), intent(in) :: GSQRT (ka,ia,ja)
842  real(RP), intent(in) :: J13G (ka,ia,ja)
843  real(RP), intent(in) :: MAPF ( ia,ja,2)
844  real(RP), intent(in) :: CDZ (ka)
845  integer, intent(in) :: IIS, IIE, JJS, JJE
846 
847  real(RP) :: vel
848  integer :: k, i, j
849  !---------------------------------------------------------------------------
850 
851  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
852  do j = jjs, jje
853  do i = iis, iie
854  do k = ks+1, ke-2
855  vel = ( f2h(k,1,i_uyz) &
856  * mom(k+1,i,j) &
857  + f2h(k,2,i_uyz) &
858  * mom(k,i,j) ) &
859  / ( f2h(k,1,i_xyz) &
860  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
861  + f2h(k,2,i_xyz) &
862  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
863  flux(k,i,j) = j13g(k,i,j) / mapf(i,j,+2) * vel &
864  * ( ( val(k,i,j) &
865  + 0.5_rp * phi(val(k+1,i,j),val(k,i,j),val(k-1,i,j)) * ( val(k,i,j)-val(k-1,i,j) ) ) &
866  * ( 0.5_rp + sign(0.5_rp,vel) ) &
867  + ( val(k+1,i,j) &
868  + 0.5_rp * phi(val(k,i,j),val(k+1,i,j),val(k+2,i,j)) * ( val(k+1,i,j)-val(k+2,i,j) ) ) &
869  * ( 0.5_rp - sign(0.5_rp,vel) ) )
870  enddo
871  enddo
872  enddo
873 
874  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
875  do j = jjs, jje
876  do i = iis, iie
877  flux(ks-1,i,j) = 0.0_rp
878 
879  vel = ( f2h(ks,1,i_uyz) &
880  * mom(ks+1,i,j) &
881  + f2h(ks,2,i_uyz) &
882  * mom(ks,i,j) ) &
883  / ( f2h(ks,1,i_xyz) &
884  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
885  + f2h(ks,2,i_xyz) &
886  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
887  flux(ks,i,j) = j13g(ks,i,j) / mapf(i,j,+2) * vel &
888  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
889  vel = ( f2h(ke-1,1,i_uyz) &
890  * mom(ke,i,j) &
891  + f2h(ke-1,2,i_uyz) &
892  * mom(ke-1,i,j) ) &
893  / ( f2h(ke-1,1,i_xyz) &
894  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
895  + f2h(ke-1,2,i_xyz) &
896  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
897  flux(ke-1,i,j) = j13g(ke-1,i,j) / mapf(i,j,+2) * vel &
898  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
899 
900  flux(ke ,i,j) = 0.0_rp
901  enddo
902  enddo
903 
904  return
906 
907  !-----------------------------------------------------------------------------
910  flux, &
911  mom, val, DENS, &
912  GSQRT, J23G, MAPF, &
913  CDZ, &
914  IIS, IIE, JJS, JJE )
915  implicit none
916 
917  real(RP), intent(out) :: flux (ka,ia,ja)
918  real(RP), intent(in) :: mom (ka,ia,ja)
919  real(RP), intent(in) :: val (ka,ia,ja)
920  real(RP), intent(in) :: DENS (ka,ia,ja)
921  real(RP), intent(in) :: GSQRT (ka,ia,ja)
922  real(RP), intent(in) :: J23G (ka,ia,ja)
923  real(RP), intent(in) :: MAPF ( ia,ja,2)
924  real(RP), intent(in) :: CDZ (ka)
925  integer, intent(in) :: IIS, IIE, JJS, JJE
926 
927  real(RP) :: vel
928  integer :: k, i, j
929  !---------------------------------------------------------------------------
930 
931  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
932  do j = jjs, jje
933  do i = iis, iie
934  do k = ks+1, ke-2
935  vel = ( f2h(k,1,i_xvz) &
936  * 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) ) &
937  + f2h(k,2,i_xvz) &
938  * 0.25_rp * ( mom(k,i,j)+mom(k,i+1,j)+mom(k,i,j-1)+mom(k,i+1,j-1) ) ) &
939  / ( f2h(k,1,i_xyz) &
940  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i+1,j) ) &
941  + f2h(k,2,i_xyz) &
942  * 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) ) )
943  flux(k,i,j) = j23g(k,i,j) / mapf(i,j,+1) * vel &
944  * ( ( val(k,i,j) &
945  + 0.5_rp * phi(val(k+1,i,j),val(k,i,j),val(k-1,i,j)) * ( val(k,i,j)-val(k-1,i,j) ) ) &
946  * ( 0.5_rp + sign(0.5_rp,vel) ) &
947  + ( val(k+1,i,j) &
948  + 0.5_rp * phi(val(k,i,j),val(k+1,i,j),val(k+2,i,j)) * ( val(k+1,i,j)-val(k+2,i,j) ) ) &
949  * ( 0.5_rp - sign(0.5_rp,vel) ) )
950  enddo
951  enddo
952  enddo
953 
954  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
955  do j = jjs, jje
956  do i = iis, iie
957  flux(ks-1,i,j) = 0.0_rp
958 
959  vel = ( f2h(ks,1,i_xvz) &
960  * 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) ) &
961  + f2h(ks,2,i_xvz) &
962  * 0.25_rp * ( mom(ks,i,j)+mom(ks,i+1,j)+mom(ks,i,j-1)+mom(ks,i+1,j-1) ) ) &
963  / ( f2h(ks,1,i_xyz) &
964  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i+1,j) ) &
965  + f2h(ks,2,i_xyz) &
966  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) ) )
967  flux(ks,i,j) = j23g(ks,i,j) / mapf(i,j,+1) * vel &
968  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
969  vel = ( f2h(ke-1,1,i_xvz) &
970  * 0.25_rp * ( mom(ke,i,j)+mom(ke,i+1,j)+mom(ke,i,j-1)+mom(ke,i+1,j-1) ) &
971  + f2h(ke-1,2,i_xvz) &
972  * 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) ) ) &
973  / ( f2h(ke-1,1,i_xyz) &
974  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
975  + f2h(ke-1,2,i_xyz) &
976  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i+1,j) ) )
977  flux(ke-1,i,j) = j23g(ke-1,i,j) / mapf(i,j,+1) * vel &
978  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
979 
980  flux(ke ,i,j) = 0.0_rp
981  enddo
982  enddo
983 
984  return
986 
987  !-----------------------------------------------------------------------------
990  flux, &
991  mom, val, DENS, &
992  GSQRT, MAPF, &
993 
994  num_diff, &
995 
996  CDZ, &
997  IIS, IIE, JJS, JJE )
998  implicit none
999 
1000  real(RP), intent(out) :: flux (ka,ia,ja)
1001  real(RP), intent(in) :: mom (ka,ia,ja)
1002  real(RP), intent(in) :: val (ka,ia,ja)
1003  real(RP), intent(in) :: DENS (ka,ia,ja)
1004  real(RP), intent(in) :: GSQRT (ka,ia,ja)
1005  real(RP), intent(in) :: MAPF ( ia,ja,2)
1006 
1007  real(RP), intent(in) :: num_diff(ka,ia,ja)
1008 
1009  real(RP), intent(in) :: CDZ (ka)
1010  integer, intent(in) :: IIS, IIE, JJS, JJE
1011 
1012  real(RP) :: vel
1013  integer :: k, i, j
1014  !---------------------------------------------------------------------------
1015 
1016  ! note that x-index is added by -1
1017 
1018  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1019  do j = jjs, jje
1020  do i = iis, iie+1
1021  do k = ks, ke
1022 #ifdef DEBUG
1023  call check( __line__, mom(k,i ,j) )
1024  call check( __line__, mom(k,i-1,j) )
1025 
1026  call check( __line__, val(k,i-1,j) )
1027  call check( __line__, val(k,i,j) )
1028 
1029  call check( __line__, val(k,i-2,j) )
1030  call check( __line__, val(k,i+1,j) )
1031 
1032 #endif
1033  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i-1,j) ) ) &
1034  / ( dens(k,i,j) )
1035  flux(k,i-1,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
1036  * ( ( val(k,i-1,j) &
1037  + 0.5_rp * phi(val(k,i,j),val(k,i-1,j),val(k,i-2,j)) * ( val(k,i-1,j)-val(k,i-2,j) ) ) &
1038  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1039  + ( val(k,i,j) &
1040  + 0.5_rp * phi(val(k,i-1,j),val(k,i,j),val(k,i+1,j)) * ( val(k,i,j)-val(k,i+1,j) ) ) &
1041  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1042  + gsqrt(k,i,j) * num_diff(k,i,j)
1043  enddo
1044  enddo
1045  enddo
1046 #ifdef DEBUG
1047  k = iundef; i = iundef; j = iundef
1048 #endif
1049 
1050  return
1052 
1053  !-----------------------------------------------------------------------------
1056  flux, &
1057  mom, val, DENS, &
1058  GSQRT, MAPF, &
1059 
1060  num_diff, &
1061 
1062  CDZ, &
1063  IIS, IIE, JJS, JJE )
1064  implicit none
1065 
1066  real(RP), intent(out) :: flux (ka,ia,ja)
1067  real(RP), intent(in) :: mom (ka,ia,ja)
1068  real(RP), intent(in) :: val (ka,ia,ja)
1069  real(RP), intent(in) :: DENS (ka,ia,ja)
1070  real(RP), intent(in) :: GSQRT (ka,ia,ja)
1071  real(RP), intent(in) :: MAPF ( ia,ja,2)
1072 
1073  real(RP), intent(in) :: num_diff(ka,ia,ja)
1074 
1075  real(RP), intent(in) :: CDZ (ka)
1076  integer, intent(in) :: IIS, IIE, JJS, JJE
1077 
1078  real(RP) :: vel
1079  integer :: k, i, j
1080  !---------------------------------------------------------------------------
1081 
1082  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1083  do j = jjs-1, jje
1084  do i = iis, iie
1085  do k = ks, ke
1086 #ifdef DEBUG
1087  call check( __line__, mom(k,i ,j) )
1088  call check( __line__, mom(k,i-1,j) )
1089 
1090  call check( __line__, val(k,i,j) )
1091  call check( __line__, val(k,i,j+1) )
1092 
1093  call check( __line__, val(k,i,j-1) )
1094  call check( __line__, val(k,i,j+2) )
1095 
1096 #endif
1097  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i+1,j) ) ) &
1098  / ( 0.25_rp * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) )
1099  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
1100  * ( ( val(k,i,j) &
1101  + 0.5_rp * phi(val(k,i,j+1),val(k,i,j),val(k,i,j-1)) * ( val(k,i,j)-val(k,i,j-1) ) ) &
1102  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1103  + ( val(k,i,j+1) &
1104  + 0.5_rp * phi(val(k,i,j),val(k,i,j+1),val(k,i,j+2)) * ( val(k,i,j+1)-val(k,i,j+2) ) ) &
1105  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1106  + gsqrt(k,i,j) * num_diff(k,i,j)
1107  enddo
1108  enddo
1109  enddo
1110 #ifdef DEBUG
1111  k = iundef; i = iundef; j = iundef
1112 #endif
1113 
1114  return
1116 
1117 
1118 
1119  !-----------------------------------------------------------------------------
1122  flux, &
1123  mom, val, DENS, &
1124  GSQRT, J33G, &
1125 
1126  num_diff, &
1127 
1128  CDZ, &
1129  IIS, IIE, JJS, JJE )
1130  implicit none
1131 
1132  real(RP), intent(out) :: flux (ka,ia,ja)
1133  real(RP), intent(in) :: mom (ka,ia,ja)
1134  real(RP), intent(in) :: val (ka,ia,ja)
1135  real(RP), intent(in) :: DENS (ka,ia,ja)
1136  real(RP), intent(in) :: GSQRT (ka,ia,ja)
1137  real(RP), intent(in) :: J33G
1138 
1139  real(RP), intent(in) :: num_diff(ka,ia,ja)
1140 
1141  real(RP), intent(in) :: CDZ (ka)
1142  integer, intent(in) :: IIS, IIE, JJS, JJE
1143 
1144  real(RP) :: vel
1145  integer :: k, i, j
1146  !---------------------------------------------------------------------------
1147 
1148  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1149  do j = jjs, jje
1150  do i = iis, iie
1151  do k = ks+1, ke-2
1152 #ifdef DEBUG
1153  call check( __line__, mom(k,i,j) )
1154  call check( __line__, mom(k,i,j+1) )
1155 
1156  call check( __line__, val(k,i,j) )
1157  call check( __line__, val(k+1,i,j) )
1158 
1159  call check( __line__, val(k-1,i,j) )
1160  call check( __line__, val(k+2,i,j) )
1161 
1162 #endif
1163  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j+1) ) ) &
1164  / ( f2h(k,1,i_xyz) &
1165  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
1166  + f2h(k,2,i_xyz) &
1167  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
1168  flux(k,i,j) = j33g * vel &
1169  * ( ( val(k,i,j) &
1170  + 0.5_rp * phi(val(k+1,i,j),val(k,i,j),val(k-1,i,j)) * ( val(k,i,j)-val(k-1,i,j) ) ) &
1171  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1172  + ( val(k+1,i,j) &
1173  + 0.5_rp * phi(val(k,i,j),val(k+1,i,j),val(k+2,i,j)) * ( val(k+1,i,j)-val(k+2,i,j) ) ) &
1174  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1175  + gsqrt(k,i,j) * num_diff(k,i,j)
1176  enddo
1177  enddo
1178  enddo
1179 #ifdef DEBUG
1180  k = iundef; i = iundef; j = iundef
1181 #endif
1182 
1183  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
1184  do j = jjs, jje
1185  do i = iis, iie
1186 #ifdef DEBUG
1187 
1188  call check( __line__, mom(ks,i ,j) )
1189  call check( __line__, mom(ks,i,j+1) )
1190  call check( __line__, val(ks+1,i,j) )
1191  call check( __line__, val(ks,i,j) )
1192 
1193 #endif
1194  flux(ks-1,i,j) = 0.0_rp
1195 
1196  vel = ( 0.5_rp * ( mom(ks,i,j)+mom(ks,i,j+1) ) ) &
1197  / ( f2h(ks,1,i_xyz) &
1198  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
1199  + f2h(ks,2,i_xyz) &
1200  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
1201  flux(ks,i,j) = j33g * vel &
1202  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) ) &
1203  + gsqrt(ks,i,j) * num_diff(ks,i,j)
1204  vel = ( 0.5_rp * ( mom(ke-1,i,j)+mom(ke-1,i,j+1) ) ) &
1205  / ( f2h(ke-1,1,i_xyz) &
1206  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
1207  + f2h(ke-1,2,i_xyz) &
1208  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
1209  flux(ke-1,i,j) = j33g * vel &
1210  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) ) &
1211  + gsqrt(ke-1,i,j) * num_diff(ke-1,i,j)
1212 
1213  flux(ke,i,j) = 0.0_rp
1214  enddo
1215  enddo
1216 #ifdef DEBUG
1217  k = iundef; i = iundef; j = iundef
1218 #endif
1219 
1220  return
1222 
1223  !-----------------------------------------------------------------------------
1226  flux, &
1227  mom, val, DENS, &
1228  GSQRT, J13G, MAPF, &
1229  CDZ, &
1230  IIS, IIE, JJS, JJE )
1231  implicit none
1232 
1233  real(RP), intent(out) :: flux (ka,ia,ja)
1234  real(RP), intent(in) :: mom (ka,ia,ja)
1235  real(RP), intent(in) :: val (ka,ia,ja)
1236  real(RP), intent(in) :: DENS (ka,ia,ja)
1237  real(RP), intent(in) :: GSQRT (ka,ia,ja)
1238  real(RP), intent(in) :: J13G (ka,ia,ja)
1239  real(RP), intent(in) :: MAPF ( ia,ja,2)
1240  real(RP), intent(in) :: CDZ (ka)
1241  integer, intent(in) :: IIS, IIE, JJS, JJE
1242 
1243  real(RP) :: vel
1244  integer :: k, i, j
1245  !---------------------------------------------------------------------------
1246 
1247  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1248  do j = jjs, jje
1249  do i = iis, iie
1250  do k = ks+1, ke-2
1251  vel = ( f2h(k,1,i_uyz) &
1252  * 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) ) &
1253  + f2h(k,2,i_uyz) &
1254  * 0.25_rp * ( mom(k,i,j)+mom(k,i-1,j)+mom(k,i,j+1)+mom(k,i-1,j+1) ) ) &
1255  / ( f2h(k,1,i_xyz) &
1256  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
1257  + f2h(k,2,i_xyz) &
1258  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
1259  flux(k,i,j) = j13g(k,i,j) / mapf(i,j,+2) * vel &
1260  * ( ( val(k,i,j) &
1261  + 0.5_rp * phi(val(k+1,i,j),val(k,i,j),val(k-1,i,j)) * ( val(k,i,j)-val(k-1,i,j) ) ) &
1262  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1263  + ( val(k+1,i,j) &
1264  + 0.5_rp * phi(val(k,i,j),val(k+1,i,j),val(k+2,i,j)) * ( val(k+1,i,j)-val(k+2,i,j) ) ) &
1265  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1266  enddo
1267  enddo
1268  enddo
1269 
1270  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
1271  do j = jjs, jje
1272  do i = iis, iie
1273  flux(ks-1,i,j) = 0.0_rp
1274 
1275  vel = ( f2h(ks,1,i_uyz) &
1276  * 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) ) &
1277  + f2h(ks,2,i_uyz) &
1278  * 0.25_rp * ( mom(ks,i,j)+mom(ks,i-1,j)+mom(ks,i,j+1)+mom(ks,i-1,j+1) ) ) &
1279  / ( f2h(ks,1,i_xyz) &
1280  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
1281  + f2h(ks,2,i_xyz) &
1282  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
1283  flux(ks,i,j) = j13g(ks,i,j) / mapf(i,j,+2) * vel &
1284  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
1285  vel = ( f2h(ke-1,1,i_uyz) &
1286  * 0.25_rp * ( mom(ke,i,j)+mom(ke,i-1,j)+mom(ke,i,j+1)+mom(ke,i-1,j+1) ) &
1287  + f2h(ke-1,2,i_uyz) &
1288  * 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) ) ) &
1289  / ( f2h(ke-1,1,i_xyz) &
1290  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
1291  + f2h(ke-1,2,i_xyz) &
1292  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
1293  flux(ke-1,i,j) = j13g(ke-1,i,j) / mapf(i,j,+2) * vel &
1294  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
1295 
1296  flux(ke ,i,j) = 0.0_rp
1297  enddo
1298  enddo
1299 
1300  return
1302 
1303  !-----------------------------------------------------------------------------
1306  flux, &
1307  mom, val, DENS, &
1308  GSQRT, J23G, MAPF, &
1309  CDZ, &
1310  IIS, IIE, JJS, JJE )
1311  implicit none
1312 
1313  real(RP), intent(out) :: flux (ka,ia,ja)
1314  real(RP), intent(in) :: mom (ka,ia,ja)
1315  real(RP), intent(in) :: val (ka,ia,ja)
1316  real(RP), intent(in) :: DENS (ka,ia,ja)
1317  real(RP), intent(in) :: GSQRT (ka,ia,ja)
1318  real(RP), intent(in) :: J23G (ka,ia,ja)
1319  real(RP), intent(in) :: MAPF ( ia,ja,2)
1320  real(RP), intent(in) :: CDZ (ka)
1321  integer, intent(in) :: IIS, IIE, JJS, JJE
1322 
1323  real(RP) :: vel
1324  integer :: k, i, j
1325  !---------------------------------------------------------------------------
1326 
1327  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1328  do j = jjs, jje
1329  do i = iis, iie
1330  do k = ks+1, ke-2
1331  vel = ( f2h(k,1,i_xvz) &
1332  * mom(k+1,i,j) &
1333  + f2h(k,2,i_xvz) &
1334  * mom(k,i,j) ) &
1335  / ( f2h(k,1,i_xyz) &
1336  * 0.5_rp * ( dens(k+1,i,j)+dens(k+1,i,j+1) ) &
1337  + f2h(k,2,i_xyz) &
1338  * 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) ) )
1339  flux(k,i,j) = j23g(k,i,j) / mapf(i,j,+1) * vel &
1340  * ( ( val(k,i,j) &
1341  + 0.5_rp * phi(val(k+1,i,j),val(k,i,j),val(k-1,i,j)) * ( val(k,i,j)-val(k-1,i,j) ) ) &
1342  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1343  + ( val(k+1,i,j) &
1344  + 0.5_rp * phi(val(k,i,j),val(k+1,i,j),val(k+2,i,j)) * ( val(k+1,i,j)-val(k+2,i,j) ) ) &
1345  * ( 0.5_rp - sign(0.5_rp,vel) ) )
1346  enddo
1347  enddo
1348  enddo
1349 
1350  !$omp parallel do private(i,j) OMP_SCHEDULE_ collapse(2)
1351  do j = jjs, jje
1352  do i = iis, iie
1353  flux(ks-1,i,j) = 0.0_rp
1354 
1355  vel = ( f2h(ks,1,i_xvz) &
1356  * mom(ks+1,i,j) &
1357  + f2h(ks,2,i_xvz) &
1358  * mom(ks,i,j) ) &
1359  / ( f2h(ks,1,i_xyz) &
1360  * 0.5_rp * ( dens(ks+1,i,j)+dens(ks+1,i,j+1) ) &
1361  + f2h(ks,2,i_xyz) &
1362  * 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) ) )
1363  flux(ks,i,j) = j23g(ks,i,j) / mapf(i,j,+1) * vel &
1364  * ( f2 * ( val(ks+1,i,j)+val(ks,i,j) ) )
1365  vel = ( f2h(ke-1,1,i_xvz) &
1366  * mom(ke,i,j) &
1367  + f2h(ke-1,2,i_xvz) &
1368  * mom(ke-1,i,j) ) &
1369  / ( f2h(ke-1,1,i_xyz) &
1370  * 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
1371  + f2h(ke-1,2,i_xyz) &
1372  * 0.5_rp * ( dens(ke-1,i,j)+dens(ke-1,i,j+1) ) )
1373  flux(ke-1,i,j) = j23g(ke-1,i,j) / mapf(i,j,+1) * vel &
1374  * ( f2 * ( val(ke,i,j)+val(ke-1,i,j) ) )
1375 
1376  flux(ke ,i,j) = 0.0_rp
1377  enddo
1378  enddo
1379 
1380  return
1382 
1383  !-----------------------------------------------------------------------------
1386  flux, &
1387  mom, val, DENS, &
1388  GSQRT, MAPF, &
1389 
1390  num_diff, &
1391 
1392  CDZ, &
1393  IIS, IIE, JJS, JJE )
1394  implicit none
1395 
1396  real(RP), intent(out) :: flux (ka,ia,ja)
1397  real(RP), intent(in) :: mom (ka,ia,ja)
1398  real(RP), intent(in) :: val (ka,ia,ja)
1399  real(RP), intent(in) :: DENS (ka,ia,ja)
1400  real(RP), intent(in) :: GSQRT (ka,ia,ja)
1401  real(RP), intent(in) :: MAPF ( ia,ja,2)
1402 
1403  real(RP), intent(in) :: num_diff(ka,ia,ja)
1404 
1405  real(RP), intent(in) :: CDZ (ka)
1406  integer, intent(in) :: IIS, IIE, JJS, JJE
1407 
1408  real(RP) :: vel
1409  integer :: k, i, j
1410  !---------------------------------------------------------------------------
1411 
1412  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1413  do j = jjs, jje
1414  do i = iis-1, iie
1415  do k = ks, ke
1416 #ifdef DEBUG
1417  call check( __line__, mom(k,i ,j) )
1418  call check( __line__, mom(k,i,j-1) )
1419 
1420  call check( __line__, val(k,i,j) )
1421  call check( __line__, val(k,i+1,j) )
1422 
1423  call check( __line__, val(k,i-1,j) )
1424  call check( __line__, val(k,i+2,j) )
1425 
1426 #endif
1427  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j+1) ) ) &
1428  / ( 0.25_rp * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) )
1429  flux(k,i,j) = gsqrt(k,i,j) / mapf(i,j,+2) * vel &
1430  * ( ( val(k,i,j) &
1431  + 0.5_rp * phi(val(k,i+1,j),val(k,i,j),val(k,i-1,j)) * ( val(k,i,j)-val(k,i-1,j) ) ) &
1432  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1433  + ( val(k,i+1,j) &
1434  + 0.5_rp * phi(val(k,i,j),val(k,i+1,j),val(k,i+2,j)) * ( val(k,i+1,j)-val(k,i+2,j) ) ) &
1435  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1436  + gsqrt(k,i,j) * num_diff(k,i,j)
1437  enddo
1438  enddo
1439  enddo
1440 #ifdef DEBUG
1441  k = iundef; i = iundef; j = iundef
1442 #endif
1443 
1444  return
1446 
1447  !-----------------------------------------------------------------------------
1450  flux, &
1451  mom, val, DENS, &
1452  GSQRT, MAPF, &
1453 
1454  num_diff, &
1455 
1456  CDZ, &
1457  IIS, IIE, JJS, JJE )
1458  implicit none
1459 
1460  real(RP), intent(out) :: flux (ka,ia,ja)
1461  real(RP), intent(in) :: mom (ka,ia,ja)
1462  real(RP), intent(in) :: val (ka,ia,ja)
1463  real(RP), intent(in) :: DENS (ka,ia,ja)
1464  real(RP), intent(in) :: GSQRT (ka,ia,ja)
1465  real(RP), intent(in) :: MAPF ( ia,ja,2)
1466 
1467  real(RP), intent(in) :: num_diff(ka,ia,ja)
1468 
1469  real(RP), intent(in) :: CDZ (ka)
1470  integer, intent(in) :: IIS, IIE, JJS, JJE
1471 
1472  real(RP) :: vel
1473  integer :: k, i, j
1474  !---------------------------------------------------------------------------
1475 
1476  ! note that y-index is added by -1
1477 
1478  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
1479  do j = jjs, jje+1
1480  do i = iis, iie
1481  do k = ks, ke
1482 #ifdef DEBUG
1483  call check( __line__, mom(k,i ,j) )
1484  call check( __line__, mom(k,i,j-1) )
1485 
1486  call check( __line__, val(k,i,j-1) )
1487  call check( __line__, val(k,i,j) )
1488 
1489  call check( __line__, val(k,i,j-2) )
1490  call check( __line__, val(k,i,j+1) )
1491 
1492 #endif
1493  vel = ( 0.5_rp * ( mom(k,i,j)+mom(k,i,j-1) ) ) &
1494  / ( dens(k,i,j) )
1495  flux(k,i,j-1) = gsqrt(k,i,j) / mapf(i,j,+1) * vel &
1496  * ( ( val(k,i,j-1) &
1497  + 0.5_rp * phi(val(k,i,j),val(k,i,j-1),val(k,i,j-2)) * ( val(k,i,j-1)-val(k,i,j-2) ) ) &
1498  * ( 0.5_rp + sign(0.5_rp,vel) ) &
1499  + ( val(k,i,j) &
1500  + 0.5_rp * phi(val(k,i,j-1),val(k,i,j),val(k,i,j+1)) * ( val(k,i,j)-val(k,i,j+1) ) ) &
1501  * ( 0.5_rp - sign(0.5_rp,vel) ) ) &
1502  + gsqrt(k,i,j) * num_diff(k,i,j)
1503  enddo
1504  enddo
1505  enddo
1506 #ifdef DEBUG
1507  k = iundef; i = iundef; j = iundef
1508 #endif
1509 
1510  return
1512 
1513 
1514 
1515 
1516 
1517 
1518  !-----------------------------------------------------------------------------
1519  function phi(v1, v2, v3)
1520  use scale_const, only: &
1521  eps => const_eps
1522  implicit none
1523 
1524  real(RP) :: phi
1525  real(RP), intent(in) :: v1
1526  real(RP), intent(in) :: v2
1527  real(RP), intent(in) :: v3
1528 
1529  real(RP) :: r2
1530  real(RP) :: zerosw1, zerosw2
1531  !---------------------------------------------------------------------------
1532 
1533  zerosw1 = eps - sign(eps, abs(v1-v2)-eps)
1534  zerosw2 = eps - sign(eps, abs(v2-v3)-eps)
1535  r2 = 2.0_rp * (v1-v2+zerosw1*zerosw2) / (v2-v3+zerosw2)
1536 
1537  phi = max(0.0_rp, min(r2, min((1.0_rp+r2)/3.0_rp, 2.0_rp) ) )
1538 
1539  end function phi
1540 
1541 
1543 
1544 !--
1545 ! vi:set readonly sw=4 ts=8
1546 !
1547 !Local Variables:
1548 !mode: f90
1549 !buffer-read-only: t
1550 !End:
1551 !
1552 !++
subroutine, public atmos_dyn_fvm_fluxx_xvz_ud3koren1993(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation X-flux at XV
subroutine, public atmos_dyn_fvm_flux_valuew_z_ud3koren1993(valW, mflx, val, GSQRT, CDZ)
value at XYW
module DEBUG
Definition: scale_debug.F90:13
subroutine, public atmos_dyn_fvm_fluxy_xyz_ud3koren1993(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation Y-flux at XYZ
subroutine, public atmos_dyn_fvm_fluxz_uyz_ud3koren1993(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation z-flux at UY
subroutine, public atmos_dyn_fvm_fluxy_xvz_ud3koren1993(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation Y-flux at XV
module STDIO
Definition: scale_stdio.F90:12
integer, public ke
end point of inner domain: z, local
subroutine, public atmos_dyn_fvm_fluxx_xyz_ud3koren1993(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation X-flux at XYZ
subroutine, public check(current_line, v)
Undefined value checker.
Definition: scale_debug.F90:58
real(rp), public const_undef
Definition: scale_const.F90:43
module grid index
module TRACER
module Index
Definition: scale_index.F90:14
integer, public ia
of x whole cells (local, with HALO)
subroutine, public atmos_dyn_fvm_fluxj13_uyz_ud3koren1993(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, IIS, IIE, JJS, JJE)
calculation J13-flux at UYZ
integer, public ka
of z whole cells (local, with HALO)
subroutine, public atmos_dyn_fvm_fluxz_xyz_ud3koren1993(flux, mflx, val, GSQRT, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation z-flux at XYZ
subroutine, public atmos_dyn_fvm_fluxj13_xyw_ud3koren1993(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, IIS, IIE, JJS, JJE)
calculation J13-flux at XYW
integer, parameter, public const_undef2
undefined value (INT2)
Definition: scale_const.F90:40
subroutine, public atmos_dyn_fvm_fluxj23_uyz_ud3koren1993(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, IIS, IIE, JJS, JJE)
calculation J23-flux at UYZ
subroutine, public atmos_dyn_fvm_fluxj23_xyw_ud3koren1993(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, IIS, IIE, JJS, JJE)
calculation J23-flux at XYW
module PROCESS
subroutine, public atmos_dyn_fvm_fluxx_uyz_ud3koren1993(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation X-flux at UY
subroutine, public atmos_dyn_fvm_fluxj23_xvz_ud3koren1993(flux, mom, val, DENS, GSQRT, J23G, MAPF, CDZ, IIS, IIE, JJS, JJE)
calculation J23-flux at XVZ
module CONSTANT
Definition: scale_const.F90:14
integer, public ks
start point of inner domain: z, local
subroutine, public atmos_dyn_fvm_fluxj13_xvz_ud3koren1993(flux, mom, val, DENS, GSQRT, J13G, MAPF, CDZ, IIS, IIE, JJS, JJE)
calculation J13-flux at XVZ
module profiler
Definition: scale_prof.F90:10
subroutine, public atmos_dyn_fvm_fluxy_uyz_ud3koren1993(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation Y-flux at UY
real(rp), public const_eps
small number
Definition: scale_const.F90:36
module scale_atmos_dyn_fvm_flux_ud3Koren1993
module PRECISION
subroutine, public atmos_dyn_fvm_fluxz_xyw_ud3koren1993(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, FDZ, dtrk, IIS, IIE, JJS, JJE)
calculation z-flux at XYW
subroutine, public atmos_dyn_fvm_fluxx_xyw_ud3koren1993(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation X-flux at XYW
subroutine, public atmos_dyn_fvm_fluxz_xvz_ud3koren1993(flux, mom, val, DENS, GSQRT, J33G, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation z-flux at XV
subroutine, public atmos_dyn_fvm_fluxy_xyw_ud3koren1993(flux, mom, val, DENS, GSQRT, MAPF, num_diff, CDZ, IIS, IIE, JJS, JJE)
calculation Y-flux at XYW
integer, public ja
of y whole cells (local, with HALO)