SCALE-RM
Functions/Subroutines
scale_atmos_phy_tb_common Module Reference

module ATMOSPHERE / Physics Turbulence More...

Functions/Subroutines

subroutine, public atmos_phy_tb_calc_strain_tensor (S33_C, S11_C, S22_C, S31_C, S12_C, S23_C, S12_Z, S23_X, S31_Y, S2, DENS, MOMZ, MOMX, MOMY, GSQRT, J13G, J23G, J33G, MAPF)
 
subroutine, public atmos_phy_tb_calc_flux_phi (qflx_phi, DENS, PHI, Kh, FACT, GSQRT, J13G, J23G, J33G, MAPF, a, b, c, dt, implicit, IIS, IIE, JJS, JJE)
 
subroutine, public atmos_phy_tb_diffusion_solver (phi, a, b, c, d, KE_TB)
 
subroutine, public atmos_phy_tb_calc_tend_momz (MOMZ_t_TB, QFLX_MOMZ, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
 
subroutine, public atmos_phy_tb_calc_tend_momx (MOMX_t_TB, QFLX_MOMX, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
 
subroutine, public atmos_phy_tb_calc_tend_momy (MOMY_t_TB, QFLX_MOMY, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
 
subroutine, public atmos_phy_tb_calc_tend_phi (phi_t_TB, QFLX_phi, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
 

Detailed Description

module ATMOSPHERE / Physics Turbulence

Description
Common routines for turbulence
Author
Team SCALE
History
  • 2014-12-13 (S.Nishizawa) [new]

Function/Subroutine Documentation

◆ atmos_phy_tb_calc_strain_tensor()

subroutine, public scale_atmos_phy_tb_common::atmos_phy_tb_calc_strain_tensor ( real(rp), dimension (ka,ia,ja), intent(out)  S33_C,
real(rp), dimension (ka,ia,ja), intent(out)  S11_C,
real(rp), dimension (ka,ia,ja), intent(out)  S22_C,
real(rp), dimension (ka,ia,ja), intent(out)  S31_C,
real(rp), dimension (ka,ia,ja), intent(out)  S12_C,
real(rp), dimension (ka,ia,ja), intent(out)  S23_C,
real(rp), dimension (ka,ia,ja), intent(out)  S12_Z,
real(rp), dimension (ka,ia,ja), intent(out)  S23_X,
real(rp), dimension (ka,ia,ja), intent(out)  S31_Y,
real(rp), dimension (ka,ia,ja), intent(out)  S2,
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,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,4), intent(in)  MAPF 
)
Parameters
[in]gsqrtvertical metrics {G}^1/2
[in]j13g(1,3) element of Jacobian matrix
[in]j23g(1,3) element of Jacobian matrix
[in]j33g(3,3) element of Jacobian matrix
[in]mapfmap factor

Definition at line 70 of file scale_atmos_phy_tb_common.F90.

References scale_debug::check(), scale_grid::grid_fdx, scale_grid::grid_fdy, scale_grid::grid_fdz, scale_grid::grid_rcdx, scale_grid::grid_rcdy, scale_grid::grid_rcdz, scale_grid::grid_rfdx, scale_grid::grid_rfdy, scale_grid::grid_rfdz, scale_gridtrans::i_uv, scale_gridtrans::i_uvz, scale_gridtrans::i_uy, scale_gridtrans::i_uyw, scale_gridtrans::i_uyz, scale_gridtrans::i_xv, scale_gridtrans::i_xvw, scale_gridtrans::i_xvz, scale_gridtrans::i_xy, scale_gridtrans::i_xyw, scale_gridtrans::i_xyz, 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, and scale_grid_index::ks.

Referenced by scale_atmos_phy_tb_d1980::atmos_phy_tb_d1980(), and scale_atmos_phy_tb_smg::atmos_phy_tb_smg().

70  use scale_gridtrans, only: &
71  i_xyz, &
72  i_xyw, &
73  i_uyw, &
74  i_xvw, &
75  i_uyz, &
76  i_xvz, &
77  i_uvz, &
78  i_xy, &
79  i_uy, &
80  i_xv, &
81  i_uv
82  use scale_grid, only: &
83  fdz => grid_fdz, &
84  fdx => grid_fdx, &
85  fdy => grid_fdy, &
86  rcdz => grid_rcdz, &
87  rcdx => grid_rcdx, &
88  rcdy => grid_rcdy, &
89  rfdz => grid_rfdz, &
90  rfdx => grid_rfdx, &
91  rfdy => grid_rfdy
92  implicit none
93 
94  real(RP), intent(out) :: S33_C (KA,IA,JA) ! (cell center)
95  real(RP), intent(out) :: S11_C (KA,IA,JA)
96  real(RP), intent(out) :: S22_C (KA,IA,JA)
97  real(RP), intent(out) :: S31_C (KA,IA,JA)
98  real(RP), intent(out) :: S12_C (KA,IA,JA)
99  real(RP), intent(out) :: S23_C (KA,IA,JA)
100  real(RP), intent(out) :: S12_Z (KA,IA,JA) ! (z edge or x-y plane)
101  real(RP), intent(out) :: S23_X (KA,IA,JA) ! (x edge or y-z plane)
102  real(RP), intent(out) :: S31_Y (KA,IA,JA) ! (y edge or z-x plane)
103  real(RP), intent(out) :: S2 (KA,IA,JA) ! |S|^2
104 
105  real(RP), intent(in) :: DENS (KA,IA,JA)
106  real(RP), intent(in) :: MOMZ (KA,IA,JA)
107  real(RP), intent(in) :: MOMX (KA,IA,JA)
108  real(RP), intent(in) :: MOMY (KA,IA,JA)
109 
110  real(RP), intent(in) :: GSQRT (KA,IA,JA,7)
111  real(RP), intent(in) :: J13G (KA,IA,JA,7)
112  real(RP), intent(in) :: J23G (KA,IA,JA,7)
113  real(RP), intent(in) :: J33G
114  real(RP), intent(in) :: MAPF (IA,JA,2,4)
115 
116  ! velocity
117  real(RP) :: VELZ_C (KA,IA,JA)
118  real(RP) :: VELZ_XY(KA,IA,JA)
119  real(RP) :: VELX_C (KA,IA,JA)
120  real(RP) :: VELX_YZ(KA,IA,JA)
121  real(RP) :: VELY_C (KA,IA,JA)
122  real(RP) :: VELY_ZX(KA,IA,JA)
123 
124  ! work space
125  real(RP) :: WORK_V(KA,IA,JA) ! work space (vertex)
126  real(RP) :: WORK_Z(KA,IA,JA) ! (z edge or x-y plane)
127  real(RP) :: WORK_X(KA,IA,JA) ! (x edge or y-z plane)
128  real(RP) :: WORK_Y(KA,IA,JA) ! (y edge or z-x plane)
129 
130  integer :: IIS, IIE, JJS, JJE
131  integer :: k, i, j
132 
133 #ifdef DEBUG
134  s33_c(:,:,:) = undef
135  s11_c(:,:,:) = undef
136  s22_c(:,:,:) = undef
137  s31_c(:,:,:) = undef
138  s12_c(:,:,:) = undef
139  s23_c(:,:,:) = undef
140  s12_z(:,:,:) = undef
141  s23_x(:,:,:) = undef
142  s31_y(:,:,:) = undef
143  s2(:,:,:) = undef
144 
145  velz_c(:,:,:) = undef
146  velz_xy(:,:,:) = undef
147  velx_c(:,:,:) = undef
148  velx_yz(:,:,:) = undef
149  vely_c(:,:,:) = undef
150  vely_zx(:,:,:) = undef
151 
152  work_v(:,:,:) = undef
153  work_z(:,:,:) = undef
154  work_x(:,:,:) = undef
155  work_y(:,:,:) = undef
156 #endif
157 
158  ! momentum -> velocity
159  do j = js-1, je+1
160  do i = is-1, ie+1
161  do k = ks, ke-1
162 #ifdef DEBUG
163  call check( __line__, momz(k,i,j) )
164  call check( __line__, dens(k+1,i,j) )
165  call check( __line__, dens(k,i,j) )
166 #endif
167  velz_xy(k,i,j) = 2.0_rp * momz(k,i,j) / ( dens(k+1,i,j)+dens(k,i,j) )
168  enddo
169  enddo
170  enddo
171 #ifdef DEBUG
172  i = iundef; j = iundef; k = iundef
173 #endif
174  do j = js-1, je+1
175  do i = is-1, ie+1
176  velz_xy(ke,i,j) = 0.0_rp
177  enddo
178  enddo
179 #ifdef DEBUG
180  i = iundef; j = iundef; k = iundef
181 #endif
182  do j = js-2, je+2
183  do i = is-2, ie+2
184  do k = ks+1, ke
185 #ifdef DEBUG
186  call check( __line__, momz(k,i,j) )
187  call check( __line__, momz(k-1,i,j) )
188  call check( __line__, dens(k,i,j) )
189 #endif
190  velz_c(k,i,j) = 0.5_rp * ( momz(k,i,j) + momz(k-1,i,j) ) / dens(k,i,j)
191  enddo
192  enddo
193  enddo
194 #ifdef DEBUG
195  i = iundef; j = iundef; k = iundef
196 #endif
197  do j = js-2, je+2
198  do i = is-2, ie+2
199 #ifdef DEBUG
200  call check( __line__, momz(ks,i,j) )
201  call check( __line__, dens(ks,i,j) )
202 #endif
203  velz_c(ks,i,j) = 0.5_rp * momz(ks,i,j) / dens(ks,i,j) ! MOMZ(KS-1,i,j) = 0
204  enddo
205 
206  enddo
207 #ifdef DEBUG
208  i = iundef; j = iundef; k = iundef
209 #endif
210 
211  do j = js-1, je+1
212  do i = is-2, ie+1
213  do k = ks, ke
214 #ifdef DEBUG
215  call check( __line__, momx(k,i,j) )
216  call check( __line__, dens(k,i+1,j) )
217  call check( __line__, dens(k,i,j) )
218 #endif
219  velx_yz(k,i,j) = 2.0_rp * momx(k,i,j) / ( dens(k,i+1,j)+dens(k,i,j) )
220  enddo
221  enddo
222  enddo
223 #ifdef DEBUG
224  i = iundef; j = iundef; k = iundef
225 #endif
226  do j = js-1, je+1
227  do i = is-2, ie+1
228  velx_yz(ke+1,i,j) = 0.0_rp
229  enddo
230  enddo
231 #ifdef DEBUG
232  i = iundef; j = iundef; k = iundef
233 #endif
234  do j = js-2, je+2
235  do i = is-1, ie+2
236  do k = ks, ke
237 #ifdef DEBUG
238  call check( __line__, momx(k,i,j) )
239  call check( __line__, momx(k,i-1,j) )
240  call check( __line__, dens(k,i,j) )
241 #endif
242  velx_c(k,i,j) = 0.5_rp * ( momx(k,i,j) + momx(k,i-1,j) ) / dens(k,i,j)
243  enddo
244  enddo
245  enddo
246 #ifdef DEBUG
247  i = iundef; j = iundef; k = iundef
248 #endif
249  !$omp parallel do default(none) &
250  !$omp shared(JS,JE,IS,IE,KS,KE,MOMY,DENS,VELY_ZX) &
251  !$omp private(i,j,k) OMP_SCHEDULE_ collapse(2)
252  do j = js-2, je+1
253  do i = is-1, ie+1
254  do k = ks, ke
255 #ifdef DEBUG
256  call check( __line__, momy(k,i,j) )
257  call check( __line__, dens(k,i,j+1) )
258  call check( __line__, dens(k,i,j) )
259 #endif
260  vely_zx(k,i,j) = 2.0_rp * momy(k,i,j) / ( dens(k,i,j+1)+dens(k,i,j) )
261  enddo
262  enddo
263  enddo
264 #ifdef DEBUG
265  i = iundef; j = iundef; k = iundef
266 #endif
267  do j = js-2, je+1
268  do i = is-1, ie+1
269  vely_zx(ke+1,i,j) = 0.0_rp
270  enddo
271  enddo
272 #ifdef DEBUG
273  i = iundef; j = iundef; k = iundef
274 #endif
275  do j = js-1, je+2
276  do i = is-2, ie+2
277  do k = ks, ke
278 #ifdef DEBUG
279  call check( __line__, momy(k,i,j) )
280  call check( __line__, momy(k,i,j-1) )
281  call check( __line__, dens(k,i,j) )
282 #endif
283  vely_c(k,i,j) = 0.5_rp * ( momy(k,i,j) + momy(k,i,j-1) ) / dens(k,i,j)
284  enddo
285  enddo
286  enddo
287 #ifdef DEBUG
288  i = iundef; j = iundef; k = iundef
289 #endif
290 
291  do jjs = js, je, jblock
292  jje = jjs+jblock-1
293  do iis = is, ie, iblock
294  iie = iis+iblock-1
295 
296 #ifdef DEBUG
297  work_z(:,:,:) = undef; work_x(:,:,:) = undef; work_y(:,:,:) = undef; work_v(:,:,:) = undef
298 #endif
299  ! w
300  ! (x-y plane; x,y,w)
301  ! WORK_Z = VELZ_XY
302  ! (y-z plane; u,y,z)
303  do j = jjs-1, jje+1
304  do i = iis-1, iie+1
305  do k = ks, ke
306 #ifdef DEBUG
307  call check( __line__, velz_c(k,i+1,j) )
308  call check( __line__, velz_c(k,i,j) )
309 #endif
310  work_x(k,i,j) = 0.5_rp * ( velz_c(k,i+1,j) + velz_c(k,i,j) )
311  enddo
312  enddo
313  enddo
314 #ifdef DEBUG
315  i = iundef; j = iundef; k = iundef
316 #endif
317  do j = jjs-1, jje+1
318  do i = iis-1, iie+1
319  work_x(ke+1,i,j) = 0.0_rp
320  enddo
321  enddo
322 #ifdef DEBUG
323  i = iundef; j = iundef; k = iundef
324 #endif
325  ! (z-x plane; x,v,z)
326  do j = jjs-1, jje+1
327  do i = iis-1, iie+1
328  do k = ks, ke
329 #ifdef DEBUG
330  call check( __line__, velz_c(k,i,j+1) )
331  call check( __line__, velz_c(k,i,j) )
332 #endif
333  work_y(k,i,j) = 0.5_rp * ( velz_c(k,i,j+1) + velz_c(k,i,j) )
334  enddo
335  enddo
336  enddo
337 #ifdef DEBUG
338  i = iundef; j = iundef; k = iundef
339 #endif
340  do j = jjs-1, jje+1
341  do i = iis-1, iie+1
342  work_y(ke+1,i,j) = 0.0_rp
343  enddo
344  enddo
345 #ifdef DEBUG
346  i = iundef; j = iundef; k = iundef
347 #endif
348 
349  ! dw/dz
350  ! (cell center; x,y,z)
351  do j = jjs-1, jje+1
352  do i = iis-1, iie+1
353  do k = ks+1, ke
354 #ifdef DEBUG
355  call check( __line__, velz_xy(k,i,j) )
356  call check( __line__, velz_xy(k-1,i,j) )
357  call check( __line__, rcdz(k) )
358 #endif
359  s33_c(k,i,j) = ( velz_xy(k,i,j) - velz_xy(k-1,i,j) ) * rcdz(k) &
360  * j33g / gsqrt(k,i,j,i_xyz)
361  enddo
362  enddo
363  enddo
364 #ifdef DEBUG
365  i = iundef; j = iundef; k = iundef
366 #endif
367  do j = jjs-1, jje+1
368  do i = iis-1, iie+1
369 #ifdef DEBUG
370  call check( __line__, velz_xy(ks,i,j) )
371  call check( __line__, gsqrt(ks,i,j,i_xyz) )
372  call check( __line__, rcdz(ks) )
373 #endif
374  s33_c(ks,i,j) = velz_xy(ks,i,j) * rcdz(ks) & ! VELZ_XY(KS-1,i,j) == 0
375  * j33g / gsqrt(ks,i,j,i_xyz)
376  enddo
377  enddo
378 #ifdef DEBUG
379  i = iundef; j = iundef; k = iundef
380 #endif
381 
382  ! 1/2 * dw/dx
383  ! (cell center; x,y,z)
384  do j = jjs-1, jje+1
385  do i = iis-1, iie+1
386  do k = ks+1, ke-1
387 #ifdef DEBUG
388  call check( __line__, velz_c(k,i+1,j) )
389  call check( __line__, velz_c(k,i-1,j) )
390  call check( __line__, gsqrt(k,i+1,j,i_xyz) )
391  call check( __line__, gsqrt(k,i-1,j,i_xyz) )
392  call check( __line__, velz_xy(k,i,j) )
393  call check( __line__, velz_xy(k-1,i,j) )
394  call check( __line__, j13g(k,i,j,i_xyw) )
395  call check( __line__, j13g(k-1,i,j,i_xyw) )
396  call check( __line__, fdx(i) )
397  call check( __line__, fdx(i-1) )
398 #endif
399  s31_c(k,i,j) = 0.5_rp * ( &
400  ( gsqrt(k,i+1,j,i_xyz)*velz_c(k,i+1,j) - gsqrt(k,i-1,j,i_xyz)*velz_c(k,i-1,j) ) / ( fdx(i) + fdx(i-1) ) &
401  + ( j13g(k,i,j,i_xyw)*velz_xy(k,i,j) - j13g(k-1,i,j,i_xyw)*velz_xy(k-1,i,j) ) * rcdz(k) &
402  ) * mapf(i,j,1,i_xy)
403 
404  enddo
405  enddo
406  enddo
407 #ifdef DEBUG
408  i = iundef; j = iundef; k = iundef
409 #endif
410  do j = jjs-1, jje+1
411  do i = iis-1, iie+1
412 #ifdef DEBUG
413  call check( __line__, velz_c(ks,i+1,j) )
414  call check( __line__, velz_c(ks,i-1,j) )
415  call check( __line__, gsqrt(ks,i+1,j,i_xyz) )
416  call check( __line__, gsqrt(ks,i-1,j,i_xyz) )
417  call check( __line__, velz_xy(ks,i,j) )
418  call check( __line__, j13g(ks,i,j,i_xyw) )
419  call check( __line__, velz_c(ke,i+1,j) )
420  call check( __line__, velz_c(ke,i-1,j) )
421  call check( __line__, gsqrt(ke,i+1,j,i_xyz) )
422  call check( __line__, gsqrt(ke,i-1,j,i_xyz) )
423  call check( __line__, velz_xy(ke,i,j) )
424  call check( __line__, j13g(ke,i,j,i_xyw) )
425  call check( __line__, fdx(i) )
426  call check( __line__, fdx(i-1) )
427 #endif
428  s31_c(ks,i,j) = 0.5_rp * ( &
429  ( gsqrt(ks,i+1,j,i_xyz)*velz_c(ks,i+1,j) - gsqrt(ks,i-1,j,i_xyz)*velz_c(ks,i-1,j) ) / ( fdx(i) + fdx(i-1) ) &
430  + ( j13g(ks,i,j,i_xyw)*velz_xy(ks,i,j) ) * rcdz(ks) &
431  ) * mapf(i,j,1,i_xy)
432  s31_c(ke,i,j) = 0.5_rp * ( &
433  ( gsqrt(ke,i+1,j,i_xyz)*velz_c(ke,i+1,j) - gsqrt(ke,i-1,j,i_xyz)*velz_c(ke,i-1,j) ) / ( fdx(i) + fdx(i-1) ) &
434  - ( j13g(ke-1,i,j,i_xyw)*velz_xy(ke-1,i,j) ) * rcdz(ke) &
435  ) * mapf(i,j,1,i_xy)
436  enddo
437  enddo
438 #ifdef DEBUG
439  i = iundef; j = iundef; k = iundef
440 #endif
441 
442  ! (y edge, u,y,w)
443  do j = jjs , jje
444  do i = iis-1, iie
445  do k = ks, ke-1
446 #ifdef DEBUG
447  call check( __line__, velz_xy(k,i+1,j) )
448  call check( __line__, velz_xy(k,i,j) )
449  call check( __line__, rfdx(i) )
450 #endif
451  s31_y(k,i,j) = 0.5_rp * ( &
452  ( gsqrt(k,i+1,j,i_xyw)*velz_xy(k,i+1,j) - gsqrt(k,i,j,i_xyw)*velz_xy(k,i,j) ) * rfdx(i) &
453  + ( j13g(k+1,i,j,i_uyz)*work_x(k+1,i,j) - j13g(k,i,j,i_uyz)*work_x(k,i,j)) * rfdz(k) &
454  ) * mapf(i,j,1,i_uy)
455  enddo
456  enddo
457  enddo
458 #ifdef DEBUG
459  i = iundef; j = iundef; k = iundef
460 #endif
461 
462  ! 1/2 * dw/dy
463  ! (cell center; x,y,z)
464  do j = jjs-1, jje+1
465  do i = iis-1, iie+1
466  do k = ks+1, ke-1
467 #ifdef DEBUG
468  call check( __line__, velz_c(k,i,j+1) )
469  call check( __line__, velz_c(k,i,j-1) )
470  call check( __line__, gsqrt(k,i,j+1,i_xyz) )
471  call check( __line__, gsqrt(k,i,j-1,i_xyz) )
472  call check( __line__, velz_xy(k,i,j) )
473  call check( __line__, velz_xy(k-1,i,j) )
474  call check( __line__, j23g(k,i,j,i_xyw) )
475  call check( __line__, j23g(k-1,i,j,i_xyw) )
476  call check( __line__, fdy(j) )
477  call check( __line__, fdy(j-1) )
478 #endif
479  s23_c(k,i,j) = 0.5_rp * ( &
480  ( gsqrt(k,i,j+1,i_xyz)*velz_c(k,i,j+1) - gsqrt(k,i,j-1,i_xyz)*velz_c(k,i,j-1) ) / ( fdy(j) + fdy(j-1) ) &
481  + ( j23g(k,i,j,i_xyw)*velz_xy(k,i,j) - j23g(k-1,i,j,i_xyw)*velz_xy(k-1,i,j) ) * rcdz(k) &
482  ) * mapf(i,j,2,i_xy)
483  enddo
484  enddo
485  enddo
486 #ifdef DEBUG
487  i = iundef; j = iundef; k = iundef
488 #endif
489  do j = jjs-1, jje+1
490  do i = iis-1, iie+1
491 #ifdef DEBUG
492  call check( __line__, velz_c(ks,i,j+1) )
493  call check( __line__, velz_c(ks,i,j-1) )
494  call check( __line__, gsqrt(ks,i,j+1,i_xyz) )
495  call check( __line__, gsqrt(ks,i,j-1,i_xyz) )
496  call check( __line__, velz_xy(ks,i,j) )
497  call check( __line__, j23g(ks,i,j,i_xyw) )
498  call check( __line__, velz_c(ke,i,j+1) )
499  call check( __line__, velz_c(ke,i,j-1) )
500  call check( __line__, gsqrt(ke,i,j+1,i_xyz) )
501  call check( __line__, gsqrt(ke,i,j-1,i_xyz) )
502  call check( __line__, velz_xy(ke,i,j) )
503  call check( __line__, j23g(ke,i,j,i_xyw) )
504  call check( __line__, fdy(j) )
505  call check( __line__, fdy(j-1) )
506 #endif
507  s23_c(ks,i,j) = 0.5_rp * ( &
508  ( gsqrt(ks,i,j+1,i_xyz)*velz_c(ks,i,j+1) - gsqrt(ks,i,j-1,i_xyz)*velz_c(ks,i,j-1) ) / ( fdy(j) + fdy(j-1) ) &
509  + ( j23g(ks,i,j,i_xyw)*velz_xy(ks,i,j) ) * rcdz(ks) &
510  ) * mapf(i,j,2,i_xy)
511  s23_c(ke,i,j) = 0.5_rp * ( &
512  ( gsqrt(ke,i,j+1,i_xyz)*velz_c(ke,i,j+1) - gsqrt(ke,i,j-1,i_xyz)*velz_c(ke,i,j-1) ) / ( fdy(j) + fdy(j-1) ) &
513  - ( j23g(ke-1,i,j,i_xyw)*velz_xy(ke-1,i,j) ) * rcdz(ke) &
514  ) * mapf(i,j,2,i_xy)
515  enddo
516  enddo
517 #ifdef DEBUG
518  i = iundef; j = iundef; k = iundef
519 #endif
520 
521  ! (x edge; x,v,w)
522  do j = jjs-1, jje
523  do i = iis , iie
524  do k = ks, ke-1
525 #ifdef DEBUG
526  call check( __line__, velz_xy(k,i,j+1) )
527  call check( __line__, velz_xy(k,i,j) )
528  call check( __line__, rfdy(j) )
529 #endif
530  s23_x(k,i,j) = 0.5_rp * ( &
531  ( gsqrt(k,i,j+1,i_xyw)*velz_xy(k,i,j+1) - gsqrt(k,i,j,i_xyw)*velz_xy(k,i,j) ) * rfdy(j) &
532  + ( j23g(k+1,i,j,i_xvz)*work_y(k+1,i,j) - j23g(k,i,j,i_xvz)*work_y(k,i,j) ) * rfdz(k) &
533  ) * mapf(i,j,2,i_xv)
534  enddo
535  enddo
536  enddo
537 #ifdef DEBUG
538  i = iundef; j = iundef; k = iundef
539 #endif
540 
541 #ifdef DEBUG
542  work_z(:,:,:) = undef; work_x(:,:,:) = undef; work_y(:,:,:) = undef; work_v(:,:,:) = undef
543 #endif
544  ! u
545  ! (x-y plane; x,y,w)
546  do j = jjs-1, jje+1
547  do i = iis-1, iie+1
548  do k = ks, ke-1
549 #ifdef DEBUG
550  call check( __line__, velx_c(k+1,i,j) )
551  call check( __line__, velx_c(k,i,j) )
552 #endif
553  work_z(k,i,j) = 0.5_rp * ( velx_c(k+1,i,j) + velx_c(k,i,j) )
554  enddo
555  enddo
556  enddo
557 #ifdef DEBUG
558  i = iundef; j = iundef; k = iundef
559 #endif
560  ! (y-z plane; u,y,z)
561  ! WORK_X = VELX_YZ
562  ! (z-x plane; x,v,z)
563  do j = jjs-1, jje
564  do i = iis-1, iie+1
565  do k = ks, ke
566 #ifdef DEBUG
567  call check( __line__, velx_c(k,i,j+1) )
568  call check( __line__, velx_c(k,i,j) )
569 #endif
570  work_y(k,i,j) = 0.5_rp * ( velx_c(k,i,j+1) + velx_c(k,i,j) )
571  enddo
572  enddo
573  enddo
574 #ifdef DEBUG
575  i = iundef; j = iundef; k = iundef
576 #endif
577  ! (vertex; u,v,w)
578  do j = jjs-1, jje
579  do i = iis-1, iie
580  do k = ks, ke-1
581 #ifdef DEBUG
582  call check( __line__, velx_yz(k,i,j) )
583  call check( __line__, velx_yz(k,i,j+1) )
584  call check( __line__, velx_yz(k+1,i,j) )
585  call check( __line__, velx_yz(k+1,i,j+1) )
586  call check( __line__, j23g(k ,i,j ,i_uvz) )
587  call check( __line__, j23g(k+1,i,j ,i_uvz) )
588  call check( __line__, j23g(k ,i,j+1,i_uvz) )
589  call check( __line__, j23g(k+1,i,j+1,i_uvz) )
590 #endif
591  work_v(k,i,j) = 0.25_rp &
592  * ( j23g(k ,i,j ,i_uyz)*velx_yz(k ,i,j ) &
593  + j23g(k+1,i,j ,i_uyz)*velx_yz(k+1,i,j ) &
594  + j23g(k ,i,j+1,i_uyz)*velx_yz(k ,i,j+1) &
595  + j23g(k+1,i,j+1,i_uyz)*velx_yz(k+1,i,j+1) )
596  enddo
597  enddo
598  enddo
599 #ifdef DEBUG
600  i = iundef; j = iundef; k = iundef
601 #endif
602 
603  ! du/dx
604  ! (cell center; x,y,z)
605  do j = jjs-1, jje+1
606  do i = iis-1, iie+1
607  do k = ks+1, ke-1
608 #ifdef DEBUG
609  call check( __line__, velx_yz(k,i,j) )
610  call check( __line__, velx_yz(k,i-1,j) )
611  call check( __line__, gsqrt(k,i,j,i_uyz) )
612  call check( __line__, gsqrt(k,i-1,j,i_uyz) )
613  call check( __line__, work_z(k,i,j) )
614  call check( __line__, work_z(k-1,i,j) )
615  call check( __line__, j13g(k,i,j,i_xyw) )
616  call check( __line__, j13g(k-1,i,j,i_xyw) )
617  call check( __line__, gsqrt(k,i,j,i_xyz) )
618  call check( __line__, rcdx(i) )
619 #endif
620  s11_c(k,i,j) = ( &
621  ( gsqrt(k,i,j,i_uyz)*velx_yz(k,i,j) - gsqrt(k,i-1,j,i_uyz)*velx_yz(k,i-1,j) ) * rcdx(i) &
622  + ( j13g(k,i,j,i_xyw)*work_z(k,i,j) - j13g(k-1,i,j,i_xyw)*work_z(k-1,i,j) ) * rcdz(k) &
623  ) * mapf(i,j,1,i_xy) / gsqrt(k,i,j,i_xyz)
624  enddo
625  enddo
626  enddo
627 #ifdef DEBUG
628  i = iundef; j = iundef; k = iundef
629 #endif
630  do j = jjs-1, jje+1
631  do i = iis-1, iie+1
632 #ifdef DEBUG
633  call check( __line__, velx_yz(ks,i,j) )
634  call check( __line__, velx_yz(ks,i-1,j) )
635  call check( __line__, gsqrt(ks,i,j,i_uyz) )
636  call check( __line__, gsqrt(ks,i-1,j,i_uyz) )
637  call check( __line__, velx_c(ks+1,i,j) )
638  call check( __line__, velx_c(ks,i,j) )
639  call check( __line__, j13g(ks+1,i,j,i_xyz) )
640  call check( __line__, j13g(ks,i,j,i_xyz) )
641  call check( __line__, gsqrt(ks,i,j,i_xyz) )
642  call check( __line__, velx_yz(ke,i,j) )
643  call check( __line__, velx_yz(ke,i-1,j) )
644  call check( __line__, gsqrt(ke,i,j,i_uyz) )
645  call check( __line__, gsqrt(ke,i-1,j,i_uyz) )
646  call check( __line__, velx_c(ke,i,j) )
647  call check( __line__, velx_c(ke-1,i,j) )
648  call check( __line__, j13g(ke,i,j,i_xyz) )
649  call check( __line__, j13g(ke-1,i,j,i_xyz) )
650  call check( __line__, gsqrt(ke,i,j,i_xyz) )
651  call check( __line__, rcdx(i) )
652 #endif
653  s11_c(ks,i,j) = ( &
654  ( gsqrt(ks,i,j,i_uyz)*velx_yz(ks,i,j) - gsqrt(ks,i-1,j,i_uyz)*velx_yz(ks,i-1,j) ) * rcdx(i) &
655  + ( j13g(ks+1,i,j,i_xyz)*velx_c(ks+1,i,j) - j13g(ks,i,j,i_xyz)*velx_c(ks,i,j) ) * rfdz(ks) &
656  ) * mapf(i,j,1,i_xy) / gsqrt(ks,i,j,i_xyz)
657  s11_c(ke,i,j) = ( &
658  ( gsqrt(ke,i,j,i_uyz)*velx_yz(ke,i,j) - gsqrt(ke,i-1,j,i_uyz)*velx_yz(ke,i-1,j) ) * rcdx(i) &
659  + ( j13g(ke,i,j,i_xyz)*velx_c(ke,i,j) - j13g(ke-1,i,j,i_xyz)*velx_c(ke-1,i,j) ) * rfdz(ke-1) &
660  ) * mapf(i,j,1,i_xy) / gsqrt(ke,i,j,i_xyz)
661  enddo
662  enddo
663 #ifdef DEBUG
664  i = iundef; j = iundef; k = iundef
665 #endif
666 
667  ! 1/2 * du/dz
668  ! (cell center; x,y,z)
669  do j = jjs-1, jje+1
670  do i = iis-1, iie+1
671  do k = ks+1, ke-1
672 #ifdef DEBUG
673  call check( __line__, s31_c(k,i,j) )
674  call check( __line__, velx_c(k+1,i,j) )
675  call check( __line__, velx_c(k-1,i,j) )
676  call check( __line__, fdz(k) )
677  call check( __line__, fdz(k-1) )
678 #endif
679  s31_c(k,i,j) = ( s31_c(k,i,j) & ! dw/dx
680  + 0.5_rp * ( velx_c(k+1,i,j) - velx_c(k-1,i,j) ) * j33g / ( fdz(k) + fdz(k-1) ) &
681  ) / gsqrt(k,i,j,i_xyz)
682  enddo
683  enddo
684  enddo
685 #ifdef DEBUG
686  i = iundef; j = iundef; k = iundef
687 #endif
688  do j = jjs-1, jje+1
689  do i = iis-1, iie+1
690 #ifdef DEBUG
691  call check( __line__, s31_c(ks,i,j) )
692  call check( __line__, velx_c(ks+1,i,j) )
693  call check( __line__, velx_c(ks,i,j) )
694  call check( __line__, rfdz(ks) )
695  call check( __line__, s31_c(ke,i,j) )
696  call check( __line__, velx_c(ke,i,j) )
697  call check( __line__, velx_c(ke-1,i,j) )
698  call check( __line__, rfdz(ke-1) )
699 #endif
700  s31_c(ks,i,j) = ( s31_c(ks,i,j) &
701  + 0.5_rp * ( velx_c(ks+1,i,j) - velx_c(ks,i,j) ) * j33g * rfdz(ks) &
702  ) / gsqrt(ks,i,j,i_xyz)
703  s31_c(ke,i,j) = ( s31_c(ke,i,j) &
704  + 0.5_rp * ( velx_c(ke,i,j) - velx_c(ke-1,i,j) ) * j33g * rfdz(ke-1) &
705  ) / gsqrt(ke,i,j,i_xyz)
706  enddo
707  enddo
708 #ifdef DEBUG
709  i = iundef; j = iundef; k = iundef
710 #endif
711  ! (y edge; u,y,w)
712  do j = jjs , jje
713  do i = iis-1, iie
714  do k = ks, ke-1
715 #ifdef DEBUG
716  call check( __line__, s31_y(k,i,j) )
717  call check( __line__, velx_yz(k+1,i,j) )
718  call check( __line__, velx_yz(k,i,j) )
719  call check( __line__, rfdz(k) )
720 #endif
721  s31_y(k,i,j) = ( s31_y(k,i,j) & ! dw/dx
722  + 0.5_rp * ( velx_yz(k+1,i,j) - velx_yz(k,i,j) ) * j33g * rfdz(k) &
723  ) / gsqrt(k,i,j,i_uyw)
724  enddo
725  enddo
726  enddo
727 #ifdef DEBUG
728  i = iundef; j = iundef; k = iundef
729 #endif
730 
731  ! 1/2 * du/dy
732  ! (cell center; x,y,z)
733  do j = jjs-1, jje+1
734  do i = iis-1, iie+1
735  do k = ks+1, ke-1
736 #ifdef DEBUG
737  call check( __line__, velx_c(k,i,j+1) )
738  call check( __line__, velx_c(k,i,j-1) )
739  call check( __line__, gsqrt(k,i,j+1,i_xyz) )
740  call check( __line__, gsqrt(k,i,j-1,i_xyz) )
741  call check( __line__, work_z(k,i,j) )
742  call check( __line__, work_z(k-1,i,j) )
743  call check( __line__, j23g(k,i,j,i_xyw) )
744  call check( __line__, j23g(k-1,i,j,i_xyw) )
745 #endif
746  s12_c(k,i,j) = 0.5_rp * ( &
747  ( gsqrt(k,i,j+1,i_xyz)*velx_c(k,i,j+1) - gsqrt(k,i,j-1,i_xyz)*velx_c(k,i,j-1) ) / ( fdy(j) + fdy(j-1) ) &
748  + ( j23g(k,i,j,i_xyw)*work_z(k,i,j) - j23g(k-1,i,j,i_xyw)*work_z(k-1,i,j) ) * rcdz(k) &
749  ) * mapf(i,j,2,i_xy) / gsqrt(k,i,j,i_xyz)
750  enddo
751  enddo
752  enddo
753 #ifdef DEBUG
754  i = iundef; j = iundef; k = iundef
755 #endif
756  do j = jjs-1, jje+1
757  do i = iis-1, iie+1
758 #ifdef DEBUG
759  call check( __line__, velx_c(ks,i,j+1) )
760  call check( __line__, velx_c(ks,i,j-1) )
761  call check( __line__, gsqrt(ks,i,j+1,i_xyz) )
762  call check( __line__, gsqrt(ks,i,j-1,i_xyz) )
763  call check( __line__, velx_c(ks+1,i,j) )
764  call check( __line__, velx_c(ks,i,j) )
765  call check( __line__, j23g(ks+1,i,j,i_xyz) )
766  call check( __line__, j23g(ks,i,j,i_xyz) )
767  call check( __line__, gsqrt(ks,i,j,i_xyz) )
768  call check( __line__, velx_c(ke,i,j+1) )
769  call check( __line__, velx_c(ke,i,j-1) )
770  call check( __line__, gsqrt(ke,i,j+1,i_xyz) )
771  call check( __line__, gsqrt(ke,i,j-1,i_xyz) )
772  call check( __line__, velx_c(ke,i,j) )
773  call check( __line__, velx_c(ke-1,i,j) )
774  call check( __line__, j23g(ke,i,j,i_xyz) )
775  call check( __line__, j23g(ke-1,i,j,i_xyz) )
776  call check( __line__, gsqrt(ke,i,j,i_xyz) )
777  call check( __line__, fdy(j) )
778  call check( __line__, fdy(j-1) )
779 #endif
780  s12_c(ks,i,j) = 0.5_rp * ( &
781  ( gsqrt(ks,i,j+1,i_xyz)*velx_c(ks,i,j+1) - gsqrt(ks,i,j-1,i_xyz)*velx_c(ks,i,j-1) ) / ( fdy(j) + fdy(j-1) ) &
782  + ( j23g(ks+1,i,j,i_xyz)*velx_c(ks+1,i,j) - j23g(ks,i,j,i_xyz)*velx_c(ks,i,j) ) * rfdz(ks) &
783  ) * mapf(i,j,2,i_xy) / gsqrt(ks,i,j,i_xyz)
784  s12_c(ke,i,j) = 0.5_rp * ( &
785  ( gsqrt(ke,i,j+1,i_xyz)*velx_c(ke,i,j+1) - gsqrt(ke,i,j-1,i_xyz)*velx_c(ke,i,j-1) ) / ( fdy(j) + fdy(j-1) ) &
786  + ( j23g(ke,i,j,i_xyz)*velx_c(ke,i,j) - j23g(ke-1,i,j,i_xyz)*velx_c(ke-1,i,j) ) * rfdz(ke-1) &
787  ) * mapf(i,j,2,i_xy) / gsqrt(ke,i,j,i_xyz)
788  enddo
789  enddo
790 #ifdef DEBUG
791  i = iundef; j = iundef; k = iundef
792 #endif
793 
794  ! (z edge; u,v,z)
795  do j = jjs-1, jje
796  do i = iis-1, iie
797  do k = ks+1, ke-1
798 #ifdef DEBUG
799  call check( __line__, velx_yz(k,i,j+1) )
800  call check( __line__, velx_yz(k,i,j) )
801  call check( __line__, work_v(k,i,j) )
802  call check( __line__, work_v(k-1,i,j) )
803  call check( __line__, rfdy(j) )
804 #endif
805  s12_z(k,i,j) = 0.5_rp * ( &
806  ( gsqrt(k,i,j+1,i_uyz)*velx_yz(k,i,j+1) - gsqrt(k,i,j,i_uyz)*velx_yz(k,i,j) ) * rfdy(j) &
807  + ( work_v(k,i,j) - work_v(k-1,i,j) ) * rcdz(k) &
808  ) * mapf(i,j,2,i_uv)
809  enddo
810  enddo
811  enddo
812 #ifdef DEBUG
813  i = iundef; j = iundef; k = iundef
814 #endif
815  do j = jjs-1, jje
816  do i = iis-1, iie
817 #ifdef DEBUG
818  call check( __line__, velx_yz(ks,i,j+1) )
819  call check( __line__, velx_yz(ks,i,j) )
820  call check( __line__, velx_yz(ks+1,i,j) )
821  call check( __line__, velx_yz(ks+1,i,j+1) )
822  call check( __line__, j23g(ks+1,i,j,i_uvz) )
823  call check( __line__, j23g(ks ,i,j,i_uvz) )
824  call check( __line__, velx_yz(ke,i,j+1) )
825  call check( __line__, velx_yz(ke,i,j) )
826  call check( __line__, velx_yz(ke-1,i,j) )
827  call check( __line__, velx_yz(ke-1,i,j+1) )
828  call check( __line__, j23g(ke ,i,j,i_uvz) )
829  call check( __line__, j23g(ke-1,i,j,i_uvz) )
830 #endif
831  s12_z(ks,i,j) = 0.25_rp * ( &
832  ( gsqrt(ks,i,j+1,i_uyz)*velx_yz(ks,i,j+1) - gsqrt(ks,i,j,i_uyz)*velx_yz(ks,i,j) ) * rfdy(j) &
833  + ( j23g(ks+1,i,j,i_uvz) * ( velx_yz(ks+1,i,j) + velx_yz(ks+1,i,j+1) ) &
834  - j23g(ks ,i,j,i_uvz) * ( velx_yz(ks ,i,j) + velx_yz(ks ,i,j+1) ) ) * rfdz(ks) &
835  ) * mapf(i,j,2,i_uv)
836  s12_z(ke,i,j) = 0.25_rp * ( &
837  ( gsqrt(ke,i,j+1,i_uyz)*velx_yz(ke,i,j+1) - gsqrt(ke,i,j,i_uyz)*velx_yz(ke,i,j) ) * rfdy(j) &
838  + ( j23g(ke ,i,j,i_uvz) * ( velx_yz(ke ,i,j) + velx_yz(ke ,i,j+1) ) &
839  - j23g(ke-1,i,j,i_uvz) * ( velx_yz(ke-1,i,j) + velx_yz(ke-1,i,j+1) ) ) * rfdz(ke-1) &
840  ) * mapf(i,j,2,i_uv)
841  enddo
842  enddo
843 #ifdef DEBUG
844  i = iundef; j = iundef; k = iundef
845 #endif
846 
847 #ifdef DEBUG
848  work_z(:,:,:) = undef; work_x(:,:,:) = undef; work_y(:,:,:) = undef; work_v(:,:,:) = undef
849 #endif
850  ! v
851  ! (x-y plane; x,y,w)
852  do j = jjs-1, jje+1
853  do i = iis-1, iie+1
854  do k = ks, ke-1
855 #ifdef DEBUG
856  call check( __line__, vely_c(k+1,i,j) )
857  call check( __line__, vely_c(k,i,j) )
858 #endif
859  work_z(k,i,j) = 0.5_rp * ( vely_c(k+1,i,j) + vely_c(k,i,j) )
860  enddo
861  enddo
862  enddo
863 #ifdef DEBUG
864  i = iundef; j = iundef; k = iundef
865 #endif
866  ! (y-z plane; u,y,z)
867  do j = jjs-1, jje+1
868  do i = iis-1, iie
869  do k = ks, ke
870 #ifdef DEBUG
871  call check( __line__, vely_c(k,i+1,j) )
872  call check( __line__, vely_c(k,i,j) )
873 #endif
874  work_x(k,i,j) = 0.5_rp * ( velz_c(k,i+1,j) + velz_c(k,i,j) )
875  enddo
876  enddo
877  enddo
878 #ifdef DEBUG
879  i = iundef; j = iundef; k = iundef
880 #endif
881  ! (z-x plane; x,v,z)
882  ! WORK_Y = VELY_ZX
883  ! (vertex; u,v,w)
884  do j = jjs-1, jje
885  do i = iis-1, iie
886  do k = ks, ke-1
887 #ifdef DEBUG
888  call check( __line__, vely_zx(k,i,j) )
889  call check( __line__, vely_zx(k+1,i,j) )
890  call check( __line__, vely_zx(k,i+1,j) )
891  call check( __line__, vely_zx(k+1,i+1,j) )
892 #endif
893  work_v(k,i,j) = 0.25_rp &
894  * ( j13g(k ,i ,j,i_xvz)*vely_zx(k ,i ,j) &
895  + j13g(k+1,i ,j,i_xvz)*vely_zx(k+1,i ,j) &
896  + j13g(k ,i+1,j,i_xvz)*vely_zx(k ,i+1,j) &
897  + j13g(k+1,i+1,j,i_xvz)*vely_zx(k+1,i+1,j) )
898  enddo
899  enddo
900  enddo
901 #ifdef DEBUG
902  i = iundef; j = iundef; k = iundef
903 #endif
904 
905  ! dv/dy
906  ! (cell center; x,y,z)
907  do j = jjs-1, jje+1
908  do i = iis-1, iie+1
909  do k = ks+1, ke-1
910 #ifdef DEBUG
911  call check( __line__, vely_zx(k,i,j) )
912  call check( __line__, vely_zx(k,i,j-1) )
913  call check( __line__, gsqrt(k,i,j,i_xvz) )
914  call check( __line__, gsqrt(k,i,j-1,i_xvz) )
915  call check( __line__, work_z(k,i,j) )
916  call check( __line__, work_z(k-1,i,j) )
917  call check( __line__, j23g(k,i,j,i_xyw) )
918  call check( __line__, j23g(k-1,i,j,i_xyw) )
919  call check( __line__, rcdy(j) )
920 #endif
921  s22_c(k,i,j) = ( &
922  ( gsqrt(k,i,j,i_xvz)*vely_zx(k,i,j) - gsqrt(k,i,j-1,i_xvz)*vely_zx(k,i,j-1) ) * rcdy(j) &
923  + ( j23g(k,i,j,i_xyw)*work_z(k,i,j) - j23g(k-1,i,j,i_xyw)*work_z(k-1,i,j) ) * rcdz(k) &
924  ) * mapf(i,j,2,i_xy) / gsqrt(k,i,j,i_xyz)
925  enddo
926  enddo
927  enddo
928 #ifdef DEBUG
929  i = iundef; j = iundef; k = iundef
930 #endif
931  do j = jjs-1, jje+1
932  do i = iis-1, iie+1
933 #ifdef DEBUG
934  call check( __line__, vely_zx(ks,i,j) )
935  call check( __line__, vely_zx(ks,i,j-1) )
936  call check( __line__, gsqrt(ks,i,j,i_xvz) )
937  call check( __line__, gsqrt(ks,i,j-1,i_xvz) )
938  call check( __line__, vely_c(ks+1,i,j) )
939  call check( __line__, vely_c(ks,i,j) )
940  call check( __line__, j23g(ks+1,i,j,i_xyz) )
941  call check( __line__, j23g(ks,i,j,i_xyz) )
942  call check( __line__, rcdy(j) )
943  call check( __line__, vely_zx(ke,i,j) )
944  call check( __line__, vely_zx(ke,i,j-1) )
945  call check( __line__, gsqrt(ke,i,j,i_xvz) )
946  call check( __line__, gsqrt(ke,i,j-1,i_xvz) )
947  call check( __line__, vely_c(ke,i,j) )
948  call check( __line__, vely_c(ke-1,i,j) )
949  call check( __line__, j23g(ke,i,j,i_xyz) )
950  call check( __line__, j23g(ke-1,i,j,i_xyz) )
951 #endif
952  s22_c(ks,i,j) = ( &
953  ( gsqrt(ks,i,j,i_xvz)*vely_zx(ks,i,j) - gsqrt(ks,i,j-1,i_xvz)*vely_zx(ks,i,j-1) ) * rcdy(j) &
954  + ( j23g(ks+1,i,j,i_xyz)*vely_c(ks+1,i,j) - j23g(ks,i,j,i_xyz)*vely_c(ks,i,j) ) * rfdz(ks) &
955  ) * mapf(i,j,2,i_xy) / gsqrt(ks,i,j,i_xyz)
956  s22_c(ke,i,j) = ( &
957  ( gsqrt(ke,i,j,i_xvz)*vely_zx(ke,i,j) - gsqrt(ke,i,j-1,i_xvz)*vely_zx(ke,i,j-1) ) * rcdy(j) &
958  + ( j23g(ke,i,j,i_xyz)*vely_c(ke,i,j) - j23g(ke-1,i,j,i_xyz)*vely_c(ke-1,i,j) ) * rfdz(ke-1) &
959  ) * mapf(i,j,2,i_xy) / gsqrt(ke,i,j,i_xyz)
960  enddo
961  enddo
962 #ifdef DEBUG
963  i = iundef; j = iundef; k = iundef
964 #endif
965 
966  ! 1/2 * dv/dx
967  ! (cell center; x,y,z)
968  do j = jjs-1, jje+1
969  do i = iis-1, iie+1
970  do k = ks+1, ke-1
971 #ifdef DEBUG
972  call check( __line__, s12_c(k,i,j) )
973  call check( __line__, vely_c(k,i+1,j) )
974  call check( __line__, vely_c(k,i-1,j) )
975  call check( __line__, gsqrt(k,i+1,j,i_xyz) )
976  call check( __line__, gsqrt(k,i-1,j,i_xyz) )
977  call check( __line__, work_z(k,i,j) )
978  call check( __line__, work_z(k-1,i,j) )
979  call check( __line__, j13g(k,i,j,i_xyw) )
980  call check( __line__, j13g(k-1,i,j,i_xyw) )
981  call check( __line__, gsqrt(k,i,j,i_xyz) )
982  call check( __line__, fdx(i) )
983  call check( __line__, fdx(i-1) )
984 #endif
985  s12_c(k,i,j) = ( s12_c(k,i,j) & ! du/dy
986  + 0.5_rp * ( &
987  ( gsqrt(k,i+1,j,i_xyz)*vely_c(k,i+1,j) - gsqrt(k,i-1,j,i_xyz)*vely_c(k,i-1,j) ) / ( fdx(i) + fdx(i-1) ) &
988  + ( j13g(k,i,j,i_xyw)*work_z(k,i,j) - j13g(k-1,i,j,i_xyw)*work_z(k-1,i,j) ) * rcdz(k) ) * mapf(i,j,1,i_xy) &
989  ) / gsqrt(k,i,j,i_xyz)
990  enddo
991  enddo
992  enddo
993 #ifdef DEBUG
994  i = iundef; j = iundef; k = iundef
995 #endif
996  do j = jjs-1, jje+1
997  do i = iis-1, iie+1
998 #ifdef DEBUG
999  call check( __line__, s12_c(ks,i,j) )
1000  call check( __line__, vely_c(ks,i+1,j) )
1001  call check( __line__, vely_c(ks,i-1,j) )
1002  call check( __line__, gsqrt(ks,i+1,j,i_xyz) )
1003  call check( __line__, gsqrt(ks,i-1,j,i_xyz) )
1004  call check( __line__, vely_c(ks+1,i,j) )
1005  call check( __line__, vely_c(ks,i,j) )
1006  call check( __line__, j13g(ks+1,i,j,i_xyz) )
1007  call check( __line__, j13g(ks,i,j,i_xyz) )
1008  call check( __line__, gsqrt(ks,i,j,i_xyz) )
1009  call check( __line__, s12_c(ke,i,j) )
1010  call check( __line__, vely_c(ke,i+1,j) )
1011  call check( __line__, vely_c(ke,i-1,j) )
1012  call check( __line__, gsqrt(ke,i+1,j,i_xyz) )
1013  call check( __line__, gsqrt(ke,i-1,j,i_xyz) )
1014  call check( __line__, vely_c(ke,i,j) )
1015  call check( __line__, vely_c(ke-1,i,j) )
1016  call check( __line__, j13g(ke,i,j,i_xyz) )
1017  call check( __line__, j13g(ke-1,i,j,i_xyz) )
1018  call check( __line__, gsqrt(ke,i,j,i_xyz) )
1019  call check( __line__, fdx(i) )
1020  call check( __line__, fdx(i-1) )
1021 #endif
1022  s12_c(ks,i,j) = ( s12_c(ks,i,j) & ! du/dy
1023  + 0.5_rp * ( &
1024  ( gsqrt(ks,i+1,j,i_xyz)*vely_c(ks,i+1,j) - gsqrt(ks,i-1,j,i_xyz)*vely_c(ks,i-1,j) ) / ( fdx(i) + fdx(i-1) ) &
1025  + ( j13g(ks+1,i,j,i_xyz)*vely_c(ks+1,i,j) - j13g(ks,i,j,i_xyz)*vely_c(ks,i,j) ) * rfdz(ks) ) &
1026  * mapf(i,j,1,i_xy) &
1027  ) / gsqrt(ks,i,j,i_xyz)
1028  s12_c(ke,i,j) = ( s12_c(ke,i,j) & ! du/dy
1029  + 0.5_rp * ( &
1030  ( gsqrt(ke,i+1,j,i_xyz)*vely_c(ke,i+1,j) - gsqrt(ke,i-1,j,i_xyz)*vely_c(ke,i-1,j) ) / ( fdx(i) + fdx(i-1) ) &
1031  + ( j13g(ke,i,j,i_xyz)*vely_c(ke,i,j) - j13g(ke-1,i,j,i_xyz)*vely_c(ke-1,i,j) ) * rfdz(ke-1) ) &
1032  * mapf(i,j,1,i_xy) &
1033  ) / gsqrt(ke,i,j,i_xyz)
1034  enddo
1035  enddo
1036 #ifdef DEBUG
1037  i = iundef; j = iundef; k = iundef
1038 #endif
1039  ! (z edge; u,v,z)
1040  do j = jjs-1, jje
1041  do i = iis-1, iie
1042  do k = ks+1, ke-1
1043 #ifdef DEBUG
1044  call check( __line__, s12_z(k,i,j) )
1045  call check( __line__, vely_zx(k,i+1,j) )
1046  call check( __line__, vely_zx(k,i,j) )
1047  call check( __line__, work_v(k,i,j) )
1048  call check( __line__, work_v(k-1,i,j) )
1049  call check( __line__, rfdx(i) )
1050 #endif
1051  s12_z(k,i,j) = ( s12_z(k,i,j) &
1052  + 0.5_rp * ( &
1053  ( gsqrt(k,i+1,j,i_xvz)*vely_zx(k,i+1,j) - gsqrt(k,i,j,i_xvz)*vely_zx(k,i,j) ) * rfdx(i) &
1054  + ( work_v(k,i,j) - work_v(k-1,i,j) ) * rcdz(k) ) * mapf(i,j,1,i_uv) &
1055  ) / gsqrt(k,i,j,i_uvz)
1056  enddo
1057  enddo
1058  enddo
1059 #ifdef DEBUG
1060  i = iundef; j = iundef; k = iundef
1061 #endif
1062  do j = jjs-1, jje
1063  do i = iis-1, iie
1064 #ifdef DEBUG
1065  call check( __line__, s12_z(ks,i,j) )
1066  call check( __line__, vely_zx(ks,i+1,j) )
1067  call check( __line__, vely_zx(ks,i,j) )
1068  call check( __line__, vely_zx(ks+1,i,j) )
1069  call check( __line__, vely_zx(ks+1,i+1,j) )
1070  call check( __line__, s12_z(ke,i,j) )
1071  call check( __line__, vely_zx(ke,i+1,j) )
1072  call check( __line__, vely_zx(ke,i,j) )
1073  call check( __line__, vely_zx(ke-1,i,j) )
1074  call check( __line__, vely_zx(ke-1,i+1,j) )
1075  call check( __line__, rfdx(i) )
1076 #endif
1077  s12_z(ks,i,j) = ( s12_z(ks,i,j) &
1078  + 0.5_rp * ( &
1079  ( gsqrt(ks,i+1,j,i_xvz)*vely_zx(ks,i+1,j) - gsqrt(ks,i,j,i_xvz)*vely_zx(ks,i,j) ) * rfdx(i) &
1080  + ( j13g(ks+1,i,j,i_uvz) * ( vely_zx(ks+1,i,j) + vely_zx(ks+1,i+1,j) ) &
1081  - j13g(ks ,i,j,i_uvz) * ( vely_zx(ks ,i,j) + vely_zx(ks ,i+1,j) ) ) * rfdz(ks) ) * mapf(i,j,1,i_uv) &
1082  ) / gsqrt(ks,i,j,i_uvz)
1083  s12_z(ke,i,j) = ( s12_z(ke,i,j) &
1084  + 0.5_rp * ( &
1085  ( gsqrt(ke,i+1,j,i_xvz)*vely_zx(ke,i+1,j) - gsqrt(ke,i,j,i_xvz)*vely_zx(ke,i,j) ) * rfdx(i) &
1086  + ( j13g(ke ,i,j,i_uvz) * ( vely_zx(ke ,i,j) + vely_zx(ke ,i+1,j) ) &
1087  - j13g(ke-1,i,j,i_uvz) * ( vely_zx(ke-1,i,j) + vely_zx(ke-1,i+1,j) ) ) * rfdz(ke-1) ) * mapf(i,j,1,i_uv) &
1088  ) / gsqrt(ke,i,j,i_uvz)
1089  enddo
1090  enddo
1091 #ifdef DEBUG
1092  i = iundef; j = iundef; k = iundef
1093 #endif
1094 
1095  ! 1/2 * dv/dz
1096  ! (cell center; x,y,z)
1097  do j = jjs-1, jje+1
1098  do i = iis-1, iie+1
1099  do k = ks+1, ke-1
1100 #ifdef DEBUG
1101  call check( __line__, s23_c(k,i,j) )
1102  call check( __line__, vely_c(k+1,i,j) )
1103  call check( __line__, vely_c(k-1,i,j) )
1104  call check( __line__, fdz(k) )
1105  call check( __line__, fdz(k-1) )
1106 #endif
1107  s23_c(k,i,j) = ( s23_c(k,i,j) & ! dw/dy
1108  + 0.5_rp * ( vely_c(k+1,i,j) - vely_c(k-1,i,j) ) * j33g / ( fdz(k) + fdz(k-1) ) &
1109  ) / gsqrt(k,i,j,i_xyz)
1110  enddo
1111  enddo
1112  enddo
1113 #ifdef DEBUG
1114  i = iundef; j = iundef; k = iundef
1115 #endif
1116  do j = jjs-1, jje+1
1117  do i = iis-1, iie+1
1118 #ifdef DEBUG
1119  call check( __line__, s23_c(ks,i,j) )
1120  call check( __line__, vely_c(ks+1,i,j) )
1121  call check( __line__, vely_c(ks,i,j) )
1122  call check( __line__, rfdz(ks) )
1123  call check( __line__, s23_c(ke,i,j) )
1124  call check( __line__, vely_c(ke,i,j) )
1125  call check( __line__, vely_c(ke-1,i,j) )
1126  call check( __line__, rfdz(ke-1) )
1127 #endif
1128  s23_c(ks,i,j) = ( s23_c(ks,i,j) &
1129  + 0.5_rp * ( vely_c(ks+1,i,j) - vely_c(ks,i,j) ) * j33g * rfdz(ks) &
1130  ) / gsqrt(ks,i,j,i_xyz)
1131  s23_c(ke,i,j) = ( s23_c(ke,i,j) &
1132  + 0.5_rp * ( vely_c(ke,i,j) - vely_c(ke-1,i,j) ) * j33g * rfdz(ke-1) &
1133  ) / gsqrt(ke,i,j,i_xyz)
1134  enddo
1135  enddo
1136 #ifdef DEBUG
1137  i = iundef; j = iundef; k = iundef
1138 #endif
1139 
1140  ! (x edge; x,v,w)
1141  do j = jjs-1, jje
1142  do i = iis , iie
1143  do k = ks, ke-1
1144 #ifdef DEBUG
1145  call check( __line__, s23_x(k,i,j) )
1146  call check( __line__, vely_zx(k+1,i,j) )
1147  call check( __line__, vely_zx(k,i,j) )
1148  call check( __line__, rfdz(k) )
1149 #endif
1150  s23_x(k,i,j) = ( s23_x(k,i,j) &
1151  + 0.5_rp * ( vely_zx(k+1,i,j) - vely_zx(k,i,j) ) * j33g * rfdz(k) &
1152  ) / gsqrt(k,i,j,i_xvw)
1153  enddo
1154  enddo
1155  enddo
1156 #ifdef DEBUG
1157  i = iundef; j = iundef; k = iundef
1158 #endif
1159 
1160 
1161  ! |S|^2 = 2*Sij*Sij
1162 #ifdef DEBUG
1163  work_z(:,:,:) = undef; work_x(:,:,:) = undef; work_y(:,:,:) = undef
1164 #endif
1165  ! (cell center)
1166  !$omp parallel do default(none) &
1167  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,S11_C,S22_C,S33_C,S31_C,S12_C,S23_C,S2) &
1168  !$omp private(i,j,k) OMP_SCHEDULE_ collapse(2)
1169  do j = jjs-1, jje+1
1170  do i = iis-1, iie+1
1171  do k = ks, ke
1172 #ifdef DEBUG
1173  call check( __line__, s11_c(k,i,j) )
1174  call check( __line__, s22_c(k,i,j) )
1175  call check( __line__, s33_c(k,i,j) )
1176  call check( __line__, s31_c(k,i,j) )
1177  call check( __line__, s12_c(k,i,j) )
1178  call check( __line__, s23_c(k,i,j) )
1179 #endif
1180  s2(k,i,j) = max( 1e-10_rp, &
1181  2.0_rp * ( s11_c(k,i,j)**2 + s22_c(k,i,j)**2 + s33_c(k,i,j)**2 ) &
1182  + 4.0_rp * ( s31_c(k,i,j)**2 + s12_c(k,i,j)**2 + s23_c(k,i,j)**2 ) )
1183  enddo
1184  enddo
1185  enddo
1186 #ifdef DEBUG
1187  i = iundef; j = iundef; k = iundef
1188 #endif
1189 
1190  enddo
1191  enddo
1192 
1193  return
integer, public is
start point of inner domain: x, local
integer, public i_xvz
integer, public je
end point of inner domain: y, local
real(rp), dimension(:), allocatable, public grid_rcdy
reciprocal of center-dy
real(rp), dimension(:), allocatable, public grid_fdy
y-length of grid(j+1) to grid(j) [m]
integer, public iblock
block size for cache blocking: x
real(rp), dimension(:), allocatable, public grid_rcdx
reciprocal of center-dx
integer, public i_xy
real(rp), dimension(:), allocatable, public grid_rfdy
reciprocal of face-dy
real(rp), dimension(:), allocatable, public grid_rcdz
reciprocal of center-dz
integer, public i_xvw
module GRIDTRANS
real(rp), dimension(:), allocatable, public grid_fdz
z-length of grid(k+1) to grid(k) [m]
integer, public i_uy
integer, public jblock
block size for cache blocking: y
integer, public i_xyw
integer, public js
start point of inner domain: y, local
integer, public i_uyw
module GRID (cartesian)
integer, public i_xv
integer, public ie
end point of inner domain: x, local
integer, public i_uyz
real(rp), dimension(:), allocatable, public grid_fdx
x-length of grid(i+1) to grid(i) [m]
real(rp), dimension(:), allocatable, public grid_rfdx
reciprocal of face-dx
real(rp), dimension(:), allocatable, public grid_rfdz
reciprocal of face-dz
integer, public i_uvz
integer, public i_xyz
integer, public i_uv
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_phy_tb_calc_flux_phi()

subroutine, public scale_atmos_phy_tb_common::atmos_phy_tb_calc_flux_phi ( real(rp), dimension(ka,ia,ja,3), intent(inout)  qflx_phi,
real(rp), dimension (ka,ia,ja), intent(in)  DENS,
real(rp), dimension (ka,ia,ja), intent(in)  PHI,
real(rp), dimension (ka,ia,ja), intent(in)  Kh,
real(rp), intent(in)  FACT,
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,4), intent(in)  MAPF,
real(rp), dimension (ka,ia,ja), intent(in)  a,
real(rp), dimension (ka,ia,ja), intent(in)  b,
real(rp), dimension (ka,ia,ja), intent(in)  c,
real(dp), intent(in)  dt,
logical, intent(in)  implicit,
integer, intent(in)  IIS,
integer, intent(in)  IIE,
integer, intent(in)  JJS,
integer, intent(in)  JJE 
)

Definition at line 1204 of file scale_atmos_phy_tb_common.F90.

References atmos_phy_tb_calc_tend_phi(), atmos_phy_tb_diffusion_solver(), scale_debug::check(), scale_grid::grid_fdz, scale_grid::grid_rfdx, scale_grid::grid_rfdy, scale_grid::grid_rfdz, scale_gridtrans::i_uy, scale_gridtrans::i_uyz, scale_gridtrans::i_xv, scale_gridtrans::i_xvz, scale_gridtrans::i_xyw, scale_gridtrans::i_xyz, 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_phy_tb_d1980::atmos_phy_tb_d1980(), and scale_atmos_phy_tb_smg::atmos_phy_tb_smg().

1204  use scale_gridtrans, only: &
1205  i_xyz, &
1206  i_xyw, &
1207  i_uyz, &
1208  i_xvz, &
1209  i_uy, &
1210  i_xv
1211  use scale_grid, only: &
1212  fdz => grid_fdz, &
1213  rfdz => grid_rfdz, &
1214  rfdx => grid_rfdx, &
1215  rfdy => grid_rfdy
1216  implicit none
1217 
1218  real(RP), intent(inout) :: qflx_phi(KA,IA,JA,3)
1219  real(RP), intent(in) :: DENS (KA,IA,JA)
1220  real(RP), intent(in) :: PHI (KA,IA,JA)
1221  real(RP), intent(in) :: Kh (KA,IA,JA)
1222  real(RP), intent(in) :: FACT
1223  real(RP), intent(in) :: GSQRT (KA,IA,JA,7)
1224  real(RP), intent(in) :: J13G (KA,IA,JA,7)
1225  real(RP), intent(in) :: J23G (KA,IA,JA,7)
1226  real(RP), intent(in) :: J33G
1227  real(RP), intent(in) :: MAPF (IA,JA,2,4)
1228  real(RP), intent(in) :: a (KA,IA,JA)
1229  real(RP), intent(in) :: b (KA,IA,JA)
1230  real(RP), intent(in) :: c (KA,IA,JA)
1231  real(DP), intent(in) :: dt
1232  logical, intent(in) :: implicit
1233  integer, intent(in) :: IIS
1234  integer, intent(in) :: IIE
1235  integer, intent(in) :: JJS
1236  integer, intent(in) :: JJE
1237 
1238  real(RP) :: TEND(KA,IA,JA)
1239  real(RP) :: d(KA)
1240 
1241  integer :: k, i, j
1242 
1243  ! (x-y plane; x,y,w)
1244  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
1245  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,DENS,Kh,PHI,qflx_phi,GSQRT,I_XYW,RFDZ,J33G) &
1246  !$omp shared(FDZ)
1247  do j = jjs, jje
1248  do i = iis, iie
1249  do k = ks, ke-1
1250 #ifdef DEBUG
1251  call check( __line__, dens(k,i,j) )
1252  call check( __line__, dens(k+1,i,j) )
1253  call check( __line__, kh(k,i,j) )
1254  call check( __line__, kh(k+1,i,j) )
1255  call check( __line__, phi(k+1,i,j) )
1256  call check( __line__, phi(k,i,j) )
1257  call check( __line__, rfdz(k) )
1258 #endif
1259  qflx_phi(k,i,j,zdir) = - 0.25_rp & ! 1/2/2
1260  * ( dens(k,i,j)+dens(k+1,i,j) ) &
1261  * ( kh(k,i,j) + kh(k+1,i,j) ) &
1262  * ( phi(k+1,i,j)-phi(k,i,j) ) * rfdz(k) * j33g &
1263  / gsqrt(k,i,j,i_xyw)
1264  enddo
1265  enddo
1266  enddo
1267 #ifdef DEBUG
1268  i = iundef; j = iundef; k = iundef
1269 #endif
1270  do j = jjs, jje
1271  do i = iis, iie
1272  qflx_phi(ks-1,i,j,zdir) = 0.0_rp
1273  qflx_phi(ke ,i,j,zdir) = 0.0_rp
1274  enddo
1275  enddo
1276 #ifdef DEBUG
1277  i = iundef; j = iundef; k = iundef
1278 #endif
1279 
1280  ! (y-z plane; u,y,z)
1281  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
1282  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,DENS,Kh,PHI,qflx_phi,GSQRT,I_XYZ,RFDX,J13G,I_UYZ) &
1283  !$omp shared(FDZ)
1284  do j = jjs, jje
1285  do i = iis-1, iie
1286  do k = ks+1, ke-1
1287 #ifdef DEBUG
1288  call check( __line__, dens(k,i,j) )
1289  call check( __line__, dens(k,i+1,j) )
1290  call check( __line__, kh(k,i,j) )
1291  call check( __line__, kh(k,i+1,j) )
1292  call check( __line__, phi(k,i+1,j) )
1293  call check( __line__, phi(k,i,j) )
1294  call check( __line__, rfdx(i) )
1295 #endif
1296  qflx_phi(k,i,j,xdir) = - 0.25_rp & ! 1/2/2
1297  * ( dens(k,i,j)+dens(k,i+1,j) ) &
1298  * ( kh(k,i,j) + kh(k,i+1,j) ) &
1299  * ( &
1300  ( gsqrt(k,i+1,j,i_xyz) * phi(k,i+1,j) &
1301  - gsqrt(k,i ,j,i_xyz) * phi(k,i ,j) ) * rfdx(i) &
1302  + ( j13g(k+1,i,j,i_uyz) * ( phi(k+1,i+1,j)+phi(k+1,i,j) ) &
1303  - j13g(k-1,i,j,i_uyz) * ( phi(k-1,i+1,j)+phi(k-1,i,j) ) &
1304  ) * 0.5_rp / ( fdz(k)+fdz(k-1) ) &
1305  ) / gsqrt(k,i,j,i_uyz)
1306  enddo
1307  enddo
1308  enddo
1309 #ifdef DEBUG
1310  i = iundef; j = iundef; k = iundef
1311 #endif
1312  do j = jjs, jje
1313  do i = iis-1, iie
1314 #ifdef DEBUG
1315  call check( __line__, dens(ks,i,j) )
1316  call check( __line__, dens(ks,i+1,j) )
1317  call check( __line__, kh(ks,i,j) )
1318  call check( __line__, kh(ks,i+1,j) )
1319  call check( __line__, phi(ks,i+1,j) )
1320  call check( __line__, phi(ks,i,j) )
1321  call check( __line__, rfdx(i) )
1322 #endif
1323  qflx_phi(ks,i,j,xdir) = - 0.25_rp & ! 1/2/2
1324  * ( dens(ks,i,j)+dens(ks,i+1,j) ) &
1325  * ( kh(ks,i,j) + kh(ks,i+1,j) ) &
1326  * ( &
1327  ( gsqrt(ks,i+1,j,i_xyz) * phi(ks,i+1,j) &
1328  - gsqrt(ks,i ,j,i_xyz) * phi(ks,i ,j) ) * rfdx(i) &
1329  + ( j13g(ks+1,i,j,i_uyz) * ( phi(ks+1,i+1,j)+phi(ks+1,i,j) ) &
1330  - j13g(ks ,i,j,i_uyz) * ( phi(ks ,i+1,j)+phi(ks ,i,j) ) &
1331  ) * 0.5_rp * rfdz(ks) &
1332  ) * mapf(i,j,1,i_uy) / gsqrt(ks,i,j,i_uyz)
1333  qflx_phi(ke,i,j,xdir) = - 0.25_rp & ! 1/2/2
1334  * ( dens(ke,i,j)+dens(ke,i+1,j) ) &
1335  * ( kh(ke,i,j) + kh(ke,i+1,j) ) &
1336  * ( &
1337  ( gsqrt(ke,i+1,j,i_xyz) * phi(ke,i+1,j) &
1338  - gsqrt(ke,i ,j,i_xyz) * phi(ke,i ,j) ) * rfdx(i) &
1339  + ( j13g(ke ,i,j,i_uyz) * ( phi(ke ,i+1,j)+phi(ke ,i,j) ) &
1340  - j13g(ke-1,i,j,i_uyz) * ( phi(ke-1,i+1,j)+phi(ke-1,i,j) ) &
1341  ) * 0.5_rp * rfdz(ke-1) &
1342  ) * mapf(i,j,1,i_uy) / gsqrt(ke,i,j,i_uyz)
1343  enddo
1344  enddo
1345 #ifdef DEBUG
1346  i = iundef; j = iundef; k = iundef
1347 #endif
1348  ! (z-x plane; x,v,z)
1349  !$omp parallel do default(none) &
1350  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,Kh,PHI,RFDY,DENS,qflx_phi,GSQRT,I_XYZ,J23G,I_XVZ,FDZ) &
1351  !$omp shared(MAPF,I_XV) &
1352  !$omp private(i,j,k) OMP_SCHEDULE_ collapse(2)
1353  do j = jjs-1, jje
1354  do i = iis, iie
1355  do k = ks+1, ke-1
1356 #ifdef DEBUG
1357  call check( __line__, dens(k,i,j) )
1358  call check( __line__, dens(k,i,j+1) )
1359  call check( __line__, kh(k,i,j) )
1360  call check( __line__, kh(k,i,j+1) )
1361  call check( __line__, phi(k,i,j+1) )
1362  call check( __line__, phi(k,i,j) )
1363  call check( __line__, rfdy(j) )
1364 #endif
1365  qflx_phi(k,i,j,ydir) = - 0.25_rp &
1366  * ( dens(k,i,j)+dens(k,i,j+1) ) &
1367  * ( kh(k,i,j) + kh(k,i,j+1) ) &
1368  * ( &
1369  ( gsqrt(k,i,j+1,i_xyz) * phi(k,i,j+1) &
1370  - gsqrt(k,i,j ,i_xyz) * phi(k,i,j ) ) * rfdy(j) &
1371  + ( j23g(k+1,i,j,i_xvz) * ( phi(k+1,i,j+1)+phi(k+1,i,j) ) &
1372  - j23g(k-1,i,j,i_xvz) * ( phi(k-1,i,j+1)+phi(k-1,i,j) ) &
1373  ) * 0.5_rp / ( fdz(k)+fdz(k-1) ) &
1374  ) * mapf(i,j,2,i_xv) / gsqrt(k,i,j,i_xvz)
1375  enddo
1376  enddo
1377  enddo
1378 #ifdef DEBUG
1379  i = iundef; j = iundef; k = iundef
1380 #endif
1381  do j = jjs-1, jje
1382  do i = iis, iie
1383 #ifdef DEBUG
1384  call check( __line__, dens(ks,i,j) )
1385  call check( __line__, dens(ks,i,j+1) )
1386  call check( __line__, kh(ks,i,j) )
1387  call check( __line__, kh(ks,i,j+1) )
1388  call check( __line__, phi(ks,i,j+1) )
1389  call check( __line__, phi(ks,i,j) )
1390  call check( __line__, rfdy(j) )
1391 #endif
1392  qflx_phi(ks,i,j,ydir) = - 0.25_rp &
1393  * ( dens(ks,i,j)+dens(ks,i,j+1) ) &
1394  * ( kh(ks,i,j) + kh(ks,i,j+1) ) &
1395  * ( &
1396  ( gsqrt(ks,i,j+1,i_xyz) * phi(ks,i,j+1) &
1397  - gsqrt(ks,i,j ,i_xyz) * phi(ks,i,j ) ) * rfdy(j) &
1398  + ( j23g(ks+1,i,j,i_xvz) * ( phi(ks+1,i,j+1)+phi(ks+1,i,j) ) &
1399  - j23g(ks ,i,j,i_xvz) * ( phi(ks ,i,j+1)+phi(ks ,i,j) ) &
1400  ) * 0.5_rp * rfdz(ks) &
1401  ) * mapf(i,j,2,i_xv) / gsqrt(ks,i,j,i_xvz)
1402  qflx_phi(ke,i,j,ydir) = - 0.25_rp &
1403  * ( dens(ke,i,j)+dens(ke,i,j+1) ) &
1404  * ( kh(ke,i,j) + kh(ke,i,j+1) ) &
1405  * ( &
1406  ( gsqrt(ke,i,j+1,i_xyz) * phi(ke,i,j+1) &
1407  - gsqrt(ke,i,j ,i_xyz) * phi(ke,i,j ) ) * rfdy(j) &
1408  + ( j23g(ke ,i,j,i_xvz) * ( phi(ke ,i,j+1)+phi(ke ,i,j) ) &
1409  - j23g(ke-1,i,j,i_xvz) * ( phi(ke-1,i,j+1)+phi(ke-1,i,j) ) &
1410  ) * 0.5_rp * rfdz(ke-1) &
1411  ) * mapf(i,j,2,i_xv) / gsqrt(ke,i,j,i_xvz)
1412  enddo
1413  enddo
1414 #ifdef DEBUG
1415  i = iundef; j = iundef; k = iundef
1416 #endif
1417 
1418  if ( implicit ) then
1419  call atmos_phy_tb_calc_tend_phi( tend, & ! (out)
1420  qflx_phi, & ! (in)
1421  gsqrt, j13g, j23g, j33g, mapf, & ! (in)
1422  iis, iie, jjs, jje ) ! (in)
1423 
1424  do j = jjs, jje
1425  do i = iis, iie
1426 
1427  do k = ks, ke
1428  d(k) = tend(k,i,j)
1429  end do
1430 
1431  call atmos_phy_tb_diffusion_solver( &
1432  tend(:,i,j), & ! (out)
1433  a(:,i,j), b(:,i,j), c(:,i,j), d, & ! (in)
1434  ke ) ! (in)
1435 
1436  do k = ks, ke-1
1437  qflx_phi(k,i,j,zdir) = qflx_phi(k,i,j,zdir) &
1438  - 0.25_rp & ! 1/2/2
1439  * ( dens(k,i,j)+dens(k+1,i,j) ) &
1440  * ( kh(k,i,j) + kh(k+1,i,j) ) &
1441  * dt * ( tend(k+1,i,j)-tend(k,i,j) ) * rfdz(k) * j33g &
1442  / gsqrt(k,i,j,i_xyw)
1443  end do
1444 
1445  end do
1446  end do
1447 
1448  end if
1449 
1450  return
integer, public i_xvz
integer, parameter, public zdir
integer, parameter, public ydir
integer, parameter, public xdir
real(rp), dimension(:), allocatable, public grid_rfdy
reciprocal of face-dy
subroutine, public atmos_phy_tb_calc_tend_phi(phi_t_TB, QFLX_phi, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
module GRIDTRANS
real(rp), dimension(:), allocatable, public grid_fdz
z-length of grid(k+1) to grid(k) [m]
integer, public i_uy
integer, public i_xyw
module GRID (cartesian)
integer, public i_xv
integer, public i_uyz
real(rp), dimension(:), allocatable, public grid_rfdx
reciprocal of face-dx
real(rp), dimension(:), allocatable, public grid_rfdz
reciprocal of face-dz
integer, public i_xyz
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_phy_tb_diffusion_solver()

subroutine, public scale_atmos_phy_tb_common::atmos_phy_tb_diffusion_solver ( real(rp), dimension(ka), intent(out)  phi,
real(rp), dimension(ka), intent(in)  a,
real(rp), dimension(ka), intent(in)  b,
real(rp), dimension(ka), intent(in)  c,
real(rp), dimension(ka), intent(in)  d,
integer, intent(in)  KE_TB 
)

Definition at line 1459 of file scale_atmos_phy_tb_common.F90.

References scale_grid_index::ka, and scale_grid_index::ks.

Referenced by atmos_phy_tb_calc_flux_phi(), scale_atmos_phy_tb_d1980::atmos_phy_tb_d1980(), scale_atmos_phy_tb_mynn::atmos_phy_tb_mynn(), and scale_atmos_phy_tb_smg::atmos_phy_tb_smg().

1459  implicit none
1460  real(RP), intent(out) :: phi(KA)
1461  real(RP), intent(in) :: a(KA)
1462  real(RP), intent(in) :: b(KA)
1463  real(RP), intent(in) :: c(KA)
1464  real(RP), intent(in) :: d(KA)
1465  integer, intent(in) :: KE_TB
1466  real(RP) :: e(KA)
1467  real(RP) :: f(KA)
1468  real(RP) :: denom
1469  integer :: k
1470 
1471  e(ks) = - a(ks) / b(ks)
1472  f(ks) = d(ks) / b(ks)
1473  do k = ks+1, ke_tb-1
1474  denom = b(k) + c(k)*e(k-1)
1475  e(k) = - a(k) / denom
1476  f(k) = ( d(k) - c(k)*f(k-1) ) / denom
1477  end do
1478 
1479  ! flux at the top boundary is zero
1480  phi(ke_tb) = ( d(ke_tb) - c(ke_tb)*f(ke_tb-1) ) / ( b(ke_tb) + c(ke_tb)*e(ke_tb-1) ) ! = f(KE_PBL)
1481 
1482  do k = ke_tb-1, ks, -1
1483  phi(k) = e(k) * phi(k+1) + f(k)
1484  end do
1485  do k = 1, ks-1
1486  phi(k) = 0.0_rp
1487  end do
1488  do k = ke_tb+1, ka
1489  phi(k) = 0.0_rp
1490  end do
1491 
1492  return
Here is the caller graph for this function:

◆ atmos_phy_tb_calc_tend_momz()

subroutine, public scale_atmos_phy_tb_common::atmos_phy_tb_calc_tend_momz ( real(rp), dimension(ka,ia,ja), intent(out)  MOMZ_t_TB,
real(rp), dimension(ka,ia,ja,3), intent(in)  QFLX_MOMZ,
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,4), intent(in)  MAPF,
integer, intent(in)  IIS,
integer, intent(in)  IIE,
integer, intent(in)  JJS,
integer, intent(in)  JJE 
)

