SCALE-RM
Functions/Subroutines
scale_atmos_dyn_common Module Reference

module Atmosphere / Dynamics common More...

Functions/Subroutines

subroutine, public atmos_dyn_filter_setup (num_diff, num_diff_q, CDZ, CDX, CDY, FDZ, FDX, FDY)
 Setup. More...
 
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. More...
 
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. More...
 
subroutine, public atmos_dyn_filter_tend (phi_t, phi, rdz, rdx, rdy, KO, IO, JO)
 
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)
 
subroutine, public atmos_dyn_divergence (DDIV, MOMZ, MOMX, MOMY, GSQRT, J13G, J23G, J33G, MAPF, RCDZ, RCDX, RCDY, RFDZ, FDZ)
 
subroutine calc_numdiff (work, iwork, data, nd_order, KO, IO, JO, KEE)
 
subroutine calc_diff3 (diff, phi, KO, IO, JO)
 
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. More...
 

Detailed Description

module Atmosphere / Dynamics common

Description
common subroutines for Atmospheric dynamical process
Author
Team SCALE
History
  • 2013-06-18 (S.Nishizawa) [new] splited from dynamical core

Function/Subroutine Documentation

◆ atmos_dyn_filter_setup()

subroutine, public scale_atmos_dyn_common::atmos_dyn_filter_setup ( real(rp), dimension(ka,ia,ja,5,3), intent(inout)  num_diff,
real(rp), dimension(ka,ia,ja,3), intent(inout)  num_diff_q,
real(rp), dimension(ka), intent(in)  CDZ,
real(rp), dimension(ia), intent(in)  CDX,
real(rp), dimension(ja), intent(in)  CDY,
real(rp), dimension(ka-1), intent(in)  FDZ,
real(rp), dimension(ia-1), intent(in)  FDX,
real(rp), dimension(ja-1), intent(in)  FDY 
)

Setup.

Definition at line 96 of file scale_atmos_dyn_common.F90.

References scale_comm::comm_vars8_init(), scale_index::i_dens, scale_index::i_momx, scale_index::i_momy, scale_index::i_momz, scale_index::i_rhot, scale_grid_index::ia, scale_grid_index::ie, scale_grid_index::ihalo, scale_grid_index::is, scale_grid_index::ja, scale_grid_index::je, scale_grid_index::jhalo, scale_grid_index::js, scale_grid_index::ka, scale_grid_index::ke, scale_grid_index::khalo, scale_grid_index::ks, scale_process::prc_mpistop(), scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by scale_atmos_dyn::atmos_dyn_setup().

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
scale_process::prc_mpistop
subroutine, public prc_mpistop
Abort MPI.
Definition: scale_process.F90:234
scale_comm::comm_vars8_init
subroutine, public comm_vars8_init(var, vid)
Register variables.
Definition: scale_comm.F90:328
scale_comm
module COMMUNICATION
Definition: scale_comm.F90:23
scale_process
module PROCESS
Definition: scale_process.F90:14
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_numfilter_coef()

subroutine, public scale_atmos_dyn_common::atmos_dyn_numfilter_coef ( real(rp), dimension(ka,ia,ja,5,3), intent(out)  num_diff,
real(rp), dimension(ka,ia,ja), intent(in)  DENS,
real(rp), dimension(ka,ia,ja), intent(in)  MOMZ,
real(rp), dimension(ka,ia,ja), intent(in)  MOMX,
real(rp), dimension(ka,ia,ja), intent(in)  MOMY,
real(rp), dimension(ka,ia,ja), intent(in)  RHOT,
real(rp), dimension(ka), intent(in)  CDZ,
real(rp), dimension(ia), intent(in)  CDX,
real(rp), dimension(ja), intent(in)  CDY,
real(rp), dimension(ka-1), intent(in)  FDZ,
real(rp), dimension(ia-1), intent(in)  FDX,
real(rp), dimension(ja-1), intent(in)  FDY,
real(rp), intent(in)  DT,
real(rp), dimension(ka,ia,ja), intent(in)  REF_dens,
real(rp), dimension(ka,ia,ja), intent(in)  REF_pott,
real(rp), intent(in)  ND_COEF,
integer, intent(in)  ND_ORDER,
real(rp), intent(in)  ND_SFC_FACT,
logical, intent(in)  ND_USE_RS 
)

