SCALE-RM
scale_atmos_phy_mp_kessler.F90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
18 !-------------------------------------------------------------------------------
19 #include "inc_openmp.h"
21  !-----------------------------------------------------------------------------
22  !
23  !++ used modules
24  !
25  use scale_precision
26  use scale_stdio
27  use scale_prof
29 
30  use scale_atmos_hydrometeor, only: &
31  n_hyd, &
32  i_qv, &
33  i_qc, &
34  i_qr, &
35  i_hc, &
36  i_hr
37  !-----------------------------------------------------------------------------
38  implicit none
39  private
40  !-----------------------------------------------------------------------------
41  !
42  !++ Public procedure
43  !
46  public :: atmos_phy_mp_kessler
50 
51  !-----------------------------------------------------------------------------
52  !
53  !++ Public parameters & variables
54  !
55  integer, public, parameter :: qa_mp = 3
56 
57  character(len=H_SHORT), public, target :: atmos_phy_mp_kessler_name(qa_mp)
58  character(len=H_MID) , public, target :: atmos_phy_mp_kessler_desc(qa_mp)
59  character(len=H_SHORT), public, target :: atmos_phy_mp_kessler_unit(qa_mp)
60 
61  real(RP), public, target :: atmos_phy_mp_kessler_dens(n_hyd) ! hydrometeor density [kg/m3]=[g/L]
62 
64  'QV', &
65  'QC', &
66  'QR' /
67 
69  'Ratio of Water Vapor mass to total mass (Specific humidity)', &
70  'Ratio of Cloud Water mass to total mass', &
71  'Ratio of Rain Water mass to total mass' /
72 
74  'kg/kg', &
75  'kg/kg', &
76  'kg/kg' /
77 
78  !-----------------------------------------------------------------------------
79  !
80  !++ Private procedure
81  !
82  private :: mp_kessler
83  private :: mp_kessler_vterm
84 
85  !-----------------------------------------------------------------------------
86  !
87  !++ Private parameters & variables
88  !
89  integer, private, parameter :: i_mp_qc = 1
90  integer, private, parameter :: i_mp_qr = 2
91  integer, private :: qs_mp
92  integer, private :: qe_mp
93 
94  logical, private :: mp_donegative_fixer = .true. ! apply negative fixer?
95  logical, private :: mp_doprecipitation = .true. ! apply sedimentation (precipitation)?
96  logical, private :: mp_couple_aerosol = .false. ! Consider CCN effect ?
97  real(RP), private :: mp_limit_negative = 1.0_rp ! Abort if abs(fixed negative vaue) > abs(MP_limit_negative)
98 
99  real(RP), private, allocatable :: factor_vterm(:) ! collection factor for terminal velocity of QR
100 
101  real(RP), private, parameter :: re_qc = 8.e-6_rp ! effective radius for cloud water
102 
103  integer, private :: mp_ntmax_sedimentation = 1 ! number of time step for sedimentation
104  integer, private :: mp_nstep_sedimentation
105  real(RP), private :: mp_rnstep_sedimentation
106  real(DP), private :: mp_dtsec_sedimentation
107 
108  logical, private :: first = .true.
109 
110  !-----------------------------------------------------------------------------
111 contains
112  !-----------------------------------------------------------------------------
114  subroutine atmos_phy_mp_kessler_config( &
115  MP_TYPE, &
116  QA, QS )
117  use scale_process, only: &
119  use scale_atmos_hydrometeor, only: &
121  implicit none
122 
123  character(len=*), intent(in) :: mp_type
124  integer, intent(out) :: qa
125  integer, intent(out) :: qs
126  !---------------------------------------------------------------------------
127 
128  if( io_l ) write(io_fid_log,*)
129  if( io_l ) write(io_fid_log,*) '++++++ Module[Cloud Microphysics Tracer] / Categ[ATMOS PHYSICS] / Origin[SCALElib]'
130  if( io_l ) write(io_fid_log,*) '*** Tracers for KESSLER-type 1-moment bulk 3 category'
131 
132  if ( mp_type /= 'KESSLER' ) then
133  write(*,*) 'xxx ATMOS_PHY_MP_TYPE is not KESSLER. Check!'
134  call prc_mpistop
135  endif
136 
137  call atmos_hydrometeor_regist( qs_mp, & ! [OUT]
138  1, 2, 0, & ! [IN]
139  atmos_phy_mp_kessler_name, & ! [IN]
140  atmos_phy_mp_kessler_desc, & ! [IN]
142 
143  qa = qa_mp
144  qs = qs_mp
145  qe_mp = qs_mp + qa_mp - 1
146 
147  i_qv = qs
148  i_qc = qs + i_mp_qc
149  i_qr = qs + i_mp_qr
150 
151  return
152  end subroutine atmos_phy_mp_kessler_config
153 
154  !-----------------------------------------------------------------------------
156  subroutine atmos_phy_mp_kessler_setup
157  use scale_process, only: &
159  use scale_const, only: &
160  const_undef, &
162  use scale_time, only: &
164  use scale_grid, only: &
165  cdz => grid_cdz
166  implicit none
167 
168  real(RP), parameter :: max_term_vel = 10.0_rp !-- terminal velocity for calculate dt of sedimentation
169  integer :: nstep_max
170 
171  namelist / param_atmos_phy_mp / &
172  mp_doprecipitation, &
173  mp_donegative_fixer, &
174  mp_limit_negative, &
175  mp_ntmax_sedimentation, &
176  mp_couple_aerosol
177 
178  integer :: ierr
179 
180  !---------------------------------------------------------------------------
181 
182  if( io_l ) write(io_fid_log,*)
183  if( io_l ) write(io_fid_log,*) '++++++ Module[Cloud Microphysics] / Categ[ATMOS PHYSICS] / Origin[SCALElib]'
184  if( io_l ) write(io_fid_log,*) '*** KESSLER-type 1-moment bulk 3 category'
185 
186  !--- read namelist
187  rewind(io_fid_conf)
188  read(io_fid_conf,nml=param_atmos_phy_mp,iostat=ierr)
189  if( ierr < 0 ) then !--- missing
190  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
191  elseif( ierr > 0 ) then !--- fatal error
192  write(*,*) 'xxx Not appropriate names in namelist PARAM_ATMOS_PHY_MP. Check!'
193  call prc_mpistop
194  endif
195  if( io_nml ) write(io_fid_nml,nml=param_atmos_phy_mp)
196 
197  if( mp_couple_aerosol ) then
198  write(*,*) 'xxx MP_aerosol_couple should be .false. for KESSLER type MP!'
199  call prc_mpistop
200  endif
201 
202  if( io_l ) write(io_fid_log,*)
203  if( io_l ) write(io_fid_log,*) '*** Enable negative fixer? : ', mp_donegative_fixer
204  if( io_l ) write(io_fid_log,*) '*** Value limit of negative fixer (abs) : ', abs(mp_limit_negative)
205  if( io_l ) write(io_fid_log,*) '*** Enable sedimentation (precipitation)? : ', mp_doprecipitation
206 
207  nstep_max = int( ( time_dtsec_atmos_phy_mp * max_term_vel ) / minval( cdz ) )
208  mp_ntmax_sedimentation = max( mp_ntmax_sedimentation, nstep_max )
209 
210  mp_nstep_sedimentation = mp_ntmax_sedimentation
211  mp_rnstep_sedimentation = 1.0_rp / real(mp_ntmax_sedimentation,kind=rp)
212  mp_dtsec_sedimentation = time_dtsec_atmos_phy_mp * mp_rnstep_sedimentation
213 
214  if( io_l ) write(io_fid_log,*) '*** Timestep of sedimentation is divided into : ', mp_ntmax_sedimentation, 'step'
215  if( io_l ) write(io_fid_log,*) '*** DT of sedimentation : ', mp_dtsec_sedimentation, '[s]'
216 
220 
221  allocate( factor_vterm(ka) )
222 
223  return
224  end subroutine atmos_phy_mp_kessler_setup
225 
226  !-----------------------------------------------------------------------------
228  subroutine atmos_phy_mp_kessler( &
229  DENS, &
230  MOMZ, &
231  MOMX, &
232  MOMY, &
233  RHOT, &
234  QTRC, &
235  CCN, &
236  EVAPORATE, &
237  SFLX_rain, &
238  SFLX_snow )
240  use scale_const, only: &
241  dwatr => const_dwatr, &
242  pi => const_pi
243  use scale_comm, only: &
245  use scale_time, only: &
246  dt_mp => time_dtsec_atmos_phy_mp
247  use scale_history, only: &
248  hist_in
249  use scale_tracer, only: &
250  qa, &
251  tracer_r, &
252  tracer_cv, &
254  use scale_atmos_thermodyn, only: &
255  thermodyn_rhoe => atmos_thermodyn_rhoe, &
256  thermodyn_rhot => atmos_thermodyn_rhot, &
257  thermodyn_temp_pres_e => atmos_thermodyn_temp_pres_e
258  use scale_atmos_phy_mp_common, only: &
259  mp_negative_fixer => atmos_phy_mp_negative_fixer, &
260  mp_precipitation => atmos_phy_mp_precipitation, &
261  mp_saturation_adjustment => atmos_phy_mp_saturation_adjustment
262  implicit none
263 
264  real(RP), intent(inout) :: dens (ka,ia,ja)
265  real(RP), intent(inout) :: momz (ka,ia,ja)
266  real(RP), intent(inout) :: momx (ka,ia,ja)
267  real(RP), intent(inout) :: momy (ka,ia,ja)
268  real(RP), intent(inout) :: rhot (ka,ia,ja)
269  real(RP), intent(inout) :: qtrc (ka,ia,ja,qa)
270  real(RP), intent(in) :: ccn (ka,ia,ja)
271  real(RP), intent(out) :: evaporate(ka,ia,ja) !--- number of evaporated cloud [/m3]
272  real(RP), intent(out) :: sflx_rain(ia,ja)
273  real(RP), intent(out) :: sflx_snow(ia,ja)
274 
275  real(RP) :: rhoe_t (ka,ia,ja)
276  real(RP) :: qtrc_t (ka,ia,ja,qa)
277  real(RP) :: rhoe (ka,ia,ja)
278  real(RP) :: temp (ka,ia,ja)
279  real(RP) :: pres (ka,ia,ja)
280 
281  real(RP) :: vterm (ka,ia,ja,qa_mp-1) ! terminal velocity of each tracer [m/s]
282  real(RP) :: pflux (ka,ia,ja,qa_mp-1) ! precipitation flux of each tracer [kg/m2/s]
283  real(RP) :: flx_hydro(ka,ia,ja,qa_mp-1)
284  real(RP) :: qc_t_sat (ka,ia,ja)
285  real(RP) :: rho_prof (ka) ! averaged profile of rho
286 
287  integer :: step
288  integer :: k, i, j
289  !---------------------------------------------------------------------------
290 
291  if( io_l ) write(io_fid_log,*) '*** Atmos physics step: Cloud microphysics(kessler)'
292 
293  if ( first ) then
294  ! Calculate collection factor for terminal velocity of QR
295  call comm_horizontal_mean( rho_prof(:), dens(:,:,:) )
296  rho_prof(:) = rho_prof(:) * 1.e-3_rp ! [kg/m3]->[g/cc]
297 
298  do k = ks, ke
299  factor_vterm(k) = sqrt( rho_prof(ks)/rho_prof(k) )
300  enddo
301  first = .false.
302  endif
303 
304  if ( mp_donegative_fixer ) then
305  call mp_negative_fixer( dens(:,:,:), & ! [INOUT]
306  rhot(:,:,:), & ! [INOUT]
307  qtrc(:,:,:,:), & ! [INOUT]
308  i_qv, & ! [IN]
309  mp_limit_negative ) ! [IN]
310  endif
311 
312  call thermodyn_rhoe( rhoe(:,:,:), & ! [OUT]
313  rhot(:,:,:), & ! [IN]
314  qtrc(:,:,:,:), & ! [IN]
315  tracer_cv(:), & ! [IN]
316  tracer_r(:), & ! [IN]
317  tracer_mass(:) ) ! [IN]
318 
319  !##### MP Main #####
320  call mp_kessler( rhoe_t(:,:,:), & ! [OUT]
321  qtrc_t(:,:,:,:), & ! [OUT]
322  rhoe(:,:,:), & ! [INOUT]
323  qtrc(:,:,:,:), & ! [INOUT]
324  dens(:,:,:) ) ! [IN]
325 
326  vterm(:,:,:,:) = 0.0_rp
327  flx_hydro(:,:,:,:) = 0.0_rp
328 
329  if ( mp_doprecipitation ) then
330 
331  do step = 1, mp_nstep_sedimentation
332 
333  call thermodyn_temp_pres_e( temp(:,:,:), & ! [OUT]
334  pres(:,:,:), & ! [OUT]
335  dens(:,:,:), & ! [IN]
336  rhoe(:,:,:), & ! [IN]
337  qtrc(:,:,:,:), & ! [IN]
338  tracer_cv(:), & ! [IN]
339  tracer_r(:), & ! [IN]
340  tracer_mass(:) ) ! [IN]
341 
342  call mp_kessler_vterm( vterm(:,:,:,:), & ! [OUT]
343  dens(:,:,:), & ! [IN]
344  qtrc(:,:,:,:) ) ! [IN]
345 
346  call mp_precipitation( qa_mp, & ! [IN]
347  qs_mp, & ! [IN]
348  pflux(:,:,:,:), & ! [OUT]
349  vterm(:,:,:,:), & ! [INOUT]
350  dens(:,:,:), & ! [INOUT]
351  momz(:,:,:), & ! [INOUT]
352  momx(:,:,:), & ! [INOUT]
353  momy(:,:,:), & ! [INOUT]
354  rhoe(:,:,:), & ! [INOUT]
355  qtrc(:,:,:,:), & ! [INOUT]
356  temp(:,:,:), & ! [IN]
357  tracer_cv(:), & ! [IN]
358  mp_dtsec_sedimentation ) ! [IN]
359 
360  do j = js, je
361  do i = is, ie
362  do k = ks-1, ke-1
363  flx_hydro(k,i,j,i_mp_qr) = flx_hydro(k,i,j,i_mp_qr) + pflux(k,i,j,i_mp_qr) * mp_rnstep_sedimentation
364  enddo
365  enddo
366  enddo
367 
368  enddo
369 
370  endif
371 
372  call hist_in( flx_hydro(:,:,:,i_mp_qr), 'pflux_QR', 'precipitation flux of QR', 'kg/m2/s', nohalo=.true. )
373 
374  call hist_in( vterm(:,:,:,i_mp_qr), 'Vterm_QR', 'terminal velocity of QR', 'm/s' )
375 
376  do j = js, je
377  do i = is, ie
378  do k = ks, ke
379  qc_t_sat(k,i,j) = qtrc(k,i,j,i_qc)
380  enddo
381  enddo
382  enddo
383 
384  call mp_saturation_adjustment( rhoe_t(:,:,:), & ! [INOUT]
385  qtrc_t(:,:,:,:), & ! [INOUT]
386  rhoe(:,:,:), & ! [INOUT]
387  qtrc(:,:,:,:), & ! [INOUT]
388  dens(:,:,:), & ! [IN]
389  i_qv, & ! [IN]
390  i_qc, & ! [IN]
391  -1, & ! [IN]
392  flag_liquid = .true. ) ! [IN]
393 
394  do j = js, je
395  do i = is, ie
396  do k = ks, ke
397  qc_t_sat(k,i,j) = ( qtrc(k,i,j,i_qc) - qc_t_sat(k,i,j) ) / dt_mp
398  enddo
399  enddo
400  enddo
401 
402  call hist_in( qc_t_sat(:,:,:), 'Pcsat', 'QC production term by satadjust', 'kg/kg/s' )
403 
404  do j = js, je
405  do i = is, ie
406  do k = ks, ke
407  evaporate(k,i,j) = max( -qc_t_sat(k,i,j), 0.0_rp ) ! if negative, condensation
408  evaporate(k,i,j) = evaporate(k,i,j) * dens(k,i,j) / (4.0_rp/3.0_rp*pi*dwatr*re_qc**3) ! mass -> number (assuming constant particle radius as re_qc)
409  enddo
410  enddo
411  enddo
412 
413  !##### END MP Main #####
414 
415  call thermodyn_rhot( rhot(:,:,:), & ! [OUT]
416  rhoe(:,:,:), & ! [IN]
417  qtrc(:,:,:,:), & ! [IN]
418  tracer_cv(:), & ! [IN]
419  tracer_r(:), & ! [IN]
420  tracer_mass(:) ) ! [IN]
421 
422  if ( mp_donegative_fixer ) then
423  call mp_negative_fixer( dens(:,:,:), & ! [INOUT]
424  rhot(:,:,:), & ! [INOUT]
425  qtrc(:,:,:,:), & ! [INOUT]
426  i_qv, & ! [IN]
427  mp_limit_negative ) ! [IN]
428  endif
429 
430  do j = js, je
431  do i = is, ie
432  sflx_rain(i,j) = - flx_hydro(ks-1,i,j,i_mp_qr)
433  sflx_snow(i,j) = 0.0_rp
434  enddo
435  enddo
436 
437  return
438  end subroutine atmos_phy_mp_kessler
439 
440  !-----------------------------------------------------------------------------
442  subroutine mp_kessler( &
443  RHOE_t, &
444  QTRC_t, &
445  RHOE0, &
446  QTRC0, &
447  DENS0 )
448  use scale_tracer, only: &
449  qa, &
450  tracer_r, &
451  tracer_cv, &
453  use scale_time, only: &
455  use scale_atmos_thermodyn, only: &
456  thermodyn_temp_pres_e => atmos_thermodyn_temp_pres_e
457  use scale_atmos_hydrometeor, only: &
458  lhv
459  use scale_atmos_saturation, only: &
460  saturation_dens2qsat_liq => atmos_saturation_dens2qsat_liq
461  implicit none
462 
463  real(RP), intent(out) :: rhoe_t(ka,ia,ja) ! tendency rhoe [J/m3/s]
464  real(RP), intent(out) :: qtrc_t(ka,ia,ja,qa) ! tendency tracer [kg/kg/s]
465  real(RP), intent(inout) :: rhoe0 (ka,ia,ja) ! density * internal energy [J/m3]
466  real(RP), intent(inout) :: qtrc0 (ka,ia,ja,qa) ! mass concentration [kg/kg]
467  real(RP), intent(in) :: dens0 (ka,ia,ja) ! density [kg/m3]
468 
469  ! working
470  real(RP) :: qsat (ka,ia,ja) ! saturated water vapor [kg/kg]
471  real(RP) :: temp0(ka,ia,ja) ! temperature [K]
472  real(RP) :: pres0(ka,ia,ja) ! pressure [Pa]
473 
474  real(RP) :: dens (ka)
475  real(RP) :: rhoe (ka)
476  real(RP) :: temp (ka)
477  real(RP) :: pres (ka)
478  real(RP) :: qv (ka)
479  real(RP) :: qc (ka)
480  real(RP) :: qr (ka)
481  real(RP) :: qsatl(ka)
482  real(RP) :: sliq (ka)
483 
484  ! tendency
485  real(RP) :: dq_evap(ka) ! tendency q (evaporation)
486  real(RP) :: dq_auto(ka) ! tendency q (autoconversion)
487  real(RP) :: dq_accr(ka) ! tendency q (accretion)
488 
489  real(RP) :: dqv, dqc, dqr
490  real(RP) :: vent_factor
491  real(RP) :: rdt
492 
493  integer :: k, i, j
494  !---------------------------------------------------------------------------
495 
496  call prof_rapstart('MP_kessler', 3)
497 
498  rdt = 1.0_rp / dt
499 
500  call thermodyn_temp_pres_e( temp0(:,:,:), & ! [OUT]
501  pres0(:,:,:), & ! [OUT]
502  dens0(:,:,:), & ! [IN]
503  rhoe0(:,:,:), & ! [IN]
504  qtrc0(:,:,:,:), & ! [IN]
505  tracer_cv(:), & ! [IN]
506  tracer_r(:), & ! [IN]
507  tracer_mass(:) ) ! [IN]
508 
509  call saturation_dens2qsat_liq( qsat(:,:,:), & ! [OUT]
510  temp0(:,:,:), & ! [IN]
511  dens0(:,:,:) ) ! [IN]
512 
513  !$omp parallel do private(i,j,k,dens,rhoe,temp,pres,qv,qc,qr,qsatl,Sliq,dq_evap,dq_auto,dq_accr,vent_factor,dqc,dqr,dqv) OMP_SCHEDULE_ collapse(2)
514  do j = js, je
515  do i = is, ie
516 
517  ! store to work
518  dens(ks:ke) = dens0(ks:ke,i,j)
519  rhoe(ks:ke) = rhoe0(ks:ke,i,j)
520  temp(ks:ke) = temp0(ks:ke,i,j)
521  pres(ks:ke) = pres0(ks:ke,i,j)
522  qv(ks:ke) = qtrc0(ks:ke,i,j,i_qv)
523  qc(ks:ke) = qtrc0(ks:ke,i,j,i_qc)
524  qr(ks:ke) = qtrc0(ks:ke,i,j,i_qr)
525  qsatl(ks:ke) = qsat(ks:ke,i,j)
526  sliq(ks:ke) = qv(ks:ke) / qsatl(ks:ke)
527 
528  ! Auto-conversion (QC->QR)
529  do k = ks, ke
530  dq_auto(k) = 1.e-3_rp * max( qc(k)-1.e-3_rp, 0.0_rp )
531  enddo
532 
533  ! Accretion (QC->QR)
534  do k = ks, ke
535  dq_accr(k) = 2.2_rp * qc(k) * qr(k)**0.875_rp
536  enddo
537 
538  ! Evaporation (QR->QV)
539  do k = ks, ke
540  vent_factor = 1.6_rp + 124.9_rp * ( dens(k)*qr(k) )**0.2046_rp
541 
542  dq_evap(k) = ( 1.0_rp-min(sliq(k),1.0_rp) ) / dens(k) * vent_factor &
543  * ( dens(k)*qr(k) )**0.525_rp / ( 5.4e5_rp + 2.55e8_rp / ( pres(k)*qsatl(k) ) )
544  enddo
545 
546  ! limiter
547  do k = ks, ke
548  dqc = (-dq_auto(k)-dq_accr(k) )
549 
550  qtrc_t(k,i,j,i_qc) = max( dqc, -qc(k)*rdt )
551  enddo
552 
553  do k = ks, ke
554  dqr = ( dq_auto(k)+dq_accr(k)-dq_evap(k) )
555 
556  qtrc_t(k,i,j,i_qr) = max( dqr, -qr(k)*rdt )
557  enddo
558 
559  do k = ks, ke
560  dqv = - ( qtrc_t(k,i,j,i_qc) &
561  + qtrc_t(k,i,j,i_qr) )
562 
563  qtrc_t(k,i,j,i_qv) = max( dqv, -qv(k)*rdt )
564  enddo
565 
566  do k = ks, ke
567  rhoe_t(k,i,j) = - dens(k) * ( lhv * qtrc_t(k,i,j,i_qv) )
568  enddo
569  enddo
570  enddo
571 
572  ! mass & energy update
573  !$omp parallel do private(i,j, k) OMP_SCHEDULE_ collapse(2)
574  do j = js, je
575  do i = is, ie
576  qtrc0(ks:ke,i,j,i_qv) = qtrc0(ks:ke,i,j,i_qv) + qtrc_t(ks:ke,i,j,i_qv) * dt
577  qtrc0(ks:ke,i,j,i_qc) = qtrc0(ks:ke,i,j,i_qc) + qtrc_t(ks:ke,i,j,i_qc) * dt
578  qtrc0(ks:ke,i,j,i_qr) = qtrc0(ks:ke,i,j,i_qr) + qtrc_t(ks:ke,i,j,i_qr) * dt
579 
580  rhoe0(ks:ke,i,j) = rhoe0(ks:ke,i,j) + rhoe_t(ks:ke,i,j) * dt
581  enddo
582  enddo
583 
584  call prof_rapend ('MP_kessler', 3)
585 
586  return
587  end subroutine mp_kessler
588 
589  !-----------------------------------------------------------------------------
591  subroutine mp_kessler_vterm( &
592  vterm, &
593  DENS0, &
594  QTRC0 )
595  use scale_tracer, only: &
596  qa
597  use scale_atmos_refstate, only: &
598  refstate_dens => atmos_refstate_dens
599  implicit none
600 
601  real(RP), intent(out) :: vterm(ka,ia,ja,qa_mp-1)
602  real(RP), intent(in) :: dens0(ka,ia,ja)
603  real(RP), intent(in) :: qtrc0(ka,ia,ja,qa)
604 
605  real(RP) :: zerosw
606 
607  integer :: k, i, j
608  !---------------------------------------------------------------------------
609 
610  !$omp parallel do private(i,j,k,zerosw) OMP_SCHEDULE_ collapse(2)
611  do j = js, je
612  do i = is, ie
613  do k = ks, ke
614  vterm(k,i,j,i_mp_qc) = 0.0_rp
615  enddo
616  enddo
617  enddo
618 
619  do j = js, je
620  do i = is, ie
621  do k = ks, ke
622  zerosw = 0.5_rp - sign(0.5_rp, qtrc0(k,i,j,i_qr) - 1.e-12_rp )
623 
624  vterm(k,i,j,i_mp_qr) = - 36.34_rp * ( dens0(k,i,j) * ( qtrc0(k,i,j,i_qr) + zerosw ) )**0.1364_rp &
625  * refstate_dens(ks,i,j) / refstate_dens(k,i,j) * ( 1.0_rp - zerosw )
626  enddo
627  enddo
628  enddo
629 
630  return
631  end subroutine mp_kessler_vterm
632 
633  !-----------------------------------------------------------------------------
636  cldfrac, &
637  QTRC, &
638  mask_criterion )
640  use scale_tracer, only: &
641  qa
642  implicit none
643 
644  real(RP), intent(out) :: cldfrac(ka,ia,ja)
645  real(RP), intent(in) :: qtrc (ka,ia,ja,qa)
646  real(RP), intent(in) :: mask_criterion
647 
648  real(RP) :: qhydro
649  integer :: k, i, j, iq
650  !---------------------------------------------------------------------------
651 
652  do j = js, je
653  do i = is, ie
654  do k = ks, ke
655  qhydro = 0.0_rp
656  do iq = qs_mp+1, qe_mp
657  qhydro = qhydro + qtrc(k,i,j,iq)
658  enddo
659  cldfrac(k,i,j) = 0.5_rp + sign(0.5_rp,qhydro-mask_criterion)
660  enddo
661  enddo
662  enddo
663 
664  return
666 
667  !-----------------------------------------------------------------------------
670  Re, &
671  QTRC0, &
672  DENS0, &
673  TEMP0 )
675  use scale_tracer, only: &
676  qa
677  use scale_atmos_hydrometeor, only: &
678  n_hyd
679  implicit none
680  real(RP), intent(out) :: re (ka,ia,ja,n_hyd) ! effective radius [cm]
681  real(RP), intent(in) :: qtrc0(ka,ia,ja,qa) ! tracer mass concentration [kg/kg]
682  real(RP), intent(in) :: dens0(ka,ia,ja) ! density [kg/m3]
683  real(RP), intent(in) :: temp0(ka,ia,ja) ! temperature [K]
684 
685  real(RP), parameter :: um2cm = 100.0_rp
686 
687  integer :: ih
688  !---------------------------------------------------------------------------
689 
690  do ih = 1, n_hyd
691  if ( ih == i_hc ) then
692  re(:,:,:,ih) = 8.e-6_rp * um2cm
693  else if ( ih == i_hr ) then
694  re(:,:,:,ih) = 100.e-6_rp * um2cm
695  else
696  re(:,:,:,ih) = 0.0_rp
697  end if
698  end do
699 
700  return
702 
703  !-----------------------------------------------------------------------------
706  Qe, &
707  QTRC0 )
709  use scale_tracer, only: &
710  qa
711  use scale_atmos_hydrometeor, only: &
712  n_hyd
713  implicit none
714  real(RP), intent(out) :: qe (ka,ia,ja,n_hyd) ! mixing ratio of each cateory [kg/kg]
715  real(RP), intent(in) :: qtrc0(ka,ia,ja,qa) ! tracer mass concentration [kg/kg]
716 
717  integer :: ih
718  !---------------------------------------------------------------------------
719 
720  do ih = 1, n_hyd
721  if ( ih == i_hc ) then
722  qe(:,:,:,ih) = qtrc0(:,:,:,i_qc)
723  else if ( ih == i_hr ) then
724  qe(:,:,:,ih) = qtrc0(:,:,:,i_qr)
725  else
726  qe(:,:,:,ih) = 0.0_rp
727  end if
728  end do
729 
730  return
731  end subroutine atmos_phy_mp_kessler_mixingratio
732 
integer, public is
start point of inner domain: x, local
subroutine, public atmos_phy_mp_kessler(DENS, MOMZ, MOMX, MOMY, RHOT, QTRC, CCN, EVAPORATE, SFLX_rain, SFLX_snow)
Cloud Microphysics.
integer, public je
end point of inner domain: y, local
module ATMOSPHERE / Saturation adjustment
subroutine, public prc_mpistop
Abort MPI.
subroutine, public atmos_phy_mp_kessler_mixingratio(Qe, QTRC0)
Calculate mixing ratio of each category.
real(rp), parameter, public const_dwatr
density of water [kg/m3]
Definition: scale_const.F90:84
subroutine, public atmos_phy_mp_negative_fixer(DENS, RHOT, QTRC, I_QV, limit_negative)
Negative fixer.
real(rp), dimension(qa_max), public tracer_r
real(dp), public time_dtsec_atmos_phy_mp
time interval of physics(microphysics) [sec]
Definition: scale_time.F90:41
logical, public io_l
output log or not? (this process)
Definition: scale_stdio.F90:61
module ATMOSPHERE / Reference state
integer, parameter, public i_hr
liquid water rain
subroutine, public atmos_phy_mp_kessler_cloudfraction(cldfrac, QTRC, mask_criterion)
Calculate Cloud Fraction.
module STDIO
Definition: scale_stdio.F90:12
integer, public ke
end point of inner domain: z, local
character(len=h_short), dimension(qa_mp), target, public atmos_phy_mp_kessler_name
integer, public qa
real(rp), dimension(qa_max), public tracer_cv
subroutine, public atmos_phy_mp_kessler_effectiveradius(Re, QTRC0, DENS0, TEMP0)
Calculate Effective Radius.
real(rp), public const_undef
Definition: scale_const.F90:43
subroutine, public atmos_phy_mp_kessler_config(MP_TYPE, QA, QS)
Setup.
module ATMOSPHERE / Physics Cloud Microphysics - Common
module grid index
character(len=h_short), dimension(qa_mp), target, public atmos_phy_mp_kessler_unit
logical, public io_nml
output log or not? (for namelist, this process)
Definition: scale_stdio.F90:62
module TRACER
integer, public ia
of whole cells: x, local, with HALO
integer, parameter, public qa_mp
subroutine, public comm_horizontal_mean(varmean, var)
calculate horizontal mean (global total with communication)
Definition: scale_comm.F90:482
module ATMOSPHERE / Physics Cloud Microphysics
real(rp), dimension(n_hyd), target, public atmos_phy_mp_kessler_dens
subroutine, public atmos_phy_mp_saturation_adjustment(RHOE_t, QTRC_t, RHOE0, QTRC0, DENS0, I_QV, I_QC, I_QI, flag_liquid)
Saturation adjustment.
integer, public ka
of whole cells: z, local, with HALO
real(rp), dimension(:,:,:), allocatable, public atmos_refstate_dens
refernce density [kg/m3]
subroutine, public atmos_phy_mp_kessler_setup
Setup.
module COMMUNICATION
Definition: scale_comm.F90:23
integer, public js
start point of inner domain: y, local
module TIME
Definition: scale_time.F90:15
module PROCESS
real(rp), public lhv
latent heat of vaporization for use [J/kg]
module CONSTANT
Definition: scale_const.F90:14
integer, parameter, public i_hc
liquid water cloud
integer, public ks
start point of inner domain: z, local
module GRID (cartesian)
subroutine, public prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:156
module profiler
Definition: scale_prof.F90:10
integer, public ie
end point of inner domain: x, local
module ATMOSPHERE / Thermodynamics
character(len=h_mid), dimension(qa_mp), target, public atmos_phy_mp_kessler_desc
module PRECISION
real(rp), dimension(:), allocatable, public grid_cdz
z-length of control volume [m]
subroutine, public atmos_hydrometeor_regist(Q0, NV, NL, NI, NAME, DESC, UNIT, ADVC)
Regist tracer.
module HISTORY
real(rp), public const_pi
pi
Definition: scale_const.F90:34
subroutine, public atmos_phy_mp_precipitation(QA_MP, QS_MP, qflx, vterm, DENS, MOMZ, MOMX, MOMY, RHOE, QTRC, temp, CVq, dt, vt_fixed)
precipitation transport
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, parameter, public n_hyd
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:204
integer, parameter, public rp
integer, public io_fid_nml
Log file ID (only for output namelist)
Definition: scale_stdio.F90:57
real(rp), dimension(qa_max), public tracer_mass
integer, public ja
of whole cells: y, local, with HALO