Definition at line 1501 of file scale_atmos_phy_tb_common.F90.

References scale_grid::grid_cdz, scale_grid::grid_rcdx, scale_grid::grid_rcdy, scale_grid::grid_rfdz, scale_gridtrans::i_uyw, scale_gridtrans::i_xvw, scale_gridtrans::i_xy, scale_gridtrans::i_xyw, scale_gridtrans::i_xyz, 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_phy_tb_d1980::atmos_phy_tb_d1980(), mod_atmos_phy_tb_driver::atmos_phy_tb_driver(), and scale_atmos_phy_tb_smg::atmos_phy_tb_smg().

1501  use scale_grid, only: &
1502  rcdx => grid_rcdx, &
1503  rcdy => grid_rcdy, &
1504  rfdz => grid_rfdz, &
1505  cdz => grid_cdz
1506  use scale_gridtrans, only: &
1507  i_xyz, &
1508  i_xyw, &
1509  i_uyw, &
1510  i_xvw, &
1511  i_xy
1512  implicit none
1513 
1514  real(RP), intent(out) :: MOMZ_t_TB(KA,IA,JA)
1515  real(RP), intent(in) :: QFLX_MOMZ(KA,IA,JA,3)
1516  real(RP), intent(in) :: GSQRT(KA,IA,JA,7)
1517  real(RP), intent(in) :: J13G(KA,IA,JA,7)
1518  real(RP), intent(in) :: J23G(KA,IA,JA,7)
1519  real(RP), intent(in) :: J33G
1520  real(RP), intent(in) :: MAPF(IA,JA,2,4)
1521  integer , intent(in) :: IIS
1522  integer , intent(in) :: IIE
1523  integer , intent(in) :: JJS
1524  integer , intent(in) :: JJE
1525 
1526  integer :: k, i, j
1527 
1528  !$omp parallel do default(none) &
1529  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,MOMZ_t_TB,GSQRT,I_UYW,I_XVW,QFLX_MOMZ,RCDX,MAPF,I_XY) &
1530  !$omp shared(RCDY,J13G,I_XYZ,J23G,CDZ,J33G,RFDZ,I_XYW) &
1531  !$omp private(i,j,k) OMP_SCHEDULE_ collapse(2)
1532  do j = jjs, jje
1533  do i = iis, iie
1534  do k = ks+1, ke-2
1535  momz_t_tb(k,i,j) = &
1536  - ( ( gsqrt(k,i ,j,i_uyw) * qflx_momz(k,i ,j,xdir) &
1537  - gsqrt(k,i-1,j,i_uyw) * qflx_momz(k,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xy) &
1538  + ( gsqrt(k,i,j ,i_xvw) * qflx_momz(k,i,j ,ydir) &
1539  - gsqrt(k,i,j-1,i_xvw) * qflx_momz(k,i,j-1,ydir) ) * rcdy(j) &
1540  + ( ( j13g(k+1,i,j,i_xyz) * ( qflx_momz(k+1,i,j,xdir) + qflx_momz(k+1,i-1,j,xdir) ) &
1541  - j13g(k-1,i,j,i_xyz) * ( qflx_momz(k-1,i,j,xdir) + qflx_momz(k-1,i-1,j,xdir) ) &
1542  ) * mapf(i,j,1,i_xy) &
1543  + ( j23g(k+1,i,j,i_xyz) * ( qflx_momz(k+1,i,j,ydir) + qflx_momz(k+1,i,j-1,ydir) ) &
1544  - j23g(k-1,i,j,i_xyz) * ( qflx_momz(k-1,i,j,ydir) + qflx_momz(k-1,i,j-1,ydir) ) &
1545  ) * mapf(i,j,2,i_xy) &
1546  ) * 0.5_rp / ( cdz(k+1)+cdz(k) ) &
1547  + j33g * ( qflx_momz(k+1,i,j,zdir) - qflx_momz(k,i,j,zdir) ) * rfdz(k) ) &
1548  / gsqrt(k,i,j,i_xyw)
1549  enddo
1550  enddo
1551  enddo
1552 
1553  do j = jjs, jje
1554  do i = iis, iie
1555  momz_t_tb(ks,i,j) = &
1556  - ( ( gsqrt(ks,i ,j,i_uyw) * qflx_momz(ks,i ,j,xdir) &
1557  - gsqrt(ks,i-1,j,i_uyw) * qflx_momz(ks,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xy) &
1558  + ( gsqrt(ks,i,j ,i_xvw) * qflx_momz(ks,i,j ,ydir) &
1559  - gsqrt(ks,i,j-1,i_xvw) * qflx_momz(ks,i,j-1,ydir) ) * rcdy(j) * mapf(i,j,2,i_xy) &
1560  + ( ( ( qflx_momz(ks+1,i,j,xdir) + qflx_momz(ks+1,i-1,j,xdir) &
1561  + qflx_momz(ks ,i,j,xdir) + qflx_momz(ks ,i-1,j,xdir) &
1562  ) * j13g(ks+1,i,j,i_xyz) * mapf(i,j,1,i_xy) &
1563  + ( qflx_momz(ks+1,i,j,ydir) + qflx_momz(ks+1,i,j-1,ydir) &
1564  + qflx_momz(ks ,i,j,ydir) + qflx_momz(ks ,i,j-1,ydir) &
1565  ) * j23g(ks+1,i,j,i_xyz) * mapf(i,j,2,i_xy) &
1566  ) * 0.25_rp &
1567  + j33g * ( qflx_momz(ks+1,i,j,zdir) ) ) * rfdz(ks) &
1568  ) / gsqrt(ks,i,j,i_xyw)
1569 
1570  momz_t_tb(ke-1,i,j) = &
1571  - ( ( gsqrt(ke-1,i ,j,i_uyw) * qflx_momz(ke-1,i ,j,xdir) &
1572  - gsqrt(ke-1,i-1,j,i_uyw) * qflx_momz(ke-1,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xy) &
1573  + ( gsqrt(ke-1,i,j ,i_xvw) * qflx_momz(ke-1,i,j ,ydir) &
1574  - gsqrt(ke-1,i,j-1,i_xvw) * qflx_momz(ke-1,i,j-1,ydir) ) * rcdy(j) * mapf(i,j,2,i_xy) &
1575  + ( ( - ( qflx_momz(ke-1,i,j,xdir) + qflx_momz(ke-1,i-1,j,xdir) &
1576  + qflx_momz(ke-2,i,j,xdir) + qflx_momz(ke-2,i-1,j,xdir) &
1577  ) * j13g(ke-1,i,j,i_xyz) * mapf(i,j,1,i_xy) &
1578  - ( qflx_momz(ke-1,i,j,ydir) + qflx_momz(ke-1,i,j-1,ydir) &
1579  + qflx_momz(ke-2,i,j,ydir) + qflx_momz(ke-2,i,j-1,ydir) &
1580  ) * j23g(ke-1,i,j,i_xyz) * mapf(i,j,2,i_xy) &
1581  ) * 0.25_rp &
1582  - j33g * ( qflx_momz(ke-1,i,j,zdir) ) ) * rfdz(ke-1) &
1583  ) / gsqrt(ke-1,i,j,i_xyw)
1584  enddo
1585  enddo
1586 
1587  return
real(rp), dimension(:), allocatable, public grid_rcdy
reciprocal of center-dy
real(rp), dimension(:), allocatable, public grid_rcdx
reciprocal of center-dx
integer, parameter, public zdir
integer, parameter, public ydir
integer, parameter, public xdir
integer, public i_xy
integer, public i_xvw
module GRIDTRANS
integer, public i_xyw
integer, public i_uyw
module GRID (cartesian)
real(rp), dimension(:), allocatable, public grid_cdz
z-length of control volume [m]
real(rp), dimension(:), allocatable, public grid_rfdz
reciprocal of face-dz
integer, public i_xyz
Here is the caller graph for this function:

◆ atmos_phy_tb_calc_tend_momx()

subroutine, public scale_atmos_phy_tb_common::atmos_phy_tb_calc_tend_momx ( real(rp), dimension(ka,ia,ja), intent(out)  MOMX_t_TB,
real(rp), dimension(ka,ia,ja,3), intent(in)  QFLX_MOMX,
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,4), intent(in)  MAPF,
integer, intent(in)  IIS,
integer, intent(in)  IIE,
integer, intent(in)  JJS,
integer, intent(in)  JJE 
)

Definition at line 1595 of file scale_atmos_phy_tb_common.F90.

References scale_grid::grid_fdz, scale_grid::grid_rcdy, scale_grid::grid_rcdz, scale_grid::grid_rfdx, scale_grid::grid_rfdz, scale_gridtrans::i_uvz, scale_gridtrans::i_uy, scale_gridtrans::i_uyw, scale_gridtrans::i_uyz, scale_gridtrans::i_xyz, 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_phy_tb_d1980::atmos_phy_tb_d1980(), mod_atmos_phy_tb_driver::atmos_phy_tb_driver(), and scale_atmos_phy_tb_smg::atmos_phy_tb_smg().

1595  use scale_grid, only: &
1596  rcdz => grid_rcdz, &
1597  rcdy => grid_rcdy, &
1598  rfdz => grid_rfdz, &
1599  rfdx => grid_rfdx, &
1600  fdz => grid_fdz
1601  use scale_gridtrans, only: &
1602  i_xyz, &
1603  i_uyw, &
1604  i_uyz, &
1605  i_uvz, &
1606  i_uy
1607  implicit none
1608  real(RP), intent(out) :: MOMX_t_TB(KA,IA,JA)
1609  real(RP), intent(in) :: QFLX_MOMX(KA,IA,JA,3)
1610  real(RP), intent(in) :: GSQRT(KA,IA,JA,7)
1611  real(RP), intent(in) :: MAPF(IA,JA,2,4)
1612  real(RP), intent(in) :: J13G(KA,IA,JA,7)
1613  real(RP), intent(in) :: J23G(KA,IA,JA,7)
1614  real(RP), intent(in) :: J33G
1615  integer , intent(in) :: IIS
1616  integer , intent(in) :: IIE
1617  integer , intent(in) :: JJS
1618  integer , intent(in) :: JJE
1619 
1620  integer :: k, i, j
1621 
1622  !$omp parallel do default(none) &
1623  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,MOMX_t_TB,GSQRT,I_UVZ,I_XYZ,QFLX_MOMX,RFDX,MAPF,I_UY) &
1624  !$omp shared(RCDY,J13G,I_UYW,J23G,FDZ,J33G,RCDZ,I_UYZ) &
1625  !$omp private(i,j,k) OMP_SCHEDULE_ collapse(2)
1626  do j = jjs, jje
1627  do i = iis, iie
1628  do k = ks+1, ke-1
1629  momx_t_tb(k,i,j) = &
1630  - ( ( gsqrt(k,i+1,j,i_xyz) * qflx_momx(k,i+1,j,xdir) &
1631  - gsqrt(k,i ,j,i_xyz) * qflx_momx(k,i ,j,xdir) ) * rfdx(i) * mapf(i,j,1,i_uy) &
1632  + ( gsqrt(k,i,j ,i_uvz) * qflx_momx(k,i,j ,ydir) &
1633  - gsqrt(k,i,j-1,i_uvz) * qflx_momx(k,i,j-1,ydir) ) * rcdy(j) * mapf(i,j,2,i_uy) &
1634  + ( ( j13g(k+1,i,j,i_uyw) * ( qflx_momx(k+1,i+1,j,xdir) + qflx_momx(k+1,i,j ,xdir) ) &
1635  - j13g(k-1,i,j,i_uyw) * ( qflx_momx(k-1,i+1,j,xdir) + qflx_momx(k-1,i,j ,xdir) ) &
1636  ) * mapf(i,j,1,i_uy) &
1637  + ( j23g(k+1,i,j,i_uyw) * ( qflx_momx(k+1,i ,j,ydir) + qflx_momx(k+1,i,j-1,ydir) ) &
1638  - j23g(k-1,i,j,i_uyw) * ( qflx_momx(k-1,i ,j,ydir) + qflx_momx(k-1,i,j-1,ydir) ) &
1639  ) * mapf(i,j,2,i_uy) &
1640  ) * 0.5_rp / ( fdz(k)+fdz(k-1) ) &
1641  + j33g * ( qflx_momx(k,i,j,zdir) - qflx_momx(k-1,i,j,zdir) ) * rcdz(k) ) &
1642  / gsqrt(k,i,j,i_uyz)
1643  enddo
1644  enddo
1645  enddo
1646  do j = jjs, jje
1647  do i = iis, iie
1648  momx_t_tb(ks,i,j) = &
1649  - ( ( gsqrt(ks,i+1,j,i_xyz) * qflx_momx(ks,i+1,j,xdir) &
1650  - gsqrt(ks,i ,j,i_xyz) * qflx_momx(ks,i ,j,xdir) ) * rfdx(i) * mapf(i,j,1,i_uy) &
1651  + ( gsqrt(ks,i,j ,i_uvz) * qflx_momx(ks,i,j ,ydir) &
1652  - gsqrt(ks,i,j-1,i_uvz) * qflx_momx(ks,i,j-1,ydir) ) * rcdy(j) &
1653  + ( ( ( qflx_momx(ks+1,i+1,j,xdir) + qflx_momx(ks+1,i,j,xdir) &
1654  + qflx_momx(ks ,i+1,j,xdir) + qflx_momx(ks ,i,j,xdir) &
1655  ) * j13g(ks,i,j,i_uyw) * mapf(i,j,1,i_uy) &
1656  + ( qflx_momx(ks+1,i,j,ydir) + qflx_momx(ks+1,i,j-1,ydir) &
1657  + qflx_momx(ks ,i,j,ydir) + qflx_momx(ks ,i,j-1,ydir) &
1658  ) * j23g(ks,i,j,i_uyw) * mapf(i,j,2,i_uy) &
1659  ) * 0.25_rp &
1660  + j33g * ( qflx_momx(ks,i,j,zdir) ) ) * rfdz(ks) &
1661  ) / gsqrt(ks,i,j,i_uyz)
1662 
1663  momx_t_tb(ke,i,j) = &
1664  - ( ( gsqrt(ke,i+1,j,i_xyz) * qflx_momx(ke,i+1,j,xdir) &
1665  - gsqrt(ke,i ,j,i_xyz) * qflx_momx(ke,i ,j,xdir) ) * rfdx(i) * mapf(i,j,1,i_uy) &
1666  + ( gsqrt(ke,i,j ,i_uvz) * qflx_momx(ke,i,j ,ydir) &
1667  - gsqrt(ke,i,j-1,i_uvz) * qflx_momx(ke,i,j-1,ydir) ) * rcdy(j) * mapf(i,j,2,i_uy)&
1668  + ( ( - ( qflx_momx(ke ,i+1,j,xdir) + qflx_momx(ke ,i,j,xdir) &
1669  + qflx_momx(ke-1,i+1,j,xdir) + qflx_momx(ke-1,i,j,xdir) &
1670  ) * j13g(ke-1,i,j,i_uyw) * mapf(i,j,1,i_uy) &
1671  - ( qflx_momx(ke ,i,j,ydir) + qflx_momx(ke ,i,j-1,ydir) &
1672  + qflx_momx(ke-1,i,j,ydir) + qflx_momx(ke-1,i,j-1,ydir) &
1673  ) * j23g(ke-1,i,j,i_uyw) * mapf(i,j,2,i_uy) &
1674  ) * 0.25_rp &
1675  - j33g * ( qflx_momx(ke-1,i,j,zdir) ) ) * rcdz(ke) &
1676  ) / gsqrt(ke,i,j,i_uyz)
1677  enddo
1678  enddo
1679 
1680  return
real(rp), dimension(:), allocatable, public grid_rcdy
reciprocal of center-dy
integer, parameter, public zdir
integer, parameter, public ydir
integer, parameter, public xdir
real(rp), dimension(:), allocatable, public grid_rcdz
reciprocal of center-dz
module GRIDTRANS
real(rp), dimension(:), allocatable, public grid_fdz
z-length of grid(k+1) to grid(k) [m]
integer, public i_uy
integer, public i_uyw
module GRID (cartesian)
integer, public i_uyz
real(rp), dimension(:), allocatable, public grid_rfdx
reciprocal of face-dx
real(rp), dimension(:), allocatable, public grid_rfdz
reciprocal of face-dz
integer, public i_uvz
integer, public i_xyz
Here is the caller graph for this function:

◆ atmos_phy_tb_calc_tend_momy()

subroutine, public scale_atmos_phy_tb_common::atmos_phy_tb_calc_tend_momy ( real(rp), dimension(ka,ia,ja), intent(out)  MOMY_t_TB,
real(rp), dimension(ka,ia,ja,3), intent(in)  QFLX_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,4), intent(in)  MAPF,
integer, intent(in)  IIS,
integer, intent(in)  IIE,
integer, intent(in)  JJS,
integer, intent(in)  JJE 
)