Calc coefficient of numerical filter.

Definition at line 330 of file scale_atmos_dyn_common.F90.

References calc_numdiff(), scale_index::i_dens, scale_index::i_momx, scale_index::i_momy, scale_index::i_momz, scale_index::i_rhot, scale_grid_index::ie, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ka, scale_grid_index::ke, scale_grid_index::ks, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by scale_atmos_dyn_tinteg_large_euler::atmos_dyn_tinteg_large_euler(), scale_atmos_dyn_tinteg_large_rk3::atmos_dyn_tinteg_large_rk3(), and scale_atmos_dyn_tstep_large_fvm_heve::atmos_dyn_tstep_large_fvm_heve().

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
scale_comm
module COMMUNICATION
Definition: scale_comm.F90:23
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_numfilter_coef_q()

subroutine, public scale_atmos_dyn_common::atmos_dyn_numfilter_coef_q ( real(rp), dimension(ka,ia,ja,3), intent(out)  num_diff_q,
real(rp), dimension(ka,ia,ja), intent(in)  DENS,
real(rp), dimension(ka,ia,ja), intent(in)  QTRC,
real(rp), dimension(ka), intent(in)  CDZ,
real(rp), dimension(ia), intent(in)  CDX,
real(rp), dimension(ja), intent(in)  CDY,
real(rp), intent(in)  dt,
real(rp), dimension(ka,ia,ja), intent(in)  REF_qv,
integer, intent(in)  iq,
real(rp), intent(in)  ND_COEF,
integer, intent(in)  ND_ORDER,
real(rp), intent(in)  ND_SFC_FACT,
logical, intent(in)  ND_USE_RS 
)

Calc coefficient of numerical filter.

Definition at line 931 of file scale_atmos_dyn_common.F90.

References calc_numdiff(), scale_tracer::i_qv, scale_grid_index::ie, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ka, scale_grid_index::ke, scale_grid_index::ks, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by scale_atmos_dyn_tinteg_large_euler::atmos_dyn_tinteg_large_euler(), scale_atmos_dyn_tinteg_large_rk3::atmos_dyn_tinteg_large_rk3(), and scale_atmos_dyn_tstep_large_fvm_heve::atmos_dyn_tstep_large_fvm_heve().

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
scale_comm
module COMMUNICATION
Definition: scale_comm.F90:23
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_filter_tend()

subroutine, public scale_atmos_dyn_common::atmos_dyn_filter_tend ( real(rp), dimension(ka,ia,ja), intent(out)  phi_t,
real(rp), dimension (ka,ia,ja), intent(in)  phi,
real(rp), dimension(:), intent(in)  rdz,
real(rp), dimension(:), intent(in)  rdx,
real(rp), dimension(:), intent(in)  rdy,
integer, intent(in)  KO,
integer, intent(in)  IO,
integer, intent(in)  JO 
)

Definition at line 1143 of file scale_atmos_dyn_common.F90.

References calc_diff3(), scale_grid_index::ie, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ke, scale_grid_index::ks, scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

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
scale_comm
module COMMUNICATION
Definition: scale_comm.F90:23
Here is the call graph for this function:

◆ atmos_dyn_copy_boundary()

subroutine, public scale_atmos_dyn_common::atmos_dyn_copy_boundary ( real(rp), dimension (ka,ia,ja), intent(inout)  DENS,
real(rp), dimension (ka,ia,ja), intent(inout)  MOMZ,
real(rp), dimension (ka,ia,ja), intent(inout)  MOMX,
real(rp), dimension (ka,ia,ja), intent(inout)  MOMY,
real(rp), dimension (ka,ia,ja), intent(inout)  RHOT,
real(rp), dimension (ka,ia,ja,va), intent(inout)  PROG,
real(rp), dimension(ka,ia,ja), intent(in)  DENS0,
real(rp), dimension(ka,ia,ja), intent(in)  MOMZ0,
real(rp), dimension(ka,ia,ja), intent(in)  MOMX0,
real(rp), dimension(ka,ia,ja), intent(in)  MOMY0,
real(rp), dimension(ka,ia,ja), intent(in)  RHOT0,
real(rp), dimension(ka,ia,ja,va), intent(in)  PROG0,
logical, intent(in)  BND_W,
logical, intent(in)  BND_E,
logical, intent(in)  BND_S,
logical, intent(in)  BND_N 
)

Definition at line 1187 of file scale_atmos_dyn_common.F90.

References scale_grid_index::ia, scale_grid_index::ie, scale_grid_index::is, scale_grid_index::ja, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ke, scale_grid_index::ks, and scale_index::va.

Referenced by scale_atmos_dyn_tinteg_short_rk3::atmos_dyn_tinteg_short_rk3(), and scale_atmos_dyn_tinteg_short_rk4::atmos_dyn_tinteg_short_rk4().

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
Here is the caller graph for this function:

◆ atmos_dyn_divergence()

subroutine, public scale_atmos_dyn_common::atmos_dyn_divergence ( real(rp), dimension(ka,ia,ja), intent(out)  DDIV,
real(rp), dimension(ka,ia,ja), intent(in)  MOMZ,
real(rp), dimension(ka,ia,ja), intent(in)  MOMX,
real(rp), dimension(ka,ia,ja), intent(in)  MOMY,
real(rp), dimension(ka,ia,ja,7), intent(in)  GSQRT,
real(rp), dimension(ka,ia,ja,7), intent(in)  J13G,
real(rp), dimension(ka,ia,ja,7), intent(in)  J23G,
real(rp), intent(in)  J33G,
real(rp), dimension(ia,ja,2,7), intent(in)  MAPF,
real(rp), dimension(ka), intent(in)  RCDZ,
real(rp), dimension(ia), intent(in)  RCDX,
real(rp), dimension(ja), intent(in)  RCDY,
real(rp), dimension(ka-1), intent(in)  RFDZ,
real(rp), dimension(ka-1), intent(in)  FDZ 
)

Definition at line 1303 of file scale_atmos_dyn_common.F90.

References scale_gridtrans::i_uy, scale_gridtrans::i_uyz, scale_gridtrans::i_xv, scale_gridtrans::i_xvz, scale_gridtrans::i_xy, scale_gridtrans::i_xyw, scale_gridtrans::i_xyz, scale_grid_index::ie, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ke, scale_grid_index::ks, scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

Referenced by scale_atmos_dyn_tstep_large_fvm_heve::atmos_dyn_tstep_large_fvm_heve().

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
scale_gridtrans::i_xvz
integer, public i_xvz
Definition: scale_gridtrans.F90:54
scale_gridtrans::i_xy
integer, public i_xy
Definition: scale_gridtrans.F90:57
scale_gridtrans
module GRIDTRANS
Definition: scale_gridtrans.F90:15
scale_gridtrans::i_uy
integer, public i_uy
Definition: scale_gridtrans.F90:58
scale_gridtrans::i_xyw
integer, public i_xyw
Definition: scale_gridtrans.F90:50
scale_gridtrans::i_xv
integer, public i_xv
Definition: scale_gridtrans.F90:59
scale_gridtrans::i_uyz
integer, public i_uyz
Definition: scale_gridtrans.F90:53
scale_gridtrans::i_xyz
integer, public i_xyz
Definition: scale_gridtrans.F90:49
Here is the call graph for this function:
Here is the caller graph for this function:

◆ calc_numdiff()

subroutine scale_atmos_dyn_common::calc_numdiff ( real(rp), dimension(ka,ia,ja,3,2), intent(out)  work,
integer, intent(out)  iwork,
real(rp), dimension(ka,ia,ja), intent(in)  data,
integer, intent(in)  nd_order,
integer, intent(in)  KO,
integer, intent(in)  IO,
integer, intent(in)  JO,
integer, intent(in)  KEE 
)

Definition at line 1392 of file scale_atmos_dyn_common.F90.

References calc_diff3(), scale_prof::prof_rapend(), scale_prof::prof_rapstart(), scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by atmos_dyn_numfilter_coef(), and atmos_dyn_numfilter_coef_q().

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
scale_comm
module COMMUNICATION
Definition: scale_comm.F90:23
Here is the call graph for this function:
Here is the caller graph for this function:

◆ calc_diff3()

subroutine scale_atmos_dyn_common::calc_diff3 ( real(rp), dimension(ka,ia,ja,3), intent(out)  diff,
real(rp), dimension(ka,ia,ja), intent(in)  phi,
integer, intent(in)  KO,
integer, intent(in)  IO,
integer, intent(in)  JO 
)

Definition at line 1465 of file scale_atmos_dyn_common.F90.

References scale_debug::check(), scale_grid_index::ie, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ka, scale_grid_index::ke, scale_grid_index::ks, scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by atmos_dyn_filter_tend(), and calc_numdiff().

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
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_dyn_fct()

subroutine, public scale_atmos_dyn_common::atmos_dyn_fct ( real(rp), dimension(ka,ia,ja,3), intent(out)  qflx_anti,
real(rp), dimension(ka,ia,ja), intent(in)  phi_in,
real(rp), dimension(ka,ia,ja), intent(in)  DENS0,
real(rp), dimension (ka,ia,ja), intent(in)  DENS,
real(rp), dimension(ka,ia,ja,3), intent(in)  qflx_hi,
real(rp), dimension(ka,ia,ja,3), intent(in)  qflx_lo,
real(rp), dimension(ka,ia,ja,3), intent(in)  mflx_hi,
real(rp), dimension(:), intent(in)  rdz,
real(rp), dimension(:), intent(in)  rdx,
real(rp), dimension(:), intent(in)  rdy,
real(rp), dimension(ka,ia,ja), intent(in)  GSQRT,
real(rp), dimension(ia,ja,2), intent(in)  MAPF,
real(rp), intent(in)  dt,
logical, intent(in)  flag_vect 
)

Flux Correction Transport Limiter.

Parameters
[in]gsqrtvertical metrics {G}^1/2
[in]mapfmap factor

Definition at line 1794 of file scale_atmos_dyn_common.F90.

References scale_debug::check(), scale_const::const_eps, scale_const::const_undef, scale_const::const_undef2, scale_grid_index::iblock, scale_grid_index::ie, scale_grid_index::is, scale_grid_index::jblock, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ke, scale_grid_index::ks, scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by scale_atmos_dyn_tinteg_large_euler::atmos_dyn_tinteg_large_euler(), scale_atmos_dyn_tinteg_large_rk3::atmos_dyn_tinteg_large_rk3(), scale_atmos_dyn_tstep_large_fvm_heve::atmos_dyn_tstep_large_fvm_heve(), scale_atmos_dyn_tstep_short_fvm_heve::atmos_dyn_tstep_short_fvm_heve(), scale_atmos_dyn_tstep_short_fvm_hevi::atmos_dyn_tstep_short_fvm_hevi(), scale_atmos_dyn_tstep_short_fvm_hivi::atmos_dyn_tstep_short_fvm_hivi(), and scale_atmos_dyn_tstep_tracer_fvm_heve::atmos_dyn_tstep_tracer_fvm_heve().

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
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:43
scale_comm
module COMMUNICATION
Definition: scale_comm.F90:23
scale_const::const_undef2
integer, parameter, public const_undef2
undefined value (INT2)
Definition: scale_const.F90:40
scale_const
module CONSTANT
Definition: scale_const.F90:14
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:36
Here is the call graph for this function:
Here is the caller graph for this function:
Generated by   doxygen 1.8.13