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