SCALE-RM
scale_atmos_phy_tb_dns.F90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
10 !-------------------------------------------------------------------------------
11 #include "scalelib.h"
13  !-----------------------------------------------------------------------------
14  !
15  !++ used modules
16  !
17  use scale_precision
18  use scale_io
19  use scale_prof
21  use scale_tracer
22 
23 #ifdef DEBUG
24  use scale_debug, only: &
25  check
26  use scale_const, only: &
27  undef => const_undef, &
28  iundef => const_undef2
29 #endif
30  !-----------------------------------------------------------------------------
31  implicit none
32  private
33  !-----------------------------------------------------------------------------
34  !
35  !++ Public procedure
36  !
37  public :: atmos_phy_tb_dns_config
38  public :: atmos_phy_tb_dns_setup
39  public :: atmos_phy_tb_dns
40 
41  !-----------------------------------------------------------------------------
42  !
43  !++ Public parameters & variables
44  !
45  !-----------------------------------------------------------------------------
46  !
47  !++ Private procedure
48  !
49  !-----------------------------------------------------------------------------
50  !
51  !++ Private parameters & variables
52  !
53  real(RP), private :: ATMOS_PHY_TB_DNS_NU = 1.512e-5_rp ! [m2/s] kinematic viscosity coefficient for air at 20degC
54 ! real(RP), private :: mu = 1.8E-5_RP ! [m2/s] molecular diffusive coefficient for air at 20degC
55  real(RP), private :: ATMOS_PHY_TB_DNS_MU = 1.512e-5_rp ! same as NU (needed based on hyposes. see Mellado 2010)
56 
57  !-----------------------------------------------------------------------------
58 contains
59  !-----------------------------------------------------------------------------
61  subroutine atmos_phy_tb_dns_config( &
62  TYPE_TB, &
63  I_TKE_out )
64  use scale_prc, only: &
65  prc_abort
66  implicit none
67 
68  character(len=*), intent(in) :: type_tb
69  integer, intent(out) :: i_tke_out
70  !---------------------------------------------------------------------------
71 
72  log_newline
73  log_info("ATMOS_PHY_TB_dns_config",*) 'Setup'
74  log_info("ATMOS_PHY_TB_dns_config",*) 'Tracers for Deardorff (1980) 1.5th TKE Model'
75 
76  if ( type_tb /= 'DNS' ) then
77  log_error("ATMOS_PHY_TB_dns_config",*) 'ATMOS_PHY_TB_TYPE is not DNS. Check!'
78  call prc_abort
79  endif
80 
81  i_tke_out = -1
82 
83  return
84  end subroutine atmos_phy_tb_dns_config
85  !-----------------------------------------------------------------------------
86  subroutine atmos_phy_tb_dns_setup( &
87  CDZ, CDX, CDY, CZ )
88  use scale_prc, only: &
89  prc_abort
90  implicit none
91 
92  real(rp), intent(in) :: cdz(ka)
93  real(rp), intent(in) :: cdx(ia)
94  real(rp), intent(in) :: cdy(ja)
95  real(rp), intent(in) :: cz (ka,ia,ja)
96 
97  namelist / param_atmos_phy_tb_dns / &
98  atmos_phy_tb_dns_nu, &
99  atmos_phy_tb_dns_mu
100 
101  integer :: ierr
102  !---------------------------------------------------------------------------
103 
104  log_newline
105  log_info("ATMOS_PHY_TB_dns_setup",*) 'Setup'
106  log_info("ATMOS_PHY_TB_dns_setup",*) 'Eddy Viscocity Model for DNS'
107 
108  !--- read namelist
109  rewind(io_fid_conf)
110  read(io_fid_conf,nml=param_atmos_phy_tb_dns,iostat=ierr)
111  if( ierr < 0 ) then !--- missing
112  log_info("ATMOS_PHY_TB_dns_setup",*) 'Not found namelist. Default used.'
113  elseif( ierr > 0 ) then !--- fatal error
114  log_error("ATMOS_PHY_TB_dns_setup",*) 'Not appropriate names in namelist PARAM_ATMOS_PHY_TB_DNS. Check!'
115  call prc_abort
116  endif
117  log_nml(param_atmos_phy_tb_dns)
118 
119  return
120  end subroutine atmos_phy_tb_dns_setup
121 
122  !-----------------------------------------------------------------------------
123  subroutine atmos_phy_tb_dns( &
124  qflx_sgs_MOMZ, qflx_sgs_MOMX, qflx_sgs_MOMY, &
125  qflx_sgs_rhot, qflx_sgs_rhoq, &
126  RHOQ_t, nu, Ri, Pr, &
127  MOMZ, MOMX, MOMY, RHOT, DENS, QTRC, N2, &
128  SFLX_MW, SFLX_MU, SFLX_MV, SFLX_SH, SFLX_Q, &
129  GSQRT, J13G, J23G, J33G, MAPF, dt )
131  use scale_tracer
132  use scale_const, only: &
133  grav => const_grav
134  use scale_atmos_grid_cartesc, only: &
135  rcdz => atmos_grid_cartesc_rcdz, &
136  rcdx => atmos_grid_cartesc_rcdx, &
137  rcdy => atmos_grid_cartesc_rcdy, &
138  rfdz => atmos_grid_cartesc_rfdz, &
139  rfdx => atmos_grid_cartesc_rfdx, &
141  implicit none
142 
143  ! SGS flux
144  real(rp), intent(out) :: qflx_sgs_momz(ka,ia,ja,3)
145  real(rp), intent(out) :: qflx_sgs_momx(ka,ia,ja,3)
146  real(rp), intent(out) :: qflx_sgs_momy(ka,ia,ja,3)
147  real(rp), intent(out) :: qflx_sgs_rhot(ka,ia,ja,3)
148  real(rp), intent(out) :: qflx_sgs_rhoq(ka,ia,ja,3,qa)
149 
150  real(rp), intent(inout) :: rhoq_t(ka,ia,ja,qa) ! tendency of rho * QTRC
151 
152  real(rp), intent(out) :: nu(ka,ia,ja) ! eddy viscosity (center)
153  real(rp), intent(out) :: ri(ka,ia,ja) ! Richardson number
154  real(rp), intent(out) :: pr(ka,ia,ja) ! Prantle number
155 
156  real(rp), intent(in) :: momz(ka,ia,ja)
157  real(rp), intent(in) :: momx(ka,ia,ja)
158  real(rp), intent(in) :: momy(ka,ia,ja)
159  real(rp), intent(in) :: rhot(ka,ia,ja)
160  real(rp), intent(in) :: dens(ka,ia,ja)
161  real(rp), intent(in) :: qtrc(ka,ia,ja,qa)
162  real(rp), intent(in) :: n2(ka,ia,ja)
163 
164  real(rp), intent(in) :: sflx_mw(ia,ja)
165  real(rp), intent(in) :: sflx_mu(ia,ja)
166  real(rp), intent(in) :: sflx_mv(ia,ja)
167  real(rp), intent(in) :: sflx_sh(ia,ja)
168  real(rp), intent(in) :: sflx_q (ia,ja,qa)
169 
170  real(rp), intent(in) :: gsqrt (ka,ia,ja,7)
171  real(rp), intent(in) :: j13g (ka,ia,ja,7)
172  real(rp), intent(in) :: j23g (ka,ia,ja,7)
173  real(rp), intent(in) :: j33g
174  real(rp), intent(in) :: mapf (ia,ja,2,4)
175  real(dp), intent(in) :: dt
176 
177  real(rp) :: pott(ka,ia,ja)
178 
179  integer :: iis, iie
180  integer :: jjs, jje
181 
182  integer :: k, i, j, iq
183  !---------------------------------------------------------------------------
184 
185  log_progress(*) 'atmosphere / physics / turbulence / DNS'
186 
187 #ifdef DEBUG
188  qflx_sgs_momz(:,:,:,:) = undef
189  qflx_sgs_momx(:,:,:,:) = undef
190  qflx_sgs_momy(:,:,:,:) = undef
191  qflx_sgs_rhot(:,:,:,:) = undef
192  qflx_sgs_rhoq(:,:,:,:,:) = undef
193 
194  pott(:,:,:) = undef
195 #endif
196 
197  nu(:,:,:) = 0.0_rp
198  ri(:,:,:) = 0.0_rp
199  pr(:,:,:) = 1.0_rp
200 
201  ! potential temperature
202  do j = js-1, je+1
203  do i = is-1, ie+1
204  do k = ks, ke
205  pott(k,i,j) = rhot(k,i,j) / dens(k,i,j)
206  enddo
207  enddo
208  enddo
209 
210  !##### Start Upadate #####
211  do jjs = js, je, jblock
212  jje = jjs+jblock-1
213  do iis = is, ie, iblock
214  iie = iis+iblock-1
215 
216  !##### momentum equation (z) #####
217  ! (cell center)
218  do j = jjs, jje
219  do i = iis, iie
220  do k = ks+1, ke-1
221  qflx_sgs_momz(k,i,j,zdir) = -atmos_phy_tb_dns_nu * ( momz(k,i,j)-momz(k-1,i,j) ) * rcdz(k)
222  enddo
223  enddo
224  enddo
225 
226  do j = jjs, jje
227  do i = iis, iie
228  qflx_sgs_momz(ks,i,j,zdir) = 0.0_rp ! bottom boundary
229  qflx_sgs_momz(ke,i,j,zdir) = 0.0_rp ! top boundary
230  enddo
231  enddo
232 
233  ! (y edge)
234  do j = jjs, jje
235  do i = iis-1, iie
236  do k = ks, ke-1
237  qflx_sgs_momz(k,i,j,xdir) = -atmos_phy_tb_dns_nu * ( momz(k,i+1,j)-momz(k,i,j) ) * rfdx(i) * mapf(i,j,1,i_xy)
238  enddo
239  enddo
240  enddo
241 
242  ! (x edge)
243  do j = jjs-1, jje
244  do i = iis, iie
245  do k = ks, ke-1
246  qflx_sgs_momz(k,i,j,ydir) = -atmos_phy_tb_dns_nu * ( momz(k,i,j+1)-momz(k,i,j) ) * rfdy(j) * mapf(i,j,2,i_xy)
247  enddo
248  enddo
249  enddo
250 
251  !##### momentum equation (x) #####
252  ! (y edge)
253  do j = jjs, jje
254  do i = iis, iie
255  do k = ks, ke-1
256  qflx_sgs_momx(k,i,j,zdir) = -atmos_phy_tb_dns_nu * ( momx(k+1,i,j)-momx(k,i,j) ) * rfdz(k)
257  enddo
258  enddo
259  enddo
260 
261  do j = jjs, jje
262  do i = iis, iie
263  qflx_sgs_momx(ks-1,i,j,zdir) = 0.0_rp ! bottom boundary
264  qflx_sgs_momx(ke ,i,j,zdir) = 0.0_rp ! top boundary
265  enddo
266  enddo
267 
268  ! (cell center)
269  do j = jjs, jje
270  do i = iis, iie+1
271  do k = ks, ke
272  qflx_sgs_momx(k,i,j,xdir) = -atmos_phy_tb_dns_nu * ( momx(k,i,j)-momx(k,i-1,j) ) * rcdx(i) * mapf(i,j,1,i_uy)
273  enddo
274  enddo
275  enddo
276 
277  ! (z edge)
278  do j = jjs-1, jje
279  do i = iis, iie
280  do k = ks, ke
281  qflx_sgs_momx(k,i,j,ydir) = -atmos_phy_tb_dns_nu * ( momx(k,i,j+1)-momx(k,i,j) ) * rfdy(j) * mapf(i,j,2,i_uy)
282  enddo
283  enddo
284  enddo
285 
286  !##### momentum equation (y) #####
287 
288  ! (x edge)
289  do j = jjs, jje
290  do i = iis, iie
291  do k = ks, ke-1
292  qflx_sgs_momy(k,i,j,zdir) = -atmos_phy_tb_dns_nu * ( momy(k+1,i,j)-momy(k,i,j) ) * rfdz(k)
293  enddo
294  enddo
295  enddo
296 
297  do j = jjs, jje
298  do i = iis, iie
299  qflx_sgs_momy(ks-1,i,j,zdir) = 0.0_rp ! bottom boundary
300  qflx_sgs_momy(ke ,i,j,zdir) = 0.0_rp ! top boundary
301  enddo
302  enddo
303 
304  ! (z edge)
305  do j = jjs, jje
306  do i = iis-1, iie
307  do k = ks, ke
308  qflx_sgs_momy(k,i,j,xdir) = -atmos_phy_tb_dns_nu * ( momy(k,i+1,j)-momy(k,i,j) ) * rfdx(i) * mapf(i,j,1,i_xv)
309  enddo
310  enddo
311  enddo
312 
313  ! (z-x plane)
314  do j = jjs, jje+1
315  do i = iis, iie
316  do k = ks, ke
317  qflx_sgs_momy(k,i,j,ydir) = -atmos_phy_tb_dns_nu * ( momy(k,i,j)-momy(k,i,j-1) ) * rcdy(j) * mapf(i,j,2,i_xv)
318  enddo
319  enddo
320  enddo
321 
322  !##### Thermodynamic Equation #####
323 
324  ! at x, y ,w
325  do j = jjs, jje
326  do i = iis, iie
327  do k = ks, ke-1
328  qflx_sgs_rhot(k,i,j,zdir) = -0.5_rp * ( dens(k+1,i,j)+dens(k,i,j) ) &
329  * atmos_phy_tb_dns_mu * ( pott(k+1,i,j)-pott(k,i,j) ) * rfdz(k)
330  enddo
331  enddo
332  enddo
333 
334  do j = jjs, jje
335  do i = iis, iie
336  qflx_sgs_rhot(ks-1,i,j,zdir) = 0.0_rp
337  qflx_sgs_rhot(ke ,i,j,zdir) = 0.0_rp
338  enddo
339  enddo
340 
341  ! at u, y, z
342  do j = jjs, jje
343  do i = iis-1, iie
344  do k = ks, ke
345  qflx_sgs_rhot(k,i,j,xdir) = -0.5_rp * ( dens(k,i+1,j)+dens(k,i,j) ) &
346  * atmos_phy_tb_dns_mu * ( pott(k,i+1,j)-pott(k,i,j) ) * rfdx(i) * mapf(i,j,1,i_xy)
347  enddo
348  enddo
349  enddo
350 
351  ! at x, v, z
352  do j = jjs-1, jje
353  do i = iis, iie
354  do k = ks, ke
355  qflx_sgs_rhot(k,i,j,ydir) = -0.5_rp * ( dens(k,i,j+1)+dens(k,i,j) ) &
356  * atmos_phy_tb_dns_mu * ( pott(k,i,j+1)-pott(k,i,j) ) * rfdy(j) * mapf(i,j,2,i_xy)
357  enddo
358  enddo
359  enddo
360 
361  enddo
362  enddo
363 
364  !##### Tracers #####
365  do iq = 1, qa
366 
367  if ( .not. tracer_advc(iq) ) cycle
368 
369  do jjs = js, je, jblock
370  jje = jjs+jblock-1
371  do iis = is, ie, iblock
372  iie = iis+iblock-1
373 
374  ! at x, y ,w
375  do j = jjs, jje
376  do i = iis, iie
377  do k = ks, ke-1
378  qflx_sgs_rhoq(k,i,j,zdir,iq) = -0.5_rp * ( dens(k+1,i,j)+dens(k,i,j) ) &
379  * atmos_phy_tb_dns_mu * ( qtrc(k+1,i,j,iq)-qtrc(k,i,j,iq) ) * rfdz(k)
380  enddo
381  enddo
382  enddo
383  do j = jjs, jje
384  do i = iis, iie
385  qflx_sgs_rhoq(ks-1,i,j,zdir,iq) = 0.0_rp
386  qflx_sgs_rhoq(ke ,i,j,zdir,iq) = 0.0_rp
387  enddo
388  enddo
389 
390  ! at u, y, z
391  do j = jjs, jje
392  do i = iis-1, iie
393  do k = ks, ke
394  qflx_sgs_rhoq(k,i,j,xdir,iq) = -0.5_rp * ( dens(k,i+1,j)+dens(k,i,j) ) &
395  * atmos_phy_tb_dns_mu * ( qtrc(k,i+1,j,iq)-qtrc(k,i,j,iq) ) * rfdx(i) * mapf(i,j,1,i_xy)
396  enddo
397  enddo
398  enddo
399 
400  ! at x, v, z
401  do j = jjs-1, jje
402  do i = iis, iie
403  do k = ks, ke
404  qflx_sgs_rhoq(k,i,j,ydir,iq) = -0.5_rp * ( dens(k,i,j+1)+dens(k,i,j) ) &
405  * atmos_phy_tb_dns_mu * ( qtrc(k,i,j+1,iq)-qtrc(k,i,j,iq) ) * rfdy(j) * mapf(i,j,2,i_xy)
406  enddo
407  enddo
408  enddo
409 
410  enddo
411  enddo
412 
413  enddo ! scalar quantities loop
414 
415  return
416  end subroutine atmos_phy_tb_dns
417 
418 end module scale_atmos_phy_tb_dns
scale_const::const_grav
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:46
scale_atmos_grid_cartesc_index::i_uy
integer, public i_uy
Definition: scale_atmos_grid_cartesC_index.F90:99
scale_atmos_grid_cartesc_index::ke
integer, public ke
end point of inner domain: z, local
Definition: scale_atmos_grid_cartesC_index.F90:52
scale_atmos_grid_cartesc_index::i_xv
integer, public i_xv
Definition: scale_atmos_grid_cartesC_index.F90:100
scale_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:342
scale_atmos_phy_tb_dns::atmos_phy_tb_dns
subroutine, public atmos_phy_tb_dns(qflx_sgs_MOMZ, qflx_sgs_MOMX, qflx_sgs_MOMY, qflx_sgs_rhot, qflx_sgs_rhoq, RHOQ_t, nu, Ri, Pr, MOMZ, MOMX, MOMY, RHOT, DENS, QTRC, N2, SFLX_MW, SFLX_MU, SFLX_MV, SFLX_SH, SFLX_Q, GSQRT, J13G, J23G, J33G, MAPF, dt)
Definition: scale_atmos_phy_tb_dns.F90:130
scale_tracer::qa
integer, public qa
Definition: scale_tracer.F90:34
scale_atmos_phy_tb_dns
module ATMOSPHERE / Physics Turbulence
Definition: scale_atmos_phy_tb_dns.F90:12
scale_atmos_grid_cartesc::atmos_grid_cartesc_rfdx
real(rp), dimension(:), allocatable, public atmos_grid_cartesc_rfdx
reciprocal of face-dx
Definition: scale_atmos_grid_cartesC.F90:67
scale_const::const_undef2
integer, parameter, public const_undef2
undefined value (INT2)
Definition: scale_const.F90:38
scale_precision
module PRECISION
Definition: scale_precision.F90:14
scale_atmos_grid_cartesc_index::ka
integer, public ka
Definition: scale_atmos_grid_cartesC_index.F90:47
scale_atmos_grid_cartesc_index::jblock
integer, public jblock
block size for cache blocking: y
Definition: scale_atmos_grid_cartesC_index.F90:41
scale_atmos_grid_cartesc::atmos_grid_cartesc_rcdx
real(rp), dimension(:), allocatable, public atmos_grid_cartesc_rcdx
reciprocal of center-dx
Definition: scale_atmos_grid_cartesC.F90:65
scale_tracer::tracer_advc
logical, dimension(qa_max), public tracer_advc
Definition: scale_tracer.F90:45
scale_atmos_grid_cartesc_index::iblock
integer, public iblock
block size for cache blocking: x
Definition: scale_atmos_grid_cartesC_index.F90:40
scale_atmos_phy_tb_dns::atmos_phy_tb_dns_setup
subroutine, public atmos_phy_tb_dns_setup(CDZ, CDX, CDY, CZ)
Definition: scale_atmos_phy_tb_dns.F90:88
scale_atmos_phy_tb_dns::atmos_phy_tb_dns_config
subroutine, public atmos_phy_tb_dns_config(TYPE_TB, I_TKE_out)
Config.
Definition: scale_atmos_phy_tb_dns.F90:64
scale_atmos_grid_cartesc::atmos_grid_cartesc_rcdz
real(rp), dimension(:), allocatable, public atmos_grid_cartesc_rcdz
reciprocal of center-dz
Definition: scale_atmos_grid_cartesC.F90:44
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_atmos_grid_cartesc_index::i_xy
integer, public i_xy
Definition: scale_atmos_grid_cartesC_index.F90:98
scale_atmos_grid_cartesc_index::ie
integer, public ie
end point of inner domain: x, local
Definition: scale_atmos_grid_cartesC_index.F90:54
scale_io
module STDIO
Definition: scale_io.F90:10
scale_atmos_grid_cartesc::atmos_grid_cartesc_rcdy
real(rp), dimension(:), allocatable, public atmos_grid_cartesc_rcdy
reciprocal of center-dy
Definition: scale_atmos_grid_cartesC.F90:66
scale_tracer::k
real(rp), public k
Definition: scale_tracer.F90:44
scale_atmos_grid_cartesc_index
module atmosphere / grid / cartesC index
Definition: scale_atmos_grid_cartesC_index.F90:12
scale_const
module CONSTANT
Definition: scale_const.F90:11
scale_atmos_grid_cartesc_index::ia
integer, public ia
Definition: scale_atmos_grid_cartesC_index.F90:48
scale_debug::check
subroutine, public check(current_line, v)
Undefined value checker.
Definition: scale_debug.F90:56
scale_atmos_grid_cartesc_index::zdir
integer, parameter, public zdir
Definition: scale_atmos_grid_cartesC_index.F90:32
scale_prof
module profiler
Definition: scale_prof.F90:11
scale_atmos_grid_cartesc_index::is
integer, public is
start point of inner domain: x, local
Definition: scale_atmos_grid_cartesC_index.F90:53
scale_precision::dp
integer, parameter, public dp
Definition: scale_precision.F90:32
scale_atmos_grid_cartesc_index::ja
integer, public ja
Definition: scale_atmos_grid_cartesC_index.F90:49
scale_tracer
module TRACER
Definition: scale_tracer.F90:12
scale_atmos_grid_cartesc_index::xdir
integer, parameter, public xdir
Definition: scale_atmos_grid_cartesC_index.F90:33
scale_atmos_grid_cartesc_index::ks
integer, public ks
start point of inner domain: z, local
Definition: scale_atmos_grid_cartesC_index.F90:51
scale_debug
module DEBUG
Definition: scale_debug.F90:11
scale_atmos_grid_cartesc_index::js
integer, public js
start point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:55
scale_atmos_grid_cartesc_index::ydir
integer, parameter, public ydir
Definition: scale_atmos_grid_cartesC_index.F90:34
scale_atmos_grid_cartesc::atmos_grid_cartesc_rfdz
real(rp), dimension(:), allocatable, public atmos_grid_cartesc_rfdz
reciprocal of face-dz
Definition: scale_atmos_grid_cartesC.F90:45
scale_atmos_grid_cartesc
module atmosphere / grid / cartesC
Definition: scale_atmos_grid_cartesC.F90:12
scale_const::const_undef
real(rp), public const_undef
Definition: scale_const.F90:41
scale_atmos_grid_cartesc::atmos_grid_cartesc_rfdy
real(rp), dimension(:), allocatable, public atmos_grid_cartesc_rfdy
reciprocal of face-dy
Definition: scale_atmos_grid_cartesC.F90:68
scale_io::io_fid_conf
integer, public io_fid_conf
Config file ID.
Definition: scale_io.F90:56
scale_atmos_grid_cartesc_index::je
integer, public je
end point of inner domain: y, local
Definition: scale_atmos_grid_cartesC_index.F90:56