Definition at line 1688 of file scale_atmos_phy_tb_common.F90.

References scale_grid::grid_fdz, scale_grid::grid_rcdx, scale_grid::grid_rcdz, scale_grid::grid_rfdy, scale_gridtrans::i_uvz, scale_gridtrans::i_xv, scale_gridtrans::i_xvw, scale_gridtrans::i_xyz, 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_phy_tb_d1980::atmos_phy_tb_d1980(), mod_atmos_phy_tb_driver::atmos_phy_tb_driver(), and scale_atmos_phy_tb_smg::atmos_phy_tb_smg().

1688  use scale_grid, only: &
1689  rcdz => grid_rcdz, &
1690  rcdx => grid_rcdx, &
1691  rfdy => grid_rfdy, &
1692  fdz => grid_fdz
1693  use scale_gridtrans, only: &
1694  i_xyz, &
1695  i_xvw, &
1696  i_uvz, &
1697  i_xv
1698  implicit none
1699 
1700  real(RP), intent(out) :: MOMY_t_TB(KA,IA,JA)
1701  real(RP), intent(in) :: QFLX_MOMY(KA,IA,JA,3)
1702  real(RP), intent(in) :: GSQRT(KA,IA,JA,7)
1703  real(RP), intent(in) :: MAPF(IA,JA,2,4)
1704  real(RP), intent(in) :: J13G(KA,IA,JA,7)
1705  real(RP), intent(in) :: J23G(KA,IA,JA,7)
1706  real(RP), intent(in) :: J33G
1707  integer , intent(in) :: IIS
1708  integer , intent(in) :: IIE
1709  integer , intent(in) :: JJS
1710  integer , intent(in) :: JJE
1711 
1712  integer :: k, i, j
1713 
1714  !$omp parallel do default(none) &
1715  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,MOMY_t_TB,GSQRT,I_UVZ,I_XYZ,QFLX_MOMY,RCDX,MAPF,I_XV) &
1716  !$omp shared(RFDY,J13G,I_XVW,J23G,FDZ,J33G,RCDZ) &
1717  !$omp private(i,j,k) OMP_SCHEDULE_ collapse(2)
1718  do j = jjs, jje
1719  do i = iis, iie
1720  do k = ks+1, ke-1
1721  momy_t_tb(k,i,j) = &
1722  - ( ( gsqrt(k,i ,j ,i_uvz) * qflx_momy(k,i ,j,xdir) &
1723  - gsqrt(k,i-1,j ,i_uvz) * qflx_momy(k,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xv) &
1724  + ( gsqrt(k,i ,j+1,i_xyz) * qflx_momy(k,i,j+1,ydir) &
1725  - gsqrt(k,i ,j ,i_xyz) * qflx_momy(k,i,j ,ydir) ) * rfdy(j) * mapf(i,j,2,i_xv) &
1726  + ( ( j13g(k+1,i,j,i_xvw) * ( qflx_momy(k+1,i,j ,xdir) + qflx_momy(k+1,i-1,j,xdir) ) &
1727  - j13g(k-1,i,j,i_xvw) * ( qflx_momy(k-1,i,j ,xdir) + qflx_momy(k-1,i-1,j,xdir) ) &
1728  ) * mapf(i,j,1,i_xv) &
1729  + ( j23g(k+1,i,j,i_xvw) * ( qflx_momy(k+1,i,j+1,ydir) + qflx_momy(k+1,i ,j,ydir) ) &
1730  - j23g(k-1,i,j,i_xvw) * ( qflx_momy(k-1,i,j+1,ydir) + qflx_momy(k-1,i ,j,ydir) ) &
1731  ) * mapf(i,j,2,i_xv) &
1732  ) * 0.5_rp / ( fdz(k)+fdz(k-1) ) &
1733  + j33g * ( qflx_momy(k,i,j,zdir) - qflx_momy(k-1,i,j,zdir) ) * rcdz(k) ) &
1734  / gsqrt(k,i,j,i_xvw)
1735  enddo
1736  enddo
1737  enddo
1738  do j = jjs, jje
1739  do i = iis, iie
1740  momy_t_tb(ks,i,j) = &
1741  - ( ( gsqrt(ks,i ,j ,i_uvz) * qflx_momy(ks,i ,j,xdir) &
1742  - gsqrt(ks,i-1,j ,i_uvz) * qflx_momy(ks,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xv) &
1743  + ( gsqrt(ks,i ,j+1,i_xyz) * qflx_momy(ks,i,j+1,ydir) &
1744  - gsqrt(ks,i ,j ,i_xyz) * qflx_momy(ks,i,j ,ydir) ) * rfdy(j) * mapf(i,j,2,i_xv) &
1745  + ( ( ( qflx_momy(ks+1,i,j ,xdir) + qflx_momy(ks+1,i-1,j,xdir) &
1746  + qflx_momy(ks ,i,j ,xdir) + qflx_momy(ks ,i-1,j,xdir) &
1747  ) * j13g(ks,i,j,i_xvw) * mapf(i,j,1,i_xv) &
1748  + ( qflx_momy(ks+1,i,j+1,ydir) + qflx_momy(ks+1,i ,j,ydir) &
1749  + qflx_momy(ks ,i,j+1,ydir) + qflx_momy(ks ,i ,j,ydir) &
1750  ) * j23g(ks,i,j,i_xvw) * mapf(i,j,2,i_xv) &
1751  ) * 0.25_rp &
1752  + j33g * ( qflx_momy(ks,i,j,zdir) ) ) * rcdz(ks) &
1753  ) / gsqrt(ks,i,j,i_xvw)
1754 
1755  momy_t_tb(ke,i,j) = &
1756  - ( ( gsqrt(ke,i ,j ,i_uvz) * qflx_momy(ke,i ,j,xdir) &
1757  - gsqrt(ke,i-1,j ,i_uvz) * qflx_momy(ke,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xv) &
1758  + ( gsqrt(ke,i ,j+1,i_xyz) * qflx_momy(ke,i,j+1,ydir) &
1759  - gsqrt(ke,i ,j ,i_xyz) * qflx_momy(ke,i,j ,ydir) ) * rfdy(j) * mapf(i,j,2,i_xv) &
1760  + ( ( - ( qflx_momy(ke ,i,j,xdir) + qflx_momy(ke ,i-1,j,xdir) &
1761  + qflx_momy(ke-1,i,j,xdir) + qflx_momy(ke-1,i-1,j,xdir) &
1762  ) * j13g(ke-1,i,j,i_xvw) * mapf(i,j,1,i_xv) &
1763  - ( qflx_momy(ke ,i,j+1,ydir) + qflx_momy(ke ,i,j,ydir) &
1764  + qflx_momy(ke-1,i,j+1,ydir) + qflx_momy(ke-1,i,j,ydir) &
1765  ) * j23g(ke-1,i,j,i_xvw) * mapf(i,j,2,i_xv) &
1766  ) * 0.25_rp &
1767  - j33g * ( qflx_momy(ke-1,i,j,zdir) ) ) * rcdz(ke) &
1768  ) / gsqrt(ke,i,j,i_xvw)
1769  end do
1770  end do
1771 
1772  return
real(rp), dimension(:), allocatable, public grid_rcdx
reciprocal of center-dx
integer, parameter, public zdir
integer, parameter, public ydir
integer, parameter, public xdir
real(rp), dimension(:), allocatable, public grid_rfdy
reciprocal of face-dy
real(rp), dimension(:), allocatable, public grid_rcdz
reciprocal of center-dz
integer, public i_xvw
module GRIDTRANS
real(rp), dimension(:), allocatable, public grid_fdz
z-length of grid(k+1) to grid(k) [m]
module GRID (cartesian)
integer, public i_xv
integer, public i_uvz
integer, public i_xyz
Here is the caller graph for this function:

◆ atmos_phy_tb_calc_tend_phi()

subroutine, public scale_atmos_phy_tb_common::atmos_phy_tb_calc_tend_phi ( real(rp), dimension(ka,ia,ja), intent(out)  phi_t_TB,
real(rp), dimension(ka,ia,ja,3), intent(in)  QFLX_phi,
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,4), intent(in)  MAPF,
integer, intent(in)  IIS,
integer, intent(in)  IIE,
integer, intent(in)  JJS,
integer, intent(in)  JJE 
)

Definition at line 1780 of file scale_atmos_phy_tb_common.F90.

References scale_grid::grid_fdz, scale_grid::grid_rcdx, scale_grid::grid_rcdy, scale_grid::grid_rcdz, scale_gridtrans::i_uvz, scale_gridtrans::i_uyz, scale_gridtrans::i_xvz, scale_gridtrans::i_xy, scale_gridtrans::i_xyw, scale_gridtrans::i_xyz, scale_grid_index::ke, scale_grid_index::ks, scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by atmos_phy_tb_calc_flux_phi(), scale_atmos_phy_tb_d1980::atmos_phy_tb_d1980(), mod_atmos_phy_tb_driver::atmos_phy_tb_driver(), and scale_atmos_phy_tb_smg::atmos_phy_tb_smg().

1780  use scale_grid, only: &
1781  rcdz => grid_rcdz, &
1782  rcdx => grid_rcdx, &
1783  rcdy => grid_rcdy, &
1784  fdz => grid_fdz
1785  use scale_gridtrans, only: &
1786  i_xyz, &
1787  i_xyw, &
1788  i_uyz, &
1789  i_xvz, &
1790  i_uvz, &
1791  i_xy
1792  implicit none
1793 
1794  real(RP), intent(out) :: phi_t_TB(KA,IA,JA)
1795  real(RP), intent(in) :: QFLX_phi(KA,IA,JA,3)
1796  real(RP), intent(in) :: GSQRT(KA,IA,JA,7)
1797  real(RP), intent(in) :: MAPF(IA,JA,2,4)
1798  real(RP), intent(in) :: J13G(KA,IA,JA,7)
1799  real(RP), intent(in) :: J23G(KA,IA,JA,7)
1800  real(RP), intent(in) :: J33G
1801  integer , intent(in) :: IIS
1802  integer , intent(in) :: IIE
1803  integer , intent(in) :: JJS
1804  integer , intent(in) :: JJE
1805 
1806  integer :: k, i, j
1807 
1808  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
1809  !$omp shared(JJS,JJE,IIS,IIE,KS,KE,phi_t_TB,GSQRT,I_UYZ,QFLX_phi,I_UVZ,RCDX,MAPF,I_XY) &
1810  !$omp shared(I_XVZ,J13G,I_XYW,J23G,FDZ,J33G,RCDZ,I_XYZ,RCDY)
1811  do j = jjs, jje
1812  do i = iis, iie
1813  do k = ks+1, ke-1
1814  phi_t_tb(k,i,j) = &
1815  - ( ( gsqrt(k,i ,j,i_uyz) * qflx_phi(k,i ,j,xdir) &
1816  - gsqrt(k,i-1,j,i_uvz) * qflx_phi(k,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xy) &
1817  + ( gsqrt(k,i,j ,i_xvz) * qflx_phi(k,i,j ,ydir) &
1818  - gsqrt(k,i,j-1,i_xvz) * qflx_phi(k,i,j-1,ydir) ) * rcdy(j) * mapf(i,j,2,i_xy) &
1819  + ( ( j13g(k+1,i,j,i_xyw) * ( qflx_phi(k+1,i,j,xdir) + qflx_phi(k+1,i-1,j,xdir) ) &
1820  - j13g(k-1,i,j,i_xyw) * ( qflx_phi(k-1,i,j,xdir) + qflx_phi(k-1,i-1,j,xdir) ) &
1821  ) * mapf(i,j,1,i_xy) &
1822  + ( j23g(k+1,i,j,i_xyw) * ( qflx_phi(k+1,i,j,ydir) + qflx_phi(k+1,i,j-1,ydir) ) &
1823  - j23g(k-1,i,j,i_xyw) * ( qflx_phi(k-1,i,j,ydir) + qflx_phi(k-1,i,j-1,ydir) ) &
1824  ) * mapf(i,j,2,i_xy) &
1825  ) * 0.5_rp / ( fdz(k) + fdz(k-1) ) &
1826  + j33g * ( qflx_phi(k,i,j,zdir) - qflx_phi(k-1,i,j,zdir) ) * rcdz(k) ) &
1827  / gsqrt(k,i,j,i_xyz)
1828  enddo
1829  enddo
1830  enddo
1831  do j = jjs, jje
1832  do i = iis, iie
1833  phi_t_tb(ks,i,j) = &
1834  - ( ( gsqrt(ks,i ,j,i_uyz) * qflx_phi(ks,i ,j,xdir) &
1835  - gsqrt(ks,i-1,j,i_uvz) * qflx_phi(ks,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xy) &
1836  + ( gsqrt(ks,i,j ,i_xvz) * qflx_phi(ks,i,j ,ydir) &
1837  - gsqrt(ks,i,j-1,i_xvz) * qflx_phi(ks,i,j-1,ydir) ) * rcdy(j) * mapf(i,j,2,i_xy) &
1838  + ( ( ( qflx_phi(ks+1,i,j,xdir) + qflx_phi(ks+1,i-1,j,xdir) &
1839  + qflx_phi(ks ,i,j,xdir) + qflx_phi(ks ,i-1,j,xdir) &
1840  ) * j13g(ks+1,i,j,i_xyw) * mapf(i,j,1,i_xy) &
1841  + ( qflx_phi(ks+1,i,j,ydir) + qflx_phi(ks+1,i,j-1,ydir) &
1842  + qflx_phi(ks ,i,j,ydir) + qflx_phi(ks ,i,j-1,ydir) &
1843  ) * j23g(ks+1,i,j,i_xyw) * mapf(i,j,2,i_xy) &
1844  ) * 0.25_rp &
1845  + j33g * ( qflx_phi(ks,i,j,zdir) ) ) * rcdz(ks) &
1846  ) / gsqrt(ks,i,j,i_xyz)
1847 
1848  phi_t_tb(ke,i,j) = &
1849  - ( ( gsqrt(ke,i ,j,i_uyz) * qflx_phi(ke,i ,j,xdir) &
1850  - gsqrt(ke,i-1,j,i_uvz) * qflx_phi(ke,i-1,j,xdir) ) * rcdx(i) * mapf(i,j,1,i_xy) &
1851  + ( gsqrt(ke,i,j ,i_xvz) * qflx_phi(ke,i,j ,ydir) &
1852  - gsqrt(ke,i,j-1,i_xvz) * qflx_phi(ke,i,j-1,ydir) ) * rcdy(j) * mapf(i,j,2,i_xy) &
1853  + ( ( - ( qflx_phi(ke ,i,j,xdir) + qflx_phi(ke ,i-1,j,xdir) &
1854  + qflx_phi(ke-1,i,j,xdir) + qflx_phi(ke-1,i-1,j,xdir) &
1855  ) * j13g(ke-1,i,j,i_xyw) * mapf(i,j,1,i_xy) &
1856  - ( qflx_phi(ke ,i,j,ydir) + qflx_phi(ke ,i,j-1,ydir) &
1857  + qflx_phi(ke-1,i,j,ydir) + qflx_phi(ke-1,i,j-1,ydir) &
1858  ) * j23g(ke-1,i,j,i_xyw) * mapf(i,j,2,i_xy) &
1859  ) * 0.25_rp &
1860  - j33g * ( qflx_phi(ke-1,i,j,zdir) ) ) * rcdz(ke) &
1861  ) / gsqrt(ke,i,j,i_xyz)
1862  end do
1863  end do
1864 
1865  return
integer, public i_xvz
real(rp), dimension(:), allocatable, public grid_rcdy
reciprocal of center-dy
real(rp), dimension(:), allocatable, public grid_rcdx
reciprocal of center-dx
integer, parameter, public zdir
integer, parameter, public ydir
integer, parameter, public xdir
integer, public i_xy
real(rp), dimension(:), allocatable, public grid_rcdz
reciprocal of center-dz
module GRIDTRANS
real(rp), dimension(:), allocatable, public grid_fdz
z-length of grid(k+1) to grid(k) [m]
integer, public i_xyw
module GRID (cartesian)
integer, public i_uyz
integer, public i_uvz
integer, public i_xyz
Here is the caller graph for this function: