SCALE-RM
Functions/Subroutines
scale_atmos_phy_tb_smg Module Reference

module ATMOSPHERE / Physics Turbulence More...

Functions/Subroutines

subroutine, public atmos_phy_tb_smg_setup (TYPE_TB, CDZ, CDX, CDY, CZ)
 
subroutine, public atmos_phy_tb_smg (qflx_sgs_momz, qflx_sgs_momx, qflx_sgs_momy, qflx_sgs_rhot, qflx_sgs_rhoq, tke, tke_t, nu, Ri, Pr, N2, MOMZ, MOMX, MOMY, RHOT, DENS, QTRC, SFLX_MW, SFLX_MU, SFLX_MV, SFLX_SH, SFLX_QV, GSQRT, J13G, J23G, J33G, MAPF, dt)
 
real(rp) function mixlen (dz, dx, dy, z, filter_fact)
 
real(rp) function fact (dz, dx, dy)
 

Detailed Description

module ATMOSPHERE / Physics Turbulence

Description
Sub-grid scale turbulence process Smagolinsky-type
Author
Team SCALE
History
  • 2011-11-29 (S.Iga) [new]
  • 2011-12-11 (H.Yashiro) [mod] integrate to SCALE-LES ver.3
  • 2012-03-23 (H.Yashiro) [mod] Explicit index parameter inclusion
  • 2012-03-27 (H.Yashiro) [mod] reconstruction
  • 2012-07-02 (S.Nishizawa) [mod] reconstruction with Brown et al. (1994)
  • 2012-10-26 (S.Nishizawa) [mod] remove surface flux
  • 2013-06-13 (S.Nishizawa) [mod] change mixing length by Brown et al. (1994) and Scotti et al. (1993)
  • 2014-04-02 (S.Nishizawa) [mod] modified for terrrain-following
  • Reference
    • Brown et al., 1994: Large-eddy simulaition of stable atmospheric boundary layers with a revised stochastic subgrid model. Roy. Meteor. Soc., 120, 1485-1512
    • Scotti et al., 1993: Generalized Smagorinsky model for anisotropic grids. Phys. Fluids A, 5, 2306-2308
    • Sullivan et al., 1994: A subgrid-scale model for large-eddy simulation of planetary boundary-layer flows. Bound.-Layer Meteor., 71, 247-276
NAMELIST
  • PARAM_ATMOS_PHY_TB_SMG
    nametypedefault valuecomment
    ATMOS_PHY_TB_SMG_CS real(RP)
    ATMOS_PHY_TB_SMG_NU_MAX real(RP) 10000.0_RP
    ATMOS_PHY_TB_SMG_FILTER_FACT real(RP) 2.0_RP
    ATMOS_PHY_TB_SMG_IMPLICIT logical .false.
    ATMOS_PHY_TB_SMG_CONSISTENT_TKE logical .true.
    ATMOS_PHY_TB_SMG_BOTTOM logical .false.

History Output
No history output

Function/Subroutine Documentation

◆ atmos_phy_tb_smg_setup()

subroutine, public scale_atmos_phy_tb_smg::atmos_phy_tb_smg_setup ( character(len=*), intent(in)  TYPE_TB,
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,ia,ja), intent(in)  CZ 
)

Definition at line 102 of file scale_atmos_phy_tb_smg.F90.

References scale_grid_index::ia, scale_grid_index::ie, scale_stdio::io_fid_conf, scale_stdio::io_fid_log, scale_stdio::io_l, scale_stdio::io_lnml, scale_grid_index::is, scale_grid_index::ja, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ka, scale_grid_index::ke, scale_grid_index::ks, mixlen(), and scale_process::prc_mpistop().

Referenced by scale_atmos_phy_tb_hybrid::atmos_phy_tb_hybrid_setup(), and scale_atmos_phy_tb::atmos_phy_tb_setup().

102  use scale_process, only: &
104  implicit none
105 
106  character(len=*), intent(in) :: type_tb
107 
108  real(RP), intent(in) :: cdz(ka)
109  real(RP), intent(in) :: cdx(ia)
110  real(RP), intent(in) :: cdy(ja)
111  real(RP), intent(in) :: cz (ka,ia,ja)
112 
113  real(RP) :: atmos_phy_tb_smg_cs
114  real(RP) :: atmos_phy_tb_smg_filter_fact = 2.0_rp
115  logical :: atmos_phy_tb_smg_consistent_tke = .true.
116 
117  namelist / param_atmos_phy_tb_smg / &
118  atmos_phy_tb_smg_cs, &
119  atmos_phy_tb_smg_nu_max, &
120  atmos_phy_tb_smg_filter_fact, &
121  atmos_phy_tb_smg_implicit, &
122  atmos_phy_tb_smg_consistent_tke, &
123  atmos_phy_tb_smg_bottom
124 
125  integer :: k, i, j
126 
127  integer :: ierr
128  !---------------------------------------------------------------------------
129 
130  if( io_l ) write(io_fid_log,*)
131  if( io_l ) write(io_fid_log,*) '++++++ Module[TURBULENCE] / Categ[ATMOS PHYSICS] / Origin[SCALElib]'
132  if( io_l ) write(io_fid_log,*) '+++ Smagorinsky-type Eddy Viscocity Model'
133 
134  if ( type_tb /= 'SMAGORINSKY' ) then
135  write(*,*) 'xxx ATMOS_PHY_TB_TYPE is not SMAGORINSKY. Check!'
136  call prc_mpistop
137  endif
138 
139  atmos_phy_tb_smg_cs = cs
140 
141  !--- read namelist
142  rewind(io_fid_conf)
143  read(io_fid_conf,nml=param_atmos_phy_tb_smg,iostat=ierr)
144  if( ierr < 0 ) then !--- missing
145  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
146  elseif( ierr > 0 ) then !--- fatal error
147  write(*,*) 'xxx Not appropriate names in namelist PARAM_ATMOS_PHY_TB_SMG. Check!'
148  call prc_mpistop
149  endif
150  if( io_lnml ) write(io_fid_log,nml=param_atmos_phy_tb_smg)
151 
152  cs = atmos_phy_tb_smg_cs
153 
154  rprn = 1.0_rp / prn
155  rric = 1.0_rp / ric
156  prnovric = (1- prn) * rric
157 
158  allocate( nu_fact(ka,ia,ja) )
159 
160 #ifdef DEBUG
161  nu_fact(:,:,:) = undef
162 #endif
163  do j = js-1, je+1
164  do i = is-1, ie+1
165  do k = ks, ke
166  nu_fact(k,i,j) = ( cs * mixlen(cdz(k),cdx(i),cdy(j),cz(k,i,j),atmos_phy_tb_smg_filter_fact) )**2
167  enddo
168  enddo
169  enddo
170 #ifdef DEBUG
171  i = iundef; j = iundef; k = iundef
172 #endif
173 
174  if ( atmos_phy_tb_smg_consistent_tke ) then
175  tke_fact = 1.0_rp
176  else
177  tke_fact = 0.0_rp
178  end if
179 
180 
181  return
integer, public is
start point of inner domain: x, local
integer, public je
end point of inner domain: y, local
subroutine, public prc_mpistop
Abort MPI.
logical, public io_l
output log or not? (this process)
Definition: scale_stdio.F90:59
integer, public ke
end point of inner domain: z, local
integer, public ia
of x whole cells (local, with HALO)
integer, public ka
of z whole cells (local, with HALO)
integer, public js
start point of inner domain: y, local
module PROCESS
integer, public ks
start point of inner domain: z, local
integer, public ie
end point of inner domain: x, local
logical, public io_lnml
output log or not? (for namelist, this process)
Definition: scale_stdio.F90:60
integer, public io_fid_conf
Config file ID.
Definition: scale_stdio.F90:55
integer, public io_fid_log
Log file ID.
Definition: scale_stdio.F90:56
integer, public ja
of y whole cells (local, with HALO)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_phy_tb_smg()

subroutine, public scale_atmos_phy_tb_smg::atmos_phy_tb_smg ( real(rp), dimension(ka,ia,ja,3), intent(out)  qflx_sgs_momz,
real(rp), dimension(ka,ia,ja,3), intent(out)  qflx_sgs_momx,
real(rp), dimension(ka,ia,ja,3), intent(out)  qflx_sgs_momy,
real(rp), dimension(ka,ia,ja,3), intent(out)  qflx_sgs_rhot,
real(rp), dimension(ka,ia,ja,3,qa), intent(out)  qflx_sgs_rhoq,
real(rp), dimension(ka,ia,ja), intent(inout)  tke,
real(rp), dimension(ka,ia,ja), intent(out)  tke_t,
real(rp), dimension (ka,ia,ja), intent(out)  nu,
real(rp), dimension (ka,ia,ja), intent(out)  Ri,
real(rp), dimension (ka,ia,ja), intent(out)  Pr,
real(rp), dimension (ka,ia,ja), intent(out)  N2,
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,ia,ja), intent(in)  DENS,
real(rp), dimension(ka,ia,ja,qa), intent(in)  QTRC,
real(rp), dimension(ia,ja), intent(in)  SFLX_MW,
real(rp), dimension(ia,ja), intent(in)  SFLX_MU,
real(rp), dimension(ia,ja), intent(in)  SFLX_MV,
real(rp), dimension(ia,ja), intent(in)  SFLX_SH,
real(rp), dimension(ia,ja), intent(in)  SFLX_QV,
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(dp), intent(in)  dt 
)
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 192 of file scale_atmos_phy_tb_smg.F90.

References scale_atmos_phy_tb_common::atmos_phy_tb_calc_flux_phi(), scale_atmos_phy_tb_common::atmos_phy_tb_calc_strain_tensor(), scale_atmos_phy_tb_common::atmos_phy_tb_calc_tend_momx(), scale_atmos_phy_tb_common::atmos_phy_tb_calc_tend_momy(), scale_atmos_phy_tb_common::atmos_phy_tb_calc_tend_momz(), scale_atmos_phy_tb_common::atmos_phy_tb_calc_tend_phi(), scale_atmos_phy_tb_common::atmos_phy_tb_diffusion_solver(), scale_debug::check(), scale_const::const_eps, scale_const::const_grav, 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_stdio::io_fid_log, scale_stdio::io_l, 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_tracer::qa, scale_grid_index::xdir, scale_grid_index::ydir, and scale_grid_index::zdir.

Referenced by scale_atmos_phy_tb_hybrid::atmos_phy_tb_hybrid_setup(), and scale_atmos_phy_tb::atmos_phy_tb_setup().

192  use scale_grid_index
193  use scale_tracer
194  use scale_const, only: &
195  eps => const_eps, &
196  grav => const_grav
197  use scale_grid, only: &
198  fdz => grid_fdz, &
199  fdx => grid_fdx, &
200  fdy => grid_fdy, &
201  rcdz => grid_rcdz, &
202  rcdx => grid_rcdx, &
203  rcdy => grid_rcdy, &
204  rfdz => grid_rfdz, &
205  rfdx => grid_rfdx, &
206  rfdy => grid_rfdy
207  use scale_gridtrans, only: &
208  i_xyz, &
209  i_xyw, &
210  i_uyw, &
211  i_xvw, &
212  i_uyz, &
213  i_xvz, &
214  i_uvz, &
215  i_xy, &
216  i_uy, &
217  i_xv, &
218  i_uv
219  use scale_atmos_phy_tb_common, only: &
220  calc_strain_tensor => atmos_phy_tb_calc_strain_tensor, &
221  diffusion_solver => atmos_phy_tb_diffusion_solver, &
222  calc_tend_momz => atmos_phy_tb_calc_tend_momz, &
223  calc_tend_momx => atmos_phy_tb_calc_tend_momx, &
224  calc_tend_momy => atmos_phy_tb_calc_tend_momy, &
225  calc_tend_phi => atmos_phy_tb_calc_tend_phi, &
226  calc_flux_phi => atmos_phy_tb_calc_flux_phi
227  implicit none
228 
229  ! SGS flux
230  real(RP), intent(out) :: qflx_sgs_momz(ka,ia,ja,3)
231  real(RP), intent(out) :: qflx_sgs_momx(ka,ia,ja,3)
232  real(RP), intent(out) :: qflx_sgs_momy(ka,ia,ja,3)
233  real(RP), intent(out) :: qflx_sgs_rhot(ka,ia,ja,3)
234  real(RP), intent(out) :: qflx_sgs_rhoq(ka,ia,ja,3,qa)
235 
236  real(RP), intent(inout) :: tke(ka,ia,ja) ! TKE
237  real(RP), intent(out) :: tke_t(ka,ia,ja) ! tendency TKE
238  real(RP), intent(out) :: nu (ka,ia,ja) ! eddy viscosity (center)
239  real(RP), intent(out) :: ri (ka,ia,ja) ! Richardson number
240  real(RP), intent(out) :: pr (ka,ia,ja) ! Prantle number
241  real(RP), intent(out) :: n2 (ka,ia,ja) ! squared Brunt-Vaisala frequency
242 
243  real(RP), intent(in) :: momz(ka,ia,ja)
244  real(RP), intent(in) :: momx(ka,ia,ja)
245  real(RP), intent(in) :: momy(ka,ia,ja)
246  real(RP), intent(in) :: rhot(ka,ia,ja)
247  real(RP), intent(in) :: dens(ka,ia,ja)
248  real(RP), intent(in) :: qtrc(ka,ia,ja,qa)
249 
250  real(RP), intent(in) :: sflx_mw(ia,ja)
251  real(RP), intent(in) :: sflx_mu(ia,ja)
252  real(RP), intent(in) :: sflx_mv(ia,ja)
253  real(RP), intent(in) :: sflx_sh(ia,ja)
254  real(RP), intent(in) :: sflx_qv(ia,ja)
255 
256  real(RP), intent(in) :: gsqrt (ka,ia,ja,7)
257  real(RP), intent(in) :: j13g (ka,ia,ja,7)
258  real(RP), intent(in) :: j23g (ka,ia,ja,7)
259  real(RP), intent(in) :: j33g
260  real(RP), intent(in) :: mapf(ia,ja,2,4)
261  real(DP), intent(in) :: dt
262 
263  ! diagnostic variables
264  real(RP) :: pott (ka,ia,ja)
265 
266  ! deformation rate tensor
267  real(RP) :: s33_c (ka,ia,ja) ! (cell center)
268  real(RP) :: s11_c (ka,ia,ja)
269  real(RP) :: s22_c (ka,ia,ja)
270  real(RP) :: s31_c (ka,ia,ja)
271  real(RP) :: s12_c (ka,ia,ja)
272  real(RP) :: s23_c (ka,ia,ja)
273  real(RP) :: s12_z (ka,ia,ja) ! (z edge or x-y plane)
274  real(RP) :: s23_x (ka,ia,ja) ! (x edge or y-z plane)
275  real(RP) :: s31_y (ka,ia,ja) ! (y edge or z-x plane)
276  real(RP) :: s2 (ka,ia,ja) ! |S|^2
277 
278  real(RP) :: kh (ka,ia,ja) ! eddy diffusion
279 
280  real(RP) :: tend(ka,ia,ja)
281  real(RP) :: a(ka,ia,ja)
282  real(RP) :: b(ka,ia,ja)
283  real(RP) :: c(ka,ia,ja)
284  real(RP) :: d(ka)
285  real(RP) :: ap
286 
287  integer :: iis, iie
288  integer :: jjs, jje
289 
290  integer :: k, i, j, iq
291  !---------------------------------------------------------------------------
292 
293  if( io_l ) write(io_fid_log,*) '*** Physics step: Turbulence(smagorinsky)'
294 
295 #ifdef DEBUG
296  qflx_sgs_momz(:,:,:,:) = undef
297  qflx_sgs_momx(:,:,:,:) = undef
298  qflx_sgs_momy(:,:,:,:) = undef
299  qflx_sgs_rhot(:,:,:,:) = undef
300  qflx_sgs_rhoq(:,:,:,:,:) = undef
301 
302  tke_t(:,:,:) = undef
303 
304  nu(:,:,:) = undef
305  tke(:,:,:) = undef
306  pr(:,:,:) = undef
307  ri(:,:,:) = undef
308  kh(:,:,:) = undef
309  n2(:,:,:) = undef
310 
311  pott(:,:,:) = undef
312 #endif
313 
314  ! potential temperature
315  do j = js-1, je+1
316  do i = is-1, ie+1
317  do k = ks, ke
318 #ifdef DEBUG
319  call check( __line__, rhot(k,i,j) )
320  call check( __line__, dens(k,i,j) )
321 #endif
322  pott(k,i,j) = rhot(k,i,j) / dens(k,i,j)
323  enddo
324  enddo
325  enddo
326 #ifdef DEBUG
327  i = iundef; j = iundef; k = iundef
328 #endif
329 
330  !##### Start Upadate #####
331 
332  call calc_strain_tensor( &
333  s33_c, s11_c, s22_c, & ! (out)
334  s31_c, s12_c, s23_c, & ! (out)
335  s12_z, s23_x, s31_y, & ! (out)
336  s2 , & ! (out)
337  dens, momz, momx, momy, & ! (in)
338  gsqrt, j13g, j23g, j33g, mapf ) ! (in)
339 
340  do jjs = js, je, jblock
341  jje = jjs+jblock-1
342  do iis = is, ie, iblock
343  iie = iis+iblock-1
344 
345  ! Ri = N^2 / |S|^2, N^2 = g / theta * dtheta/dz
346  do j = jjs-1, jje+1
347  do i = iis-1, iie+1
348  do k = ks+1, ke-1
349 #ifdef DEBUG
350  call check( __line__, pott(k+1,i,j) )
351  call check( __line__, pott(k,i,j) )
352  call check( __line__, pott(k-1,i,j) )
353  call check( __line__, fdz(k) )
354  call check( __line__, fdz(k-1) )
355  call check( __line__, s2(k,i,j) )
356 #endif
357  n2(k,i,j) = grav * ( pott(k+1,i,j) - pott(k-1,i,j) ) * j33g &
358  / ( ( fdz(k) + fdz(k-1) ) * gsqrt(k,i,j,i_xyz) * pott(k,i,j) )
359  ri(k,i,j) = n2(k,i,j) / s2(k,i,j)
360  enddo
361  enddo
362  enddo
363 #ifdef DEBUG
364  i = iundef; j = iundef; k = iundef
365 #endif
366  do j = jjs-1, jje+1
367  do i = iis-1, iie+1
368 #ifdef DEBUG
369  call check( __line__, pott(ks+1,i,j) )
370  call check( __line__, pott(ks,i,j) )
371  call check( __line__, rfdz(ks) )
372  call check( __line__, s2(ks,i,j) )
373 #endif
374  n2(ks,i,j) = grav * ( pott(ks+1,i,j) - pott(ks,i,j) ) * j33g &
375  / ( fdz(ks) * gsqrt(ks,i,j,i_xyz) * pott(ks,i,j) )
376  ri(ks,i,j) = grav * ( pott(ks+1,i,j) - pott(ks,i,j) ) * j33g * rfdz(ks) &
377  / ( gsqrt(ks,i,j,i_xyz) * pott(ks,i,j) * s2(ks,i,j) )
378  enddo
379  enddo
380 #ifdef DEBUG
381  i = iundef; j = iundef; k = iundef
382 #endif
383  do j = jjs-1, jje+1
384  do i = iis-1, iie+1
385 #ifdef DEBUG
386  call check( __line__, pott(ke,i,j) )
387  call check( __line__, pott(ke-1,i,j) )
388  call check( __line__, rfdz(ke-1) )
389  call check( __line__, s2(ke,i,j) )
390 #endif
391  n2(ke,i,j) = grav * ( pott(ke,i,j) - pott(ke-1,i,j) ) * j33g &
392  / ( fdz(ke-1) * gsqrt(ke,i,j,i_xyz) * pott(ke,i,j) )
393  ri(ke,i,j) = grav * ( pott(ke,i,j) - pott(ke-1,i,j) ) * j33g * rfdz(ke-1) &
394  / ( gsqrt(ke,i,j,i_xyz) * pott(ke,i,j) * s2(ke,i,j) )
395  enddo
396  enddo
397 #ifdef DEBUG
398  i = iundef; j = iundef; k = iundef
399 #endif
400  do j = jjs-1, jje+1
401  do i = iis-1, iie+1
402  do k = ks, ke
403 #ifdef DEBUG
404  call check( __line__, ri(k,i,j) )
405  call check( __line__, nu_fact(k,i,j) )
406  call check( __line__, s2(k,i,j) )
407 #endif
408  if ( ri(k,i,j) < 0.0_rp ) then ! stable
409  nu(k,i,j) = nu_fact(k,i,j) &
410  * sqrt( s2(k,i,j) * (1.0_rp - fmc*ri(k,i,j)) )
411  else if ( ri(k,i,j) < ric ) then ! weakly stable
412  nu(k,i,j) = nu_fact(k,i,j) &
413  * sqrt( s2(k,i,j) ) * ( 1.0_rp - ri(k,i,j)*rric )**4
414  else ! strongly stable
415  nu(k,i,j) = 0.0_rp
416  endif
417  enddo
418  enddo
419  enddo
420 #ifdef DEBUG
421  i = iundef; j = iundef; k = iundef
422 #endif
423  ! Pr = nu_m / nu_h = fm / fh
424  do j = jjs-1, jje+1
425  do i = iis-1, iie+1
426  do k = ks, ke
427 #ifdef DEBUG
428  call check( __line__, ri(k,i,j) )
429 #endif
430  if ( ri(k,i,j) < 0.0_rp ) then ! stable
431  pr(k,i,j) = sqrt( ( 1.0_rp - fmc*ri(k,i,j) ) &
432  / ( 1.0_rp - fhb*ri(k,i,j) ) ) * prn
433  else if ( ri(k,i,j) < ric ) then ! weakly stable
434  pr(k,i,j) = prn / ( 1.0_rp - prnovric * ri(k,i,j) )
435  else ! strongly stable
436  pr(k,i,j) = 1.0_rp
437  endif
438  kh(k,i,j) = min( nu(k,i,j)/pr(k,i,j), atmos_phy_tb_smg_nu_max )
439  nu(k,i,j) = min( nu(k,i,j), atmos_phy_tb_smg_nu_max )
440  pr(k,i,j) = nu(k,i,j) / ( kh(k,i,j) + ( 0.5_rp - sign(0.5_rp, kh(k,i,j)-eps) ) )
441  enddo
442  enddo
443  enddo
444 #ifdef DEBUG
445  i = iundef; j = iundef; k = iundef
446 #endif
447 
448  ! tke = (nu/(Ck * Delta))^2 = ( nu * Cs / Ck )^2 / ( Cs * Delta )^2
449  ! Sullivan et al. (1994)
450  do j = jjs, jje+1
451  do i = iis, iie+1
452  do k = ks, ke
453 #ifdef DEBUG
454  call check( __line__, nu(k,i,j) )
455  call check( __line__, nu_fact(k,i,j) )
456 #endif
457  tke(k,i,j) = ( nu(k,i,j) * cs / ck )**2 / nu_fact(k,i,j)
458  enddo
459  enddo
460  enddo
461 #ifdef DEBUG
462  i = iundef; j = iundef; k = iundef
463 #endif
464 
465  !##### momentum equation (z) #####
466  ! (cell center)
467  do j = jjs, jje
468  do i = iis, iie
469  do k = ks+1, ke-1
470 #ifdef DEBUG
471  call check( __line__, dens(k,i,j) )
472  call check( __line__, nu(k,i,j) )
473  call check( __line__, s33_c(k,i,j) )
474  call check( __line__, s11_c(k,i,j) )
475  call check( __line__, s22_c(k,i,j) )
476  call check( __line__, tke(k,i,j) )
477 #endif
478  qflx_sgs_momz(k,i,j,zdir) = dens(k,i,j) * ( &
479  - 2.0_rp * nu(k,i,j) &
480  * ( s33_c(k,i,j) - ( s11_c(k,i,j) + s22_c(k,i,j) + s33_c(k,i,j) ) * oneoverthree * tke_fact ) &
481  + twooverthree * tke(k,i,j) * tke_fact )
482  enddo
483  enddo
484  enddo
485 #ifdef DEBUG
486  i = iundef; j = iundef; k = iundef
487 #endif
488  do j = jjs, jje
489  do i = iis, iie
490  qflx_sgs_momz(ks,i,j,zdir) = 0.0_rp ! just above bottom boundary
491  qflx_sgs_momz(ke,i,j,zdir) = 0.0_rp ! just below top boundary
492  enddo
493  enddo
494 #ifdef DEBUG
495  i = iundef; j = iundef; k = iundef
496 #endif
497  ! (y edge)
498  do j = jjs, jje
499  do i = iis-1, iie
500  do k = ks, ke-1
501 #ifdef DEBUG
502  call check( __line__, dens(k,i,j) )
503  call check( __line__, dens(k,i+1,j) )
504  call check( __line__, dens(k+1,i,j) )
505  call check( __line__, dens(k+1,i+1,j) )
506  call check( __line__, nu(k,i,j) )
507  call check( __line__, nu(k,i+1,j) )
508  call check( __line__, nu(k+1,i,j) )
509  call check( __line__, nu(k+1,i+1,j) )
510  call check( __line__, s31_y(k,i,j) )
511 #endif
512  qflx_sgs_momz(k,i,j,xdir) = - 0.125_rp & ! 2.0 / 4 / 4
513  * ( dens(k,i,j)+dens(k+1,i,j)+dens(k,i+1,j)+dens(k+1,i+1,j) ) &
514  * ( nu(k,i,j)+nu(k+1,i,j)+nu(k,i+1,j)+nu(k+1,i+1,j)) &
515  * s31_y(k,i,j)
516  enddo
517  enddo
518  enddo
519 #ifdef DEBUG
520  i = iundef; j = iundef; k = iundef
521 #endif
522  ! (x edge)
523  do j = jjs-1, jje
524  do i = iis, iie
525  do k = ks, ke-1
526 #ifdef DEBUG
527  call check( __line__, dens(k,i,j) )
528  call check( __line__, dens(k,i,j+1) )
529  call check( __line__, dens(k+1,i,j) )
530  call check( __line__, dens(k+1,i,j+1) )
531  call check( __line__, nu(k,i,j) )
532  call check( __line__, nu(k,i,j+1) )
533  call check( __line__, nu(k+1,i,j) )
534  call check( __line__, nu(k+1,i,j+1) )
535  call check( __line__, s23_x(k,i,j) )
536 #endif
537  qflx_sgs_momz(k,i,j,ydir) = - 0.125_rp & ! 2/4/4
538  * ( dens(k,i,j)+dens(k+1,i,j)+dens(k,i,j+1)+dens(k+1,i,j+1) ) &
539  * ( nu(k,i,j)+nu(k+1,i,j)+nu(k,i,j+1)+nu(k+1,i,j+1) ) &
540  * s23_x(k,i,j)
541  enddo
542  enddo
543  enddo
544 #ifdef DEBUG
545  i = iundef; j = iundef; k = iundef
546 #endif
547 
548  if ( atmos_phy_tb_smg_implicit ) then
549  call calc_tend_momz( tend, & ! (out)
550  qflx_sgs_momz, & ! (in)
551  gsqrt, j13g, j23g, j33g, mapf, & ! (in)
552  iis, iie, jjs, jje ) ! (in)
553 
554  do j = jjs, jje
555  do i = iis, iie
556 
557  ap = - fouroverthree * dt &
558  * dens(ks+1,i,j)*nu(ks+1,i,j) &
559  * rcdz(ks+1) / gsqrt(ks+1,i,j,i_xyz)
560  a(ks,i,j) = ap * rfdz(ks) / gsqrt(ks,i,j,i_xyw)
561  c(ks,i,j) = 0.0_rp
562  b(ks,i,j) = - a(ks,i,j) + 0.5_rp * ( dens(ks,i,j)+dens(ks+1,i,j) )
563  do k = ks+1, ke-2
564  c(k,i,j) = ap * rfdz(k+1) / gsqrt(k+1,i,j,i_xyw)
565  ap = - fouroverthree * dt &
566  * dens(k+1,i,j)*nu(k+1,i,j) &
567  * rcdz(k+1) / gsqrt(k+1,i,j,i_xyz)
568  a(k,i,j) = ap * rfdz(k) / gsqrt(k,i,j,i_xyw)
569  b(k,i,j) = - a(k,i,j) - c(k,i,j) + 0.5_rp * ( dens(k,i,j)+dens(k+1,i,j) )
570  end do
571  a(ke-1,i,j) = 0.0_rp
572  c(ke-1,i,j) = ap * rfdz(ke) / gsqrt(ke,i,j,i_xyw)
573  b(ke-1,i,j) = - c(ke-1,i,j) + 0.5_rp * ( dens(ke-1,i,j)+dens(ke,i,j) )
574 
575  do k = ks, ke-1
576  d(k) = tend(k,i,j)
577  end do
578 
579  call diffusion_solver( &
580  tend(:,i,j), & ! (out)
581  a(:,i,j), b(:,i,j), c(:,i,j), d, & ! (in)
582  ke-1 ) ! (in)
583 
584  do k = ks+1, ke-1
585  qflx_sgs_momz(k,i,j,zdir) = qflx_sgs_momz(k,i,j,zdir) &
586  - fouroverthree * dens(k,i,j) * nu(k,i,j) * dt &
587  * ( tend(k,i,j) - tend(k-1,i,j) ) * rcdz(k) / gsqrt(k,i,j,i_xyz)
588  end do
589 
590  end do
591  end do
592 
593  end if
594 
595  !##### momentum equation (x) #####
596  ! (y edge)
597  do j = jjs, jje
598  do i = iis, iie
599  do k = ks, ke-1
600 #ifdef DEBUG
601  call check( __line__, dens(k,i,j) )
602  call check( __line__, dens(k,i+1,j) )
603  call check( __line__, dens(k+1,i,j) )
604  call check( __line__, dens(k+1,i+1,j) )
605  call check( __line__, nu(k,i,j) )
606  call check( __line__, nu(k,i+1,j) )
607  call check( __line__, nu(k+1,i,j) )
608  call check( __line__, nu(k+1,i+1,j) )
609  call check( __line__, s31_y(k,i,j) )
610 #endif
611  qflx_sgs_momx(k,i,j,zdir) = - 0.125_rp & ! 2/4/4
612  * ( dens(k,i,j)+dens(k+1,i,j)+dens(k,i+1,j)+dens(k+1,i+1,j) ) &
613  * ( nu(k,i,j)+nu(k+1,i,j)+nu(k,i+1,j)+nu(k+1,i+1,j) ) &
614  * s31_y(k,i,j)
615  enddo
616  enddo
617  enddo
618 #ifdef DEBUG
619  i = iundef; j = iundef; k = iundef
620 #endif
621  do j = jjs, jje
622  do i = iis, iie
623  qflx_sgs_momx(ks-1,i,j,zdir) = 0.0_rp ! bottom boundary
624  qflx_sgs_momx(ke ,i,j,zdir) = 0.0_rp ! top boundary
625  enddo
626  enddo
627 #ifdef DEBUG
628  i = iundef; j = iundef; k = iundef
629 #endif
630  ! (cell center)
631  do j = jjs, jje
632  do i = iis, iie+1
633  do k = ks, ke
634 #ifdef DEBUG
635  call check( __line__, dens(k,i,j) )
636  call check( __line__, nu(k,i,j) )
637  call check( __line__, s11_c(k,i,j) )
638  call check( __line__, s22_c(k,i,j) )
639  call check( __line__, s33_c(k,i,j) )
640  call check( __line__, tke(k,i,j) )
641 #endif
642  qflx_sgs_momx(k,i,j,xdir) = dens(k,i,j) * ( &
643  - 2.0_rp * nu(k,i,j) &
644  * ( s11_c(k,i,j) - ( s11_c(k,i,j) + s22_c(k,i,j) + s33_c(k,i,j) ) * oneoverthree * tke_fact ) &
645  + twooverthree * tke(k,i,j) * tke_fact )
646  enddo
647  enddo
648  enddo
649 #ifdef DEBUG
650  i = iundef; j = iundef; k = iundef
651 #endif
652  ! (z edge)
653  do j = jjs-1, jje
654  do i = iis, iie
655  do k = ks, ke
656 #ifdef DEBUG
657  call check( __line__, dens(k,i,j) )
658  call check( __line__, dens(k,i+1,j) )
659  call check( __line__, dens(k,i,j+1) )
660  call check( __line__, dens(k,i+1,j+1) )
661  call check( __line__, nu(k,i,j) )
662  call check( __line__, nu(k,i+1,j) )
663  call check( __line__, nu(k,i,j+1) )
664  call check( __line__, nu(k,i+1,j+1) )
665  call check( __line__, s12_z(k,i,j) )
666 #endif
667  qflx_sgs_momx(k,i,j,ydir) = - 0.125_rp & ! 2/4/4
668  * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) &
669  * ( nu(k,i,j)+nu(k,i+1,j)+nu(k,i,j+1)+nu(k,i+1,j+1) ) &
670  * s12_z(k,i,j)
671  enddo
672  enddo
673  enddo
674 #ifdef DEBUG
675  i = iundef; j = iundef; k = iundef
676 #endif
677 
678  if ( atmos_phy_tb_smg_implicit ) then
679  call calc_tend_momx( tend, & ! (out)
680  qflx_sgs_momx, & ! (in)
681  gsqrt, j13g, j23g, j33g, mapf, & ! (in)
682  iis, iie, jjs, jje ) ! (in)
683 
684  do j = jjs, jje
685  do i = iis, iie
686 
687  ap = - dt * 0.25_rp * ( dens(ks ,i ,j)*nu(ks ,i ,j) &
688  + dens(ks+1,i ,j)*nu(ks+1,i ,j) &
689  + dens(ks ,i+1,j)*nu(ks ,i+1,j) &
690  + dens(ks+1,i+1,j)*nu(ks+1,i+1,j) ) &
691  * rfdz(ks) / gsqrt(ks,i,j,i_uyw)
692  a(ks,i,j) = ap * rcdz(ks) / gsqrt(ks,i,j,i_uyz)
693  c(ks,i,j) = 0.0_rp
694  b(ks,i,j) = - a(ks,i,j) + 0.5_rp * ( dens(ks,i,j)+dens(ks,i+1,j) )
695  do k = ks+1, ke-1
696  c(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,i_uyz)
697  ap = - dt * 0.25_rp * ( dens(k ,i ,j)*nu(k ,i ,j) &
698  + dens(k+1,i ,j)*nu(k+1,i ,j) &
699  + dens(k ,i+1,j)*nu(k ,i+1,j) &
700  + dens(k+1,i+1,j)*nu(k+1,i+1,j) ) &
701  * rfdz(k) / gsqrt(k,i,j,i_uyw)
702  a(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,i_uyz)
703  b(k,i,j) = - a(k,i,j) - c(k,i,j) + 0.5_rp * ( dens(k,i,j)+dens(k,i+1,j) )
704  end do
705  a(ke,i,j) = 0.0_rp
706  c(ke,i,j) = ap * rcdz(ke) / gsqrt(ke,i,j,i_uyz)
707  b(ke,i,j) = - c(ke,i,j) + 0.5_rp * ( dens(ke,i,j)+dens(ke,i+1,j) )
708 
709  do k = ks, ke
710  d(k) = tend(k,i,j)
711  end do
712 
713  call diffusion_solver( &
714  tend(:,i,j), & ! (out)
715  a(:,i,j), b(:,i,j), c(:,i,j), d, & ! (in)
716  ke ) ! (in)
717 
718  do k = ks, ke-1
719  qflx_sgs_momx(k,i,j,zdir) = qflx_sgs_momx(k,i,j,zdir) &
720  - 0.25_rp * ( dens(k ,i ,j)*nu(k ,i ,j) &
721  + dens(k+1,i ,j)*nu(k+1,i ,j) &
722  + dens(k ,i+1,j)*nu(k ,i+1,j) &
723  + dens(k+1,i+1,j)*nu(k+1,i+1,j) ) &
724  * dt * ( tend(k+1,i,j) - tend(k,i,j) ) * rfdz(k) / gsqrt(k,i,j,i_uyw)
725  end do
726 
727  end do
728  end do
729 
730  end if
731 
732  !##### momentum equation (y) #####
733  ! (x edge)
734  do j = jjs, jje
735  do i = iis, iie
736  do k = ks, ke-1
737 #ifdef DEBUG
738  call check( __line__, dens(k,i,j) )
739  call check( __line__, dens(k,i,j+1) )
740  call check( __line__, dens(k+1,i,j) )
741  call check( __line__, dens(k+1,i,j+1) )
742  call check( __line__, nu(k,i,j) )
743  call check( __line__, nu(k,i,j+1) )
744  call check( __line__, nu(k+1,i,j) )
745  call check( __line__, nu(k+1,i,j+1) )
746  call check( __line__, s23_x(k,i,j) )
747 #endif
748  qflx_sgs_momy(k,i,j,zdir) = - 0.125_rp & ! 2/4/4
749  * ( dens(k,i,j)+dens(k+1,i,j)+dens(k,i,j+1)+dens(k+1,i,j+1) ) &
750  * ( nu(k,i,j)+nu(k+1,i,j)+nu(k,i,j+1)+nu(k+1,i,j+1) ) &
751  * s23_x(k,i,j)
752  enddo
753  enddo
754  enddo
755 #ifdef DEBUG
756  i = iundef; j = iundef; k = iundef
757 #endif
758  do j = jjs, jje
759  do i = iis, iie
760  qflx_sgs_momy(ks-1,i,j,zdir) = 0.0_rp ! bottom boundary
761  qflx_sgs_momy(ke ,i,j,zdir) = 0.0_rp ! top boundary
762  enddo
763  enddo
764 #ifdef DEBUG
765  i = iundef; j = iundef; k = iundef
766 #endif
767 
768  ! (z edge)
769  do j = jjs, jje
770  do i = iis-1, iie
771  do k = ks, ke
772 #ifdef DEBUG
773  call check( __line__, dens(k,i,j) )
774  call check( __line__, dens(k,i+1,j) )
775  call check( __line__, dens(k,i,j+1) )
776  call check( __line__, dens(k,i+1,j+1) )
777  call check( __line__, nu(k,i,j) )
778  call check( __line__, nu(k,i+1,j) )
779  call check( __line__, nu(k,i,j+1) )
780  call check( __line__, nu(k,i+1,j+1) )
781  call check( __line__, s12_z(k,i,j) )
782 #endif
783  qflx_sgs_momy(k,i,j,xdir) = - 0.125_rp & !
784  * ( dens(k,i,j)+dens(k,i+1,j)+dens(k,i,j+1)+dens(k,i+1,j+1) ) &
785  * ( nu(k,i,j)+nu(k,i+1,j)+nu(k,i,j+1)+nu(k,i+1,j+1) ) &
786  * s12_z(k,i,j)
787  enddo
788  enddo
789  enddo
790 #ifdef DEBUG
791  i = iundef; j = iundef; k = iundef
792 #endif
793 
794  ! (z-x plane)
795  do j = jjs, jje+1
796  do i = iis, iie
797  do k = ks, ke
798 #ifdef DEBUG
799  call check( __line__, dens(k,i,j) )
800  call check( __line__, nu(k,i,j) )
801  call check( __line__, s11_c(k,i,j) )
802  call check( __line__, s22_c(k,i,j) )
803  call check( __line__, s33_c(k,i,j) )
804  call check( __line__, tke(k,i,j) )
805 #endif
806  qflx_sgs_momy(k,i,j,ydir) = dens(k,i,j) * ( &
807  - 2.0_rp * nu(k,i,j) &
808  * ( s22_c(k,i,j) - ( s11_c(k,i,j) + s22_c(k,i,j) + s33_c(k,i,j) ) * oneoverthree * tke_fact ) &
809  + twooverthree * tke(k,i,j) * tke_fact)
810  enddo
811  enddo
812  enddo
813 #ifdef DEBUG
814  i = iundef; j = iundef; k = iundef
815 #endif
816 
817  if ( atmos_phy_tb_smg_implicit ) then
818  call calc_tend_momy( tend, & ! (out)
819  qflx_sgs_momy, & ! (in)
820  gsqrt, j13g, j23g, j33g, mapf, & ! (in)
821  iis, iie, jjs, jje ) ! (in)
822 
823  do j = jjs, jje
824  do i = iis, iie
825 
826  ap = - dt * 0.25_rp * ( dens(ks ,i,j )*nu(ks ,i,j ) &
827  + dens(ks+1,i,j )*nu(ks+1,i,j ) &
828  + dens(ks ,i,j+1)*nu(ks ,i,j+1) &
829  + dens(ks+1,i,j+1)*nu(ks+1,i,j+1) ) &
830  * rfdz(ks) / gsqrt(ks,i,j,i_xvw)
831  a(ks,i,j) = ap * rcdz(ks) / gsqrt(ks,i,j,i_xvz)
832  c(ks,i,j) = 0.0_rp
833  b(ks,i,j) = - a(ks,i,j) + 0.5_rp * ( dens(ks,i,j)+dens(ks,i,j+1) )
834  do k = ks+1, ke-1
835  c(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,i_xvz)
836  ap = - dt * 0.25_rp * ( dens(k ,i,j )*nu(k ,i,j ) &
837  + dens(k+1,i,j )*nu(k+1,i,j ) &
838  + dens(k ,i,j+1)*nu(k ,i,j+1) &
839  + dens(k+1,i,j+1)*nu(k+1,i,j+1) ) &
840  * rfdz(k) / gsqrt(k,i,j,i_xvw)
841  a(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,i_xvz)
842  b(k,i,j) = - a(k,i,j) - c(k,i,j) + 0.5_rp * ( dens(k,i,j)+dens(k,i,j+1) )
843  end do
844  a(ke,i,j) = 0.0_rp
845  c(ke,i,j) = ap * rcdz(ke) / gsqrt(ke,i,j,i_xvz)
846  b(ke,i,j) = - c(ke,i,j) + 0.5_rp * ( dens(ke,i,j)+dens(ke,i,j+1) )
847 
848  do k = ks, ke
849  d(k) = tend(k,i,j)
850  end do
851 
852  call diffusion_solver( &
853  tend(:,i,j), & ! (out)
854  a(:,i,j), b(:,i,j), c(:,i,j), d, & ! (in)
855  ke ) ! (in)
856 
857  do k = ks, ke-1
858  qflx_sgs_momy(k,i,j,zdir) = qflx_sgs_momy(k,i,j,zdir) &
859  - 0.25_rp * ( dens(k ,i,j )*nu(k ,i,j ) &
860  + dens(k+1,i,j )*nu(k+1,i,j ) &
861  + dens(k ,i,j+1)*nu(k ,i,j+1) &
862  + dens(k+1,i,j+1)*nu(k+1,i,j+1) ) &
863  * dt * ( tend(k+1,i,j) - tend(k,i,j) ) * rfdz(k) / gsqrt(k,i,j,i_xvw)
864  end do
865 
866  end do
867  end do
868 
869  end if
870 
871  !##### Thermodynamic Equation #####
872 
873  if ( atmos_phy_tb_smg_implicit ) then
874 
875  do j = jjs, jje
876  do i = iis, iie
877 
878  ap = - dt * 0.25_rp * ( dens(ks,i,j)+dens(ks+1,i,j) ) &
879  * ( kh(ks,i,j)+kh(ks+1,i,j) ) &
880  * rfdz(ks) / gsqrt(ks,i,j,i_xyw)
881  a(ks,i,j) = ap * rcdz(ks) / gsqrt(ks,i,j,i_xyz)
882  c(ks,i,j) = 0.0_rp
883  b(ks,i,j) = - a(ks,i,j) + dens(ks,i,j)
884  do k = ks+1, ke-1
885  c(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,i_xyz)
886  ap = - dt * 0.25_rp * ( dens(k,i,j)+dens(k+1,i,j) ) &
887  * ( kh(k,i,j)+kh(k+1,i,j) ) &
888  * rfdz(k) / gsqrt(k,i,j,i_xyw)
889  a(k,i,j) = ap * rcdz(k) / gsqrt(k,i,j,i_xyz)
890  b(k,i,j) = - a(k,i,j) - c(k,i,j) + dens(k,i,j)
891  end do
892  a(ke,i,j) = 0.0_rp
893  c(ke,i,j) = ap * rcdz(ke) / gsqrt(ke,i,j,i_xyz)
894  b(ke,i,j) = - c(ke,i,j) + dens(ke,i,j)
895 
896  end do
897  end do
898 
899  end if
900 
901  call calc_flux_phi( &
902  qflx_sgs_rhot, &
903  dens, pott, kh, 1.0_rp, &
904  gsqrt, j13g, j23g, j33g, mapf, &
905  a, b, c, dt, &
906  atmos_phy_tb_smg_implicit, &
907  iis, iie, jjs, jje )
908 
909  enddo
910  enddo
911 
912 
913  !##### Tracers #####
914  do iq = 1, qa
915 
916  do jjs = js, je, jblock
917  jje = jjs+jblock-1
918  do iis = is, ie, iblock
919  iie = iis+iblock-1
920 
921  call calc_flux_phi( &
922  qflx_sgs_rhoq(:,:,:,:,iq), &
923  dens, qtrc(:,:,:,iq), kh, 1.0_rp, &
924  gsqrt, j13g, j23g, j33g, mapf, &
925  a, b, c, dt, &
926  atmos_phy_tb_smg_implicit, &
927  iis, iie, jjs, jje )
928 
929 
930  enddo
931  enddo
932 #ifdef DEBUG
933  iis = iundef; iie = iundef; jjs = iundef; jje = iundef
934 #endif
935 
936  enddo ! scalar quantities loop
937 #ifdef DEBUG
938  iq = iundef
939 #endif
940 
941  tke_t = 0.0_rp
942 
943  return
integer, public is
start point of inner domain: x, local
integer, public i_xvz
subroutine, public atmos_phy_tb_calc_tend_momz(MOMZ_t_TB, QFLX_MOMZ, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
integer, public je
end point of inner domain: y, local
real(rp), dimension(:), allocatable, public grid_rcdy
reciprocal of center-dy
subroutine, public atmos_phy_tb_diffusion_solver(phi, a, b, c, d, KE_TB)
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)
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
subroutine, public atmos_phy_tb_calc_tend_momy(MOMY_t_TB, QFLX_MOMY, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
logical, public io_l
output log or not? (this process)
Definition: scale_stdio.F90:59
real(rp), dimension(:), allocatable, public grid_rcdx
reciprocal of center-dx
integer, parameter, public zdir
integer, parameter, public ydir
integer, public ke
end point of inner domain: z, local
integer, parameter, public xdir
integer, public qa
integer, public i_xy
real(rp), dimension(:), allocatable, public grid_rfdy
reciprocal of face-dy
subroutine, public check(current_line, v)
Undefined value checker.
Definition: scale_debug.F90:58
real(rp), dimension(:), allocatable, public grid_rcdz
reciprocal of center-dz
integer, public i_xvw
module grid index
subroutine, public atmos_phy_tb_calc_tend_phi(phi_t_TB, QFLX_phi, GSQRT, J13G, J23G, J33G, MAPF, IIS, IIE, JJS, JJE)
module TRACER
integer, public ia
of x whole cells (local, with HALO)
subroutine, public atmos_phy_tb_calc_tend_momx(MOMX_t_TB, QFLX_MOMX, 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 ka
of z whole cells (local, with HALO)
integer, public i_uy
integer, public jblock
block size for cache blocking: y
integer, public i_xyw
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:48
integer, public js
start point of inner domain: y, local
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)
integer, public i_uyw
module CONSTANT
Definition: scale_const.F90:14
integer, public ks
start point of inner domain: z, local
module ATMOSPHERE / Physics Turbulence
module GRID (cartesian)
integer, public i_xv
integer, public ie
end point of inner domain: x, local
real(rp), public const_eps
small number
Definition: scale_const.F90:36
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 io_fid_log
Log file ID.
Definition: scale_stdio.F90:56
integer, public i_uv
integer, public ja
of y whole cells (local, with HALO)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ mixlen()

real(rp) function scale_atmos_phy_tb_smg::mixlen ( real(rp), intent(in)  dz,
real(rp), intent(in)  dx,
real(rp), intent(in)  dy,
real(rp), intent(in)  z,
real(rp), intent(in)  filter_fact 
)

Definition at line 948 of file scale_atmos_phy_tb_smg.F90.

References scale_const::const_karman, and fact().

Referenced by atmos_phy_tb_smg_setup().

948  use scale_const, only: &
949  karman => const_karman
950  implicit none
951  real(RP), intent(in) :: dz
952  real(RP), intent(in) :: dx
953  real(RP), intent(in) :: dy
954  real(RP), intent(in) :: z
955  real(RP), intent(in) :: filter_fact
956  real(RP) :: mixlen ! (out)
957 
958  real(RP) :: d0
959 
960  d0 = fact(dz, dx, dy) * filter_fact * ( dz * dx * dy )**oneoverthree ! Scotti et al. (1993)
961  if ( atmos_phy_tb_smg_bottom ) then
962  mixlen = sqrt( 1.0_rp / ( 1.0_rp/d0**2 + 1.0_rp/(karman*z)**2 ) ) ! Brown et al. (1994)
963  else
964  mixlen = d0
965  end if
966 
967  return
real(rp), parameter, public const_karman
von Karman constant
Definition: scale_const.F90:52
module CONSTANT
Definition: scale_const.F90:14
Here is the call graph for this function:
Here is the caller graph for this function:

◆ fact()

real(rp) function scale_atmos_phy_tb_smg::fact ( real(rp), intent(in)  dz,
real(rp), intent(in)  dx,
real(rp), intent(in)  dy 
)

Definition at line 971 of file scale_atmos_phy_tb_smg.F90.

Referenced by mixlen().

971  real(RP), intent(in) :: dz
972  real(RP), intent(in) :: dx
973  real(RP), intent(in) :: dy
974  real(RP) :: fact ! (out)
975 
976  real(RP), parameter :: oot = -1.0_rp/3.0_rp
977  real(RP), parameter :: fot = 5.0_rp/3.0_rp
978  real(RP), parameter :: eot = 11.0_rp/3.0_rp
979  real(RP), parameter :: tof = -3.0_rp/4.0_rp
980  real(RP) :: a1, a2, b1, b2, dmax
981 
982 
983  dmax = max(dz, dx, dy)
984  if ( dz == dmax ) then
985  a1 = dx / dmax
986  a2 = dy / dmax
987  else if ( dx == dmax ) then
988  a1 = dz / dmax
989  a2 = dy / dmax
990  else ! dy == dmax
991  a1 = dz / dmax
992  a2 = dx / dmax
993  end if
994  b1 = atan( a1/a2 )
995  b2 = atan( a2/a1 )
996 
997  fact = 1.736_rp * (a1*a2)**oot &
998  * ( 4.0_rp*p1(b1)*a1**oot + 0.222_rp*p2(b1)*a1**fot + 0.077*p3(b1)*a1**eot - 3.0_rp*b1 &
999  + 4.0_rp*p1(b2)*a2**oot + 0.222_rp*p2(b2)*a2**fot + 0.077*p3(b2)*a2**eot - 3.0_rp*b2 &
1000  )**tof
1001  return
Here is the caller graph for this function: