SCALE-RM
mod_mkinit.f90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
30 !-------------------------------------------------------------------------------
31 module mod_mkinit
32  !-----------------------------------------------------------------------------
33  !
34  !++ used modules
35  !
36  use scale_precision
37  use scale_stdio
38  use scale_prof
39  use scale_grid_index
40  use scale_tracer
41 
42  use scale_process, only: &
43  prc_mpistop
44  use scale_const, only: &
45  pi => const_pi, &
46  grav => const_grav, &
47  pstd => const_pstd, &
48  rdry => const_rdry, &
49  cpdry => const_cpdry, &
50  p00 => const_pre00, &
51  i_sw => const_i_sw, &
52  i_lw => const_i_lw
53  use scale_random, only: &
54  random_get
55  use scale_comm, only: &
56  comm_vars8, &
57  comm_wait
58  use scale_grid, only: &
59  grid_cz, &
60  grid_cx, &
61  grid_cy, &
62  grid_fz, &
63  grid_fx, &
64  grid_fy, &
65  grid_cxg, &
66  grid_fxg, &
67  grid_fyg
68  use scale_grid_real, only: &
69  real_cz, &
70  real_fz
71  use scale_atmos_profile, only: &
72  profile_isa => atmos_profile_isa
73  use scale_atmos_hydrometeor, only: &
74  hydrometeor_lhv => atmos_hydrometeor_lhv
75  use scale_atmos_hydrostatic, only: &
76  hydrostatic_buildrho => atmos_hydrostatic_buildrho, &
77  hydrostatic_buildrho_atmos => atmos_hydrostatic_buildrho_atmos, &
78  hydrostatic_buildrho_bytemp => atmos_hydrostatic_buildrho_bytemp
79  use scale_atmos_saturation, only: &
80  saturation_pres2qsat_all => atmos_saturation_pres2qsat_all, &
81  saturation_pres2qsat_liq => atmos_saturation_pres2qsat_liq
82  use mod_atmos_vars, only: &
83  dens, &
84  momx, &
85  momy, &
86  momz, &
87  rhot, &
88  qtrc
89  use mod_atmos_phy_ae_vars, only: &
90  ccn => atmos_phy_ae_ccn
91  !-----------------------------------------------------------------------------
92  implicit none
93  private
94  !-----------------------------------------------------------------------------
95  !
96  !++ Public procedure
97  !
98  public :: mkinit_setup
99  public :: mkinit
100 
101  !-----------------------------------------------------------------------------
102  !
103  !++ Public parameters & variables
104  !
105  integer, public :: mkinit_type = -1
106  integer, public, parameter :: i_ignore = 0
107 
108  integer, public, parameter :: i_planestate = 1
109  integer, public, parameter :: i_tracerbubble = 2
110  integer, public, parameter :: i_coldbubble = 3
111 
112  integer, public, parameter :: i_lambwave = 4
113  integer, public, parameter :: i_gravitywave = 5
114  integer, public, parameter :: i_khwave = 6
115  integer, public, parameter :: i_turbulence = 7
116  integer, public, parameter :: i_mountainwave = 8
117 
118  integer, public, parameter :: i_warmbubble = 9
119  integer, public, parameter :: i_supercell = 10
120  integer, public, parameter :: i_squallline = 11
121  integer, public, parameter :: i_wk1982 = 12
122  integer, public, parameter :: i_dycoms2_rf01 = 13
123  integer, public, parameter :: i_dycoms2_rf02 = 14
124  integer, public, parameter :: i_rico = 15
125 
126  integer, public, parameter :: i_interporation = 16
127 
128  integer, public, parameter :: i_landcouple = 17
129  integer, public, parameter :: i_oceancouple = 18
130  integer, public, parameter :: i_urbancouple = 19
131  integer, public, parameter :: i_triplecouple = 20
132  integer, public, parameter :: i_bubblecouple = 21
133 
134  integer, public, parameter :: i_seabreeze = 22
135  integer, public, parameter :: i_heatisland = 23
136 
137  integer, public, parameter :: i_dycoms2_rf02_dns = 24
138 
139  integer, public, parameter :: i_real = 25
140 
141  integer, public, parameter :: i_grayzone = 26
142  integer, public, parameter :: i_boxaero = 27
143  integer, public, parameter :: i_warmbubbleaero = 28
144 
145  integer, public, parameter :: i_cavityflow = 29
146  integer, public, parameter :: i_barocwave = 30
147  integer, public, parameter :: i_bomex = 31
148 
149  !-----------------------------------------------------------------------------
150  !
151  !++ Private procedure
152  !
153  private :: bubble_setup
154  private :: sbmaero_setup
155  private :: aerosol_setup
156 
157  private :: mkinit_planestate
158  private :: mkinit_tracerbubble
159  private :: mkinit_coldbubble
160  private :: mkinit_lambwave
161  private :: mkinit_gravitywave
162  private :: mkinit_khwave
163  private :: mkinit_turbulence
164  private :: mkinit_cavityflow
165  private :: mkinit_mountainwave
166  private :: mkinit_barocwave
167 
168  private :: mkinit_warmbubble
169  private :: mkinit_supercell
170  private :: mkinit_squallline
171  private :: mkinit_wk1982
172  private :: mkinit_dycoms2_rf01
173  private :: mkinit_dycoms2_rf02
174  private :: mkinit_rico
175  private :: mkinit_bomex
176 
177  private :: mkinit_interporation
178 
179  private :: mkinit_landcouple
180  private :: mkinit_oceancouple
181  private :: mkinit_urbancouple
182  private :: mkinit_seabreeze
183  private :: mkinit_heatisland
184 
185  private :: mkinit_dycoms2_rf02_dns
186 
187  private :: mkinit_real
188 
189  private :: mkinit_grayzone
190 
191  private :: mkinit_boxaero
192  private :: mkinit_warmbubbleaero
193 
194  !-----------------------------------------------------------------------------
195  !
196  !++ Private parameters & variables
197  !
198  real(RP), private, parameter :: thetastd = 300.0_rp ! [K]
199 
200  real(RP), private, allocatable :: pres (:,:,:) ! pressure [Pa]
201  real(RP), private, allocatable :: temp (:,:,:) ! temperature [K]
202  real(RP), private, allocatable :: pott (:,:,:) ! potential temperature [K]
203  real(RP), private, allocatable :: qsat (:,:,:) ! satulated water vapor [kg/kg]
204  real(RP), private, allocatable :: qv (:,:,:) ! water vapor [kg/kg]
205  real(RP), private, allocatable :: qc (:,:,:) ! cloud water [kg/kg]
206  real(RP), private, allocatable :: velx (:,:,:) ! velocity u [m/s]
207  real(RP), private, allocatable :: vely (:,:,:) ! velocity v [m/s]
208 
209  real(RP), private, allocatable :: pres_sfc(:,:,:) ! surface pressure [Pa]
210  real(RP), private, allocatable :: temp_sfc(:,:,:) ! surface temperature [K]
211  real(RP), private, allocatable :: pott_sfc(:,:,:) ! surface potential temperature [K]
212  real(RP), private, allocatable :: qsat_sfc(:,:,:) ! surface satulated water vapor [kg/kg]
213  real(RP), private, allocatable :: qv_sfc (:,:,:) ! surface water vapor [kg/kg]
214  real(RP), private, allocatable :: qc_sfc (:,:,:) ! surface cloud water [kg/kg]
215 
216  real(RP), private, allocatable :: rndm (:,:,:) ! random number (0-1)
217  real(RP), private, allocatable, target :: bubble (:,:,:) ! bubble factor (0-1)
218  real(RP), private, allocatable, target :: rect (:,:,:) ! rectangle factor (0-1)
219  real(RP), private, allocatable :: gan (:) ! gamma factor (0-1)
220 
221  logical, private :: flg_intrp = .false.
222 
223  !-----------------------------------------------------------------------------
224 contains
225  !-----------------------------------------------------------------------------
227  subroutine mkinit_setup
228  implicit none
229 
230  character(len=H_SHORT) :: mkinit_initname = 'NONE'
231 
232  namelist / param_mkinit / &
233  mkinit_initname, &
234  flg_intrp
235 
236  integer :: ierr
237  !---------------------------------------------------------------------------
238 
239  if( io_l ) write(io_fid_log,*)
240  if( io_l ) write(io_fid_log,*) '++++++ Module[make init] / Categ[preprocess] / Origin[SCALE-RM]'
241 
242  !--- read namelist
243  rewind(io_fid_conf)
244  read(io_fid_conf,nml=param_mkinit,iostat=ierr)
245  if( ierr < 0 ) then !--- missing
246  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
247  elseif( ierr > 0 ) then !--- fatal error
248  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT. Check!'
249  call prc_mpistop
250  endif
251  if( io_nml ) write(io_fid_nml,nml=param_mkinit)
252 
253  allocate( pres(ka,ia,ja) )
254  allocate( temp(ka,ia,ja) )
255  allocate( pott(ka,ia,ja) )
256  allocate( qsat(ka,ia,ja) )
257  allocate( qv(ka,ia,ja) )
258  allocate( qc(ka,ia,ja) )
259  allocate( velx(ka,ia,ja) )
260  allocate( vely(ka,ia,ja) )
261 
262  allocate( pres_sfc(1,ia,ja) )
263  allocate( temp_sfc(1,ia,ja) )
264  allocate( pott_sfc(1,ia,ja) )
265  allocate( qsat_sfc(1,ia,ja) )
266  allocate( qv_sfc(1,ia,ja) )
267  allocate( qc_sfc(1,ia,ja) )
268 
269  allocate( rndm(ka,ia,ja) )
270  allocate( bubble(ka,ia,ja) )
271  allocate( rect(ka,ia,ja) )
272 
273  select case(trim(mkinit_initname))
274  case('NONE')
275  mkinit_type = i_ignore
276  case('PLANESTATE')
277  mkinit_type = i_planestate
278  case('TRACERBUBBLE')
279  mkinit_type = i_tracerbubble
280  case('COLDBUBBLE')
281  mkinit_type = i_coldbubble
282  call bubble_setup
283  case('LAMBWAVE')
284  mkinit_type = i_lambwave
285  call bubble_setup
286  case('GRAVITYWAVE')
287  mkinit_type = i_gravitywave
288  call bubble_setup
289  case('KHWAVE')
290  mkinit_type = i_khwave
291  case('TURBULENCE')
292  mkinit_type = i_turbulence
293  case('MOUNTAINWAVE')
294  mkinit_type = i_mountainwave
295  call bubble_setup
296  case('WARMBUBBLE')
297  mkinit_type = i_warmbubble
298  call bubble_setup
299  case('SUPERCELL')
300  mkinit_type = i_supercell
301  call bubble_setup
302  case('SQUALLLINE')
303  mkinit_type = i_squallline
304  case('WK1982')
305  mkinit_type = i_wk1982
306  call bubble_setup
307  case('DYCOMS2_RF01')
308  mkinit_type = i_dycoms2_rf01
309  case('DYCOMS2_RF02')
310  mkinit_type = i_dycoms2_rf02
311  case('RICO')
312  mkinit_type = i_rico
313  case('BOMEX')
314  mkinit_type = i_bomex
315  case('INTERPORATION')
316  mkinit_type = i_interporation
317  case('LANDCOUPLE')
318  mkinit_type = i_landcouple
319  case('OCEANCOUPLE')
320  mkinit_type = i_oceancouple
321  case('URBANCOUPLE')
322  mkinit_type = i_urbancouple
323  case('TRIPLECOUPLE')
324  mkinit_type = i_triplecouple
325  case('BUBBLECOUPLE')
326  mkinit_type = i_bubblecouple
327  call bubble_setup
328  case('SEABREEZE')
329  mkinit_type = i_seabreeze
330  case('HEATISLAND')
331  mkinit_type = i_heatisland
332  case('DYCOMS2_RF02_DNS')
333  mkinit_type = i_dycoms2_rf02_dns
334  case('REAL')
335  mkinit_type = i_real
336  case('GRAYZONE')
337  mkinit_type = i_grayzone
338  case('BOXAERO')
339  mkinit_type = i_boxaero
340  case('WARMBUBBLEAERO')
341  mkinit_type = i_warmbubbleaero
342  call bubble_setup
343  case('CAVITYFLOW')
344  mkinit_type = i_cavityflow
345  case('BAROCWAVE')
346  mkinit_type = i_barocwave
347  case default
348  write(*,*) 'xxx Unsupported TYPE:', trim(mkinit_initname)
349  call prc_mpistop
350  endselect
351 
352  return
353  end subroutine mkinit_setup
354 
355  !-----------------------------------------------------------------------------
357  subroutine mkinit
358  use scale_const, only: &
359  const_undef8
360  use scale_landuse, only: &
361  landuse_write
362  use mod_atmos_admin, only: &
363  atmos_phy_mp_type
364  use mod_atmos_driver, only: &
365  atmos_surface_get
366  use mod_ocean_driver, only: &
367  ocean_surface_set
368  use mod_land_driver, only: &
369  land_surface_set
370  use mod_urban_driver, only: &
371  urban_surface_set
372  use mod_admin_restart, only: &
373  admin_restart_write
374  use mod_admin_time, only: &
375  time_doatmos_restart, &
376  time_doland_restart, &
377  time_dourban_restart, &
378  time_doocean_restart
379  implicit none
380 
381  integer :: iq
382  !---------------------------------------------------------------------------
383 
384  if ( mkinit_type == i_ignore ) then
385  if( io_l ) write(io_fid_log,*)
386  if( io_l ) write(io_fid_log,*) '++++++ SKIP MAKING INITIAL DATA ++++++'
387  else
388  if( io_l ) write(io_fid_log,*)
389  if( io_l ) write(io_fid_log,*) '++++++ START MAKING INITIAL DATA ++++++'
390 
391  !--- Initialize variables
392 #ifndef DRY
393  do iq = 2, qa
394  qtrc(:,:,:,iq) = 0.0_rp
395  enddo
396 #endif
397 
398  pres(:,:,:) = const_undef8
399  temp(:,:,:) = const_undef8
400  pott(:,:,:) = const_undef8
401  qsat(:,:,:) = const_undef8
402  qv(:,:,:) = const_undef8
403  qc(:,:,:) = const_undef8
404  velx(:,:,:) = const_undef8
405  vely(:,:,:) = const_undef8
406 
407  rndm(:,:,:) = const_undef8
408 
409  pres_sfc(:,:,:) = const_undef8
410  temp_sfc(:,:,:) = const_undef8
411  pott_sfc(:,:,:) = const_undef8
412  qsat_sfc(:,:,:) = const_undef8
413  qv_sfc(:,:,:) = const_undef8
414  qc_sfc(:,:,:) = const_undef8
415 
416  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
417  if( io_l ) write(io_fid_log,*) '*** Aerosols for SBM are included ***'
418  call sbmaero_setup
419  endif
420 
421  select case(mkinit_type)
422  case(i_planestate)
423  call mkinit_planestate
424  case(i_tracerbubble)
425  call mkinit_tracerbubble
426  case(i_coldbubble)
427  call mkinit_coldbubble
428  case(i_lambwave)
429  call mkinit_lambwave
430  case(i_gravitywave)
431  call mkinit_gravitywave
432  case(i_khwave)
433  call mkinit_khwave
434  case(i_turbulence)
435  call mkinit_turbulence
436  case(i_mountainwave)
437  call mkinit_mountainwave
438  case(i_warmbubble)
439  call mkinit_warmbubble
440  case(i_supercell)
441  call mkinit_supercell
442  case(i_squallline)
443  call mkinit_squallline
444  case(i_wk1982)
445  call mkinit_wk1982
446  case(i_dycoms2_rf01)
447  call mkinit_dycoms2_rf01
448  case(i_dycoms2_rf02)
449  call mkinit_dycoms2_rf02
450  case(i_rico)
451  call mkinit_rico
452  case(i_bomex)
453  call mkinit_bomex
454  case(i_interporation)
455  call mkinit_interporation
456  case(i_oceancouple)
457  call mkinit_planestate
458  call mkinit_oceancouple
459  case(i_landcouple)
460  call mkinit_planestate
461  call mkinit_landcouple
462  case(i_urbancouple)
463  call mkinit_planestate
464  call mkinit_landcouple ! tentative
465  call mkinit_urbancouple
466  case(i_triplecouple)
467  call mkinit_planestate
468  call mkinit_oceancouple
469  call mkinit_landcouple
470  call mkinit_urbancouple
471  case(i_bubblecouple)
472  call mkinit_planestate
473  call mkinit_warmbubble
474  call mkinit_oceancouple
475  call mkinit_landcouple
476  call mkinit_urbancouple
477  case(i_seabreeze)
478  call mkinit_planestate
479  call mkinit_seabreeze
480  case(i_heatisland)
481  call mkinit_planestate
482  call mkinit_heatisland
483  case(i_dycoms2_rf02_dns)
484  call mkinit_dycoms2_rf02_dns
485  case(i_real)
486  call mkinit_real
487  case(i_grayzone)
488  call mkinit_grayzone
489  case(i_boxaero)
490  call mkinit_boxaero
491  case(i_warmbubbleaero)
492  call mkinit_warmbubbleaero
493  case(i_cavityflow)
494  call mkinit_cavityflow
495  case(i_barocwave)
496  call mkinit_barocwave
497  case default
498  write(*,*) 'xxx Unsupported TYPE:', mkinit_type
499  call prc_mpistop
500  endselect
501 
502  call tke_setup
503 
504  if( io_l ) write(io_fid_log,*) '++++++ END MAKING INITIAL DATA ++++++'
505 
506  ! setup surface condition
507  call ocean_surface_set( countup = .false. )
508  call land_surface_set ( countup = .false. )
509  call urban_surface_set( countup = .false. )
510  call atmos_surface_get
511 
512  ! output boundary file
513  call landuse_write
514 
515  ! output restart file
516  time_doocean_restart = .true.
517  time_doland_restart = .true.
518  time_dourban_restart = .true.
519  time_doatmos_restart = .true.
520  call admin_restart_write
521 
522  endif
523 
524  return
525  end subroutine mkinit
526 
527  !-----------------------------------------------------------------------------
529  subroutine bubble_setup
530  use scale_const, only: &
531  const_undef8
532  implicit none
533 
534  ! Bubble
535  logical :: bbl_eachnode = .false. ! Arrange bubble at each node? [kg/kg]
536  real(RP) :: bbl_cz = 2.e3_rp ! center location [m]: z
537  real(RP) :: bbl_cx = 2.e3_rp ! center location [m]: x
538  real(RP) :: bbl_cy = 2.e3_rp ! center location [m]: y
539  real(RP) :: bbl_rz = 0.0_rp ! bubble radius [m]: z
540  real(RP) :: bbl_rx = 0.0_rp ! bubble radius [m]: x
541  real(RP) :: bbl_ry = 0.0_rp ! bubble radius [m]: y
542 
543  namelist / param_bubble / &
544  bbl_eachnode, &
545  bbl_cz, &
546  bbl_cx, &
547  bbl_cy, &
548  bbl_rz, &
549  bbl_rx, &
550  bbl_ry
551 
552  real(RP) :: cz_offset
553  real(RP) :: cx_offset
554  real(RP) :: cy_offset
555  real(RP) :: distx, disty, distz
556 
557  real(RP) :: domain_rx, domain_ry
558 
559  integer :: ierr
560  integer :: k, i, j
561  !---------------------------------------------------------------------------
562 
563  if( io_l ) write(io_fid_log,*)
564  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit bubble] / Categ[preprocess] / Origin[SCALE-RM]'
565 
566  !--- read namelist
567  rewind(io_fid_conf)
568  read(io_fid_conf,nml=param_bubble,iostat=ierr)
569  if( ierr < 0 ) then !--- missing
570  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
571  elseif( ierr > 0 ) then !--- fatal error
572  write(*,*) 'xxx Not appropriate names in namelist PARAM_BUBBLE. Check!'
573  call prc_mpistop
574  endif
575  if( io_nml ) write(io_fid_nml,nml=param_bubble)
576 
577  if ( abs(bbl_rz*bbl_rx*bbl_ry) <= 0.0_rp ) then
578  if( io_l ) write(io_fid_log,*) '*** no bubble'
579  bubble(:,:,:) = 0.0_rp
580  else
581 
582  bubble(:,:,:) = const_undef8
583 
584  if ( bbl_eachnode ) then
585  cz_offset = grid_cz(ks)
586  cx_offset = grid_cx(is)
587  cy_offset = grid_cy(js)
588  domain_rx = grid_fx(ie) - grid_fx(is-1)
589  domain_ry = grid_fy(je) - grid_fy(js-1)
590  else
591  cz_offset = 0.0_rp
592  cx_offset = 0.0_rp
593  cy_offset = 0.0_rp
594  domain_rx = grid_fxg(iag-ihalo) - grid_fxg(ihalo)
595  domain_ry = grid_fyg(jag-jhalo) - grid_fyg(jhalo)
596  endif
597 
598  ! make bubble coefficient
599  do j = 1, ja
600  do i = 1, ia
601  do k = ks, ke
602 
603  distz = ( (grid_cz(k)-cz_offset-bbl_cz)/bbl_rz )**2
604 
605  distx = min( ( (grid_cx(i)-cx_offset-bbl_cx )/bbl_rx )**2, &
606  ( (grid_cx(i)-cx_offset-bbl_cx-domain_rx)/bbl_rx )**2, &
607  ( (grid_cx(i)-cx_offset-bbl_cx+domain_rx)/bbl_rx )**2 )
608 
609  disty = min( ( (grid_cy(j)-cy_offset-bbl_cy )/bbl_ry )**2, &
610  ( (grid_cy(j)-cy_offset-bbl_cy-domain_ry)/bbl_ry )**2, &
611  ( (grid_cy(j)-cy_offset-bbl_cy+domain_ry)/bbl_ry )**2 )
612 
613  bubble(k,i,j) = cos( 0.5_rp*pi*sqrt( min(distz+distx+disty,1.0_rp) ) )**2
614 
615  enddo
616  enddo
617  enddo
618  endif
619 
620  return
621  end subroutine bubble_setup
622 
623  !-----------------------------------------------------------------------------
625  subroutine rect_setup
626  use scale_const, only: &
627  const_undef8
628  implicit none
629 
630  ! Bubble
631  logical :: RCT_eachnode = .false. ! Arrange rectangle at each node? [kg/kg]
632  real(RP) :: RCT_CZ = 2.e3_rp ! center location [m]: z
633  real(RP) :: RCT_CX = 2.e3_rp ! center location [m]: x
634  real(RP) :: RCT_CY = 2.e3_rp ! center location [m]: y
635  real(RP) :: RCT_RZ = 2.e3_rp ! rectangle z width [m]: z
636  real(RP) :: RCT_RX = 2.e3_rp ! rectangle x width [m]: x
637  real(RP) :: RCT_RY = 2.e3_rp ! rectangle y width [m]: y
638 
639  namelist / param_rect / &
640  rct_eachnode, &
641  rct_cz, &
642  rct_cx, &
643  rct_cy, &
644  rct_rz, &
645  rct_rx, &
646  rct_ry
647 
648  real(RP) :: CZ_offset
649  real(RP) :: CX_offset
650  real(RP) :: CY_offset
651  real(RP) :: dist
652 
653  integer :: ierr
654  integer :: k, i, j
655  !---------------------------------------------------------------------------
656 
657  if( io_l ) write(io_fid_log,*)
658  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit rectangle] / Categ[preprocess] / Origin[SCALE-RM]'
659 
660  !--- read namelist
661  rewind(io_fid_conf)
662  read(io_fid_conf,nml=param_rect,iostat=ierr)
663  if( ierr < 0 ) then !--- missing
664  write(*,*) 'xxx Not found namelist. Check!'
665  call prc_mpistop
666  elseif( ierr > 0 ) then !--- fatal error
667  write(*,*) 'xxx Not appropriate names in namelist PARAM_RECT. Check!'
668  call prc_mpistop
669  endif
670  if( io_nml ) write(io_fid_nml,nml=param_rect)
671 
672  rect(:,:,:) = const_undef8
673 
674  if ( rct_eachnode ) then
675  cz_offset = grid_cz(ks)
676  cx_offset = grid_cx(is)
677  cy_offset = grid_cy(js)
678  else
679  cz_offset = 0.0_rp
680  cx_offset = 0.0_rp
681  cy_offset = 0.0_rp
682  endif
683 
684  do j = 1, ja
685  do i = 1, ia
686  do k = ks, ke
687 
688  ! make tracer rectangle
689  dist = 2.0_rp * max( &
690  abs(grid_cz(k) - cz_offset - rct_cz)/rct_rz, &
691  abs(grid_cx(i) - cx_offset - rct_cx)/rct_rx, &
692  abs(grid_cy(j) - cy_offset - rct_cy)/rct_ry &
693  & )
694  if ( dist <= 1.0_rp ) then
695  rect(k,i,j) = 1.0_rp
696  else
697  rect(k,i,j) = 0.0_rp
698  end if
699  enddo
700  enddo
701  enddo
702 
703  return
704  end subroutine rect_setup
705 
706  !-----------------------------------------------------------------------------
708  subroutine aerosol_setup
709  use scale_atmos_phy_ae_kajino13, only: &
710  atmos_phy_ae_kajino13_mkinit
711  implicit none
712 
713  real(RP), parameter :: d_min_def = 1.e-9_rp ! default lower bound of 1st size bin
714  real(RP), parameter :: d_max_def = 1.e-5_rp ! upper bound of last size bin
715  integer, parameter :: n_kap_def = 1 ! number of kappa bins
716  real(RP), parameter :: k_min_def = 0.e0_rp ! lower bound of 1st kappa bin
717  real(RP), parameter :: k_max_def = 1.e0_rp ! upper bound of last kappa bin
718 
719  real(RP) :: m0_init = 0.0_rp ! initial total num. conc. of modes (Atk,Acm,Cor) [#/m3]
720  real(RP) :: dg_init = 80.e-9_rp ! initial number equivalen diameters of modes [m]
721  real(RP) :: sg_init = 1.6_rp ! initial standard deviation [-]
722 
723  real(RP) :: d_min_inp(3) = d_min_def
724  real(RP) :: d_max_inp(3) = d_max_def
725  real(RP) :: k_min_inp(3) = k_min_def
726  real(RP) :: k_max_inp(3) = k_max_def
727  integer :: n_kap_inp(3) = n_kap_def
728 
729  namelist / param_aero / &
730  m0_init, &
731  dg_init, &
732  sg_init, &
733  d_min_inp, &
734  d_max_inp, &
735  k_min_inp, &
736  k_max_inp, &
737  n_kap_inp
738 
739  integer :: ierr
740  !---------------------------------------------------------------------------
741 
742  if( io_l ) write(io_fid_log,*)
743  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit aerosol] / Categ[preprocess] / Origin[SCALE-RM]'
744 
745  !--- read namelist
746  rewind(io_fid_conf)
747  read(io_fid_conf,nml=param_aero,iostat=ierr)
748  if( ierr < 0 ) then !--- missing
749  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used!'
750  elseif( ierr > 0 ) then !--- fatal error
751  write(*,*) 'xxx Not appropriate names in namelist PARAM_AERO. Check!'
752  call prc_mpistop
753  endif
754  if( io_nml ) write(io_fid_nml,nml=param_aero)
755 
756  call atmos_phy_ae_kajino13_mkinit( qtrc, ccn, & ! (out)
757  dens, rhot, & ! (in)
758  m0_init, dg_init, sg_init, & ! (in)
759  d_min_inp, d_max_inp, & ! (in)
760  k_min_inp, k_max_inp, & ! (in)
761  n_kap_inp ) ! (in)
762 
763  return
764  end subroutine aerosol_setup
765 
766  !-----------------------------------------------------------------------------
768  subroutine sbmaero_setup
769  use scale_const, only: &
770  pi => const_pi
771  use scale_atmos_hydrometeor, only: &
772  qhs, &
773  qhe
774  use scale_atmos_phy_mp_suzuki10, only: &
775  nccn, nbin
776  implicit none
777 
778 #ifndef DRY
779  real(RP) :: xasta, xaend, dxaer
780  real(RP), allocatable :: xabnd( : ), xactr( : )
781 
782  real(RP) :: F0_AERO = 1.e+7_rp
783  real(RP) :: R0_AERO = 1.e-7_rp
784  real(RP) :: R_MAX = 1.e-06_rp
785  real(RP) :: R_MIN = 1.e-08_rp
786  real(RP) :: A_ALPHA = 3.0_rp
787  real(RP) :: RHO_AERO = 2.25e+03_rp
788 
789  namelist / param_sbmaero / &
790  f0_aero, &
791  r0_aero, &
792  r_max, &
793  r_min, &
794  a_alpha, &
795  rho_aero
796 
797  integer :: ierr
798  integer :: iq, i, j, k
799  !---------------------------------------------------------------------------
800 
801  if( io_l ) write(io_fid_log,*)
802  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit aerobin] / Categ[preprocess] / Origin[SCALE-RM]'
803 
804  !--- read namelist
805  rewind(io_fid_conf)
806  read(io_fid_conf,nml=param_sbmaero,iostat=ierr)
807 
808  if( ierr < 0 ) then !--- missing
809  if( io_l ) write(io_fid_log,*) '*** Not found namelist. default value used'
810  elseif( ierr > 0 ) then !--- fatal error
811  write(*,*) 'xxx Not appropriate names in namelist SBMAERO. Check!'
812  call prc_mpistop
813  endif
814  if( io_nml ) write(io_fid_nml,nml=param_sbmaero)
815 
816  if( nccn /= 0 ) then
817  allocate( gan( nccn ) )
818  allocate( xactr(nccn) )
819  allocate( xabnd(nccn+1) )
820 
821  xasta = log( rho_aero*4.0_rp/3.0_rp*pi * ( r_min )**3 )
822  xaend = log( rho_aero*4.0_rp/3.0_rp*pi * ( r_max )**3 )
823  dxaer = ( xaend-xasta )/nccn
824  do iq = 1, nccn+1
825  xabnd( iq ) = xasta + dxaer*( iq-1 )
826  enddo
827  do iq = 1, nccn
828  xactr( iq ) = ( xabnd( iq )+xabnd( iq+1 ) )*0.5_rp
829  enddo
830  do iq = 1, nccn
831  gan( iq ) = faero( f0_aero,r0_aero,xactr( iq ), a_alpha, rho_aero )*exp( xactr(iq) )
832  enddo
833  endif
834 
835  !--- Hydrometeor is zero at initial time for Bin method
836  do iq = qhs, qhe
837  do j = js, je
838  do i = is, ie
839  do k = ks, ke
840  qtrc(k,i,j,iq) = 0.0_rp
841  enddo
842  enddo
843  enddo
844  enddo
845 
846  !-- Aerosol distribution
847  if( nccn /= 0 ) then
848  do iq = 1, nccn
849  do j = js, je
850  do i = is, ie
851  do k = ks, ke
852  qtrc(k,i,j,qhe+iq) = gan(iq) !/ DENS(k,i,j)
853  enddo
854  enddo
855  enddo
856  enddo
857 
858  deallocate( xactr )
859  deallocate( xabnd )
860  endif
861 
862 #endif
863  return
864  end subroutine sbmaero_setup
865 
866  !-----------------------------------------------------------------------------
867  function faero( f0,r0,x,alpha,rhoa )
868  use scale_const, only: &
869  pi => const_pi
870  implicit none
871 
872  real(RP), intent(in) :: x, f0, r0, alpha, rhoa
873  real(RP) :: faero
874  real(RP) :: rad
875  !---------------------------------------------------------------------------
876 
877  rad = ( exp(x) * 3.0_rp / 4.0_rp / pi / rhoa )**(1.0_rp/3.0_rp)
878 
879  faero = f0 * (rad/r0)**(-alpha)
880 
881  return
882  end function faero
883 
884  !-----------------------------------------------------------------------------
886  subroutine flux_setup
887  use mod_atmos_phy_mp_vars, only: &
888  sflx_rain => atmos_phy_mp_sflx_rain, &
889  sflx_snow => atmos_phy_mp_sflx_snow
890  use mod_atmos_phy_rd_vars, only: &
891  sflx_lw_up => atmos_phy_rd_sflx_lw_up, &
892  sflx_lw_dn => atmos_phy_rd_sflx_lw_dn, &
893  sflx_sw_up => atmos_phy_rd_sflx_sw_up, &
894  sflx_sw_dn => atmos_phy_rd_sflx_sw_dn, &
895  toaflx_lw_up => atmos_phy_rd_toaflx_lw_up, &
896  toaflx_lw_dn => atmos_phy_rd_toaflx_lw_dn, &
897  toaflx_sw_up => atmos_phy_rd_toaflx_sw_up, &
898  toaflx_sw_dn => atmos_phy_rd_toaflx_sw_dn, &
899  sflx_rad_dn => atmos_phy_rd_sflx_downall
900  implicit none
901 
902  ! Flux from Atmosphere
903  real(RP) :: FLX_rain = 0.0_rp ! surface rain flux [kg/m2/s]
904  real(RP) :: FLX_snow = 0.0_rp ! surface snow flux [kg/m2/s]
905  real(RP) :: FLX_LW_dn = 0.0_rp ! surface downwad long-wave radiation flux [J/m2/s]
906  real(RP) :: FLX_SW_dn = 0.0_rp ! surface downwad short-wave radiation flux [J/m2/s]
907 
908  namelist / param_mkinit_flux / &
909  flx_rain, &
910  flx_snow, &
911  flx_lw_dn, &
912  flx_sw_dn
913 
914  integer :: i, j
915  integer :: ierr
916  !---------------------------------------------------------------------------
917 
918  !--- read namelist
919  rewind(io_fid_conf)
920  read(io_fid_conf,nml=param_mkinit_flux,iostat=ierr)
921  if( ierr < 0 ) then !--- missing
922  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
923  elseif( ierr > 0 ) then !--- fatal error
924  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_FLUX. Check!'
925  call prc_mpistop
926  endif
927  if( io_nml ) write(io_fid_nml,nml=param_mkinit_flux)
928 
929  do j = js, je
930  do i = is, ie
931  sflx_rain(i,j) = flx_rain
932  sflx_snow(i,j) = flx_snow
933 
934  sflx_lw_up(i,j) = 0.0_rp
935  sflx_lw_dn(i,j) = flx_lw_dn
936  sflx_sw_up(i,j) = 0.0_rp
937  sflx_sw_dn(i,j) = flx_sw_dn
938 
939  toaflx_lw_up(i,j) = 0.0_rp
940  toaflx_lw_dn(i,j) = 0.0_rp
941  toaflx_sw_up(i,j) = 0.0_rp
942  toaflx_sw_dn(i,j) = 0.0_rp
943 
944  sflx_rad_dn(i,j,1,1) = flx_sw_dn
945  sflx_rad_dn(i,j,1,2) = 0.0_rp
946  sflx_rad_dn(i,j,2,1) = flx_lw_dn
947  sflx_rad_dn(i,j,2,2) = 0.0_rp
948  enddo
949  enddo
950 
951  return
952  end subroutine flux_setup
953 
954  !-----------------------------------------------------------------------------
956  subroutine land_setup
957  use mod_land_vars, only: &
958  land_temp, &
959  land_water, &
960  land_sfc_temp, &
961  land_sfc_albedo
962  implicit none
963 
964  ! Land state
965  real(RP) :: LND_TEMP ! soil temperature [K]
966  real(RP) :: LND_WATER = 0.15_rp ! soil moisture [m3/m3]
967  real(RP) :: SFC_TEMP ! land skin temperature [K]
968  real(RP) :: SFC_albedo_LW = 0.01_rp ! land surface albedo for LW (0-1)
969  real(RP) :: SFC_albedo_SW = 0.20_rp ! land surface albedo for SW (0-1)
970 
971  integer :: i, j
972  integer :: ierr
973 
974  namelist /param_mkinit_land/ &
975  lnd_temp, &
976  lnd_water, &
977  sfc_temp, &
978  sfc_albedo_lw, &
979  sfc_albedo_sw
980 
981  lnd_temp = thetastd
982  sfc_temp = thetastd
983 
984  !--- read namelist
985  rewind(io_fid_conf)
986  read(io_fid_conf,nml=param_mkinit_land,iostat=ierr)
987  if( ierr < 0 ) then !--- missing
988  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
989  elseif( ierr > 0 ) then !--- fatal error
990  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_LAND. Check!'
991  call prc_mpistop
992  endif
993  if( io_nml ) write(io_fid_nml,nml=param_mkinit_land)
994 
995  do j = js, je
996  do i = is, ie
997  land_temp(:,i,j) = lnd_temp
998  land_water(:,i,j) = lnd_water
999  land_sfc_temp(i,j) = sfc_temp
1000  land_sfc_albedo(i,j,i_lw) = sfc_albedo_lw
1001  land_sfc_albedo(i,j,i_sw) = sfc_albedo_sw
1002  end do
1003  end do
1004 
1005  return
1006  end subroutine land_setup
1007 
1008  !-----------------------------------------------------------------------------
1010  subroutine ocean_setup
1011  use mod_ocean_vars, only: &
1012  ocean_temp, &
1013  ocean_sfc_temp, &
1014  ocean_sfc_albedo, &
1015  ocean_sfc_z0m, &
1016  ocean_sfc_z0h, &
1017  ocean_sfc_z0e
1018  implicit none
1019 
1020  ! Ocean state
1021  real(RP) :: OCN_TEMP ! ocean temperature [K]
1022  real(RP) :: SFC_TEMP ! ocean skin temperature [K]
1023  real(RP) :: SFC_albedo_LW = 0.04_rp ! ocean surface albedo for LW (0-1)
1024  real(RP) :: SFC_albedo_SW = 0.05_rp ! ocean surface albedo for SW (0-1)
1025  real(RP) :: SFC_Z0M = 1.0e-4_rp ! ocean surface roughness length (momentum) [m]
1026  real(RP) :: SFC_Z0H = 1.0e-4_rp ! ocean surface roughness length (heat) [m]
1027  real(RP) :: SFC_Z0E = 1.0e-4_rp ! ocean surface roughness length (vapor) [m]
1028 
1029  integer :: i, j
1030  integer :: ierr
1031 
1032  namelist /param_mkinit_ocean/ &
1033  ocn_temp, &
1034  sfc_temp, &
1035  sfc_albedo_lw, &
1036  sfc_albedo_sw, &
1037  sfc_z0m, &
1038  sfc_z0h, &
1039  sfc_z0e
1040 
1041  ocn_temp = thetastd
1042  sfc_temp = thetastd
1043 
1044  !--- read namelist
1045  rewind(io_fid_conf)
1046  read(io_fid_conf,nml=param_mkinit_ocean,iostat=ierr)
1047  if( ierr < 0 ) then !--- missing
1048  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1049  elseif( ierr > 0 ) then !--- fatal error
1050  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_OCEAN. Check!'
1051  call prc_mpistop
1052  endif
1053  if( io_nml ) write(io_fid_nml,nml=param_mkinit_ocean)
1054 
1055 
1056  do j = js, je
1057  do i = is, ie
1058  ocean_temp(i,j) = ocn_temp
1059  ocean_sfc_temp(i,j) = sfc_temp
1060  ocean_sfc_albedo(i,j,i_lw) = sfc_albedo_lw
1061  ocean_sfc_albedo(i,j,i_sw) = sfc_albedo_sw
1062  ocean_sfc_z0m(i,j) = sfc_z0m
1063  ocean_sfc_z0h(i,j) = sfc_z0h
1064  ocean_sfc_z0e(i,j) = sfc_z0e
1065  end do
1066  end do
1067 
1068  return
1069  end subroutine ocean_setup
1070 
1071  !-----------------------------------------------------------------------------
1073  subroutine urban_setup
1074  use mod_urban_vars, only: &
1075  urban_tr, &
1076  urban_tb, &
1077  urban_tg, &
1078  urban_tc, &
1079  urban_qc, &
1080  urban_uc, &
1081  urban_trl, &
1082  urban_tbl, &
1083  urban_tgl, &
1084  urban_rainr, &
1085  urban_rainb, &
1086  urban_raing, &
1087  urban_roff, &
1088  urban_sfc_temp, &
1089  urban_sfc_albedo
1090  implicit none
1091 
1092  ! urban state
1093  real(RP) :: URB_ROOF_TEMP ! Surface temperature of roof [K]
1094  real(RP) :: URB_BLDG_TEMP ! Surface temperature of building [K
1095  real(RP) :: URB_GRND_TEMP ! Surface temperature of ground [K]
1096  real(RP) :: URB_CNPY_TEMP ! Diagnostic canopy air temperature
1097  real(RP) :: URB_CNPY_HMDT = 0.0_rp ! Diagnostic canopy humidity [-]
1098  real(RP) :: URB_CNPY_WIND = 0.0_rp ! Diagnostic canopy wind [m/s]
1099  real(RP) :: URB_ROOF_LAYER_TEMP ! temperature in layer of roof [K]
1100  real(RP) :: URB_BLDG_LAYER_TEMP ! temperature in layer of building [
1101  real(RP) :: URB_GRND_LAYER_TEMP ! temperature in layer of ground [K]
1102  real(RP) :: URB_ROOF_RAIN = 0.0_rp ! temperature in layer of roof [K]
1103  real(RP) :: URB_BLDG_RAIN = 0.0_rp ! temperature in layer of building [
1104  real(RP) :: URB_GRND_RAIN = 0.0_rp ! temperature in layer of ground [K]
1105  real(RP) :: URB_RUNOFF = 0.0_rp ! temperature in layer of ground [K]
1106  real(RP) :: URB_SFC_TEMP ! Grid average of surface temperature [K]
1107  real(RP) :: URB_ALB_LW = 0.0_rp ! Grid average of surface albedo for LW (0-1)
1108  real(RP) :: URB_ALB_SW = 0.0_rp ! Grid average of surface albedo for SW (0-1)
1109 
1110  integer :: i, j
1111  integer :: ierr
1112 
1113  namelist /param_mkinit_urban/ &
1114  urb_roof_temp, &
1115  urb_bldg_temp, &
1116  urb_grnd_temp, &
1117  urb_cnpy_temp, &
1118  urb_cnpy_hmdt, &
1119  urb_cnpy_wind, &
1120  urb_roof_layer_temp, &
1121  urb_bldg_layer_temp, &
1122  urb_grnd_layer_temp, &
1123  urb_roof_rain, &
1124  urb_bldg_rain, &
1125  urb_grnd_rain, &
1126  urb_runoff, &
1127  urb_sfc_temp, &
1128  urb_alb_lw, &
1129  urb_alb_sw
1130 
1131  urb_roof_temp = thetastd
1132  urb_bldg_temp = thetastd
1133  urb_grnd_temp = thetastd
1134  urb_cnpy_temp = thetastd
1135  urb_roof_layer_temp = thetastd
1136  urb_bldg_layer_temp = thetastd
1137  urb_grnd_layer_temp = thetastd
1138 
1139  urb_sfc_temp = thetastd
1140 
1141  !--- read namelist
1142  rewind(io_fid_conf)
1143  read(io_fid_conf,nml=param_mkinit_urban,iostat=ierr)
1144  if( ierr < 0 ) then !--- missing
1145  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1146  elseif( ierr > 0 ) then !--- fatal error
1147  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_URBAN. Check!'
1148  call prc_mpistop
1149  endif
1150  if( io_nml ) write(io_fid_nml,nml=param_mkinit_urban)
1151 
1152 
1153  do j = js, je
1154  do i = is, ie
1155  urban_tr(i,j) = urb_roof_temp
1156  urban_tb(i,j) = urb_bldg_temp
1157  urban_tg(i,j) = urb_grnd_temp
1158  urban_tc(i,j) = urb_cnpy_temp
1159  urban_qc(i,j) = urb_cnpy_hmdt
1160  urban_uc(i,j) = urb_cnpy_wind
1161  urban_trl(:,i,j) = urb_roof_layer_temp
1162  urban_tbl(:,i,j) = urb_bldg_layer_temp
1163  urban_tgl(:,i,j) = urb_grnd_layer_temp
1164  urban_rainr(i,j) = urb_roof_rain
1165  urban_rainb(i,j) = urb_bldg_rain
1166  urban_raing(i,j) = urb_grnd_rain
1167  urban_roff(i,j) = urb_runoff
1168  urban_sfc_temp(i,j) = urb_sfc_temp
1169  urban_sfc_albedo(i,j,i_lw) = urb_alb_lw
1170  urban_sfc_albedo(i,j,i_sw) = urb_alb_sw
1171  end do
1172  end do
1173 
1174  return
1175  end subroutine urban_setup
1176 
1177  !-----------------------------------------------------------------------------
1179  subroutine tke_setup
1180  use scale_const, only: &
1181  eps => const_eps
1182  use scale_atmos_phy_tb, only: &
1183  i_tke
1184  implicit none
1185 
1186  real(RP) :: TKE_CONST
1187 
1188  namelist / param_mkinit_tke / &
1189  tke_const
1190 
1191  integer :: k, i, j
1192  integer :: ierr
1193  !---------------------------------------------------------------------------
1194 
1195  tke_const = eps
1196 
1197  !--- read namelist
1198  rewind(io_fid_conf)
1199  read(io_fid_conf,nml=param_mkinit_tke,iostat=ierr)
1200  if( ierr < 0 ) then !--- missing
1201  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1202  elseif( ierr > 0 ) then !--- fatal error
1203  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_TKE. Check!'
1204  call prc_mpistop
1205  endif
1206  if( io_nml ) write(io_fid_nml,nml=param_mkinit_tke)
1207 
1208  if ( i_tke > 0 ) then
1209  do j = 1, ja
1210  do i = 1, ia
1211  do k = 1, ka
1212  qtrc(k,i,j,i_tke) = tke_const
1213  enddo
1214  enddo
1215  enddo
1216  end if
1217 
1218  return
1219  end subroutine tke_setup
1220 
1221  !-----------------------------------------------------------------------------
1223  subroutine read_sounding( &
1224  DENS, VELX, VELY, POTT, QV )
1225  use scale_atmos_hydrometeor, only: &
1226  i_qv
1227  implicit none
1228 
1229  real(RP), intent(out) :: DENS(KA)
1230  real(RP), intent(out) :: VELX(KA)
1231  real(RP), intent(out) :: VELY(KA)
1232  real(RP), intent(out) :: POTT(KA)
1233  real(RP), intent(out) :: QV (KA)
1234 
1235  real(RP) :: TEMP(KA)
1236  real(RP) :: PRES(KA)
1237  real(RP) :: QC (KA)
1238 
1239  character(len=H_LONG) :: ENV_IN_SOUNDING_file = ''
1240 
1241  integer, parameter :: EXP_klim = 100
1242  integer :: EXP_kmax
1243 
1244  real(RP) :: SFC_THETA ! surface potential temperature [K]
1245  real(RP) :: SFC_PRES ! surface pressure [hPa]
1246  real(RP) :: SFC_QV ! surface watervapor [g/kg]
1247 
1248  real(RP) :: EXP_z (EXP_klim+1) ! height [m]
1249  real(RP) :: EXP_pott(EXP_klim+1) ! potential temperature [K]
1250  real(RP) :: EXP_qv (EXP_klim+1) ! water vapor [g/kg]
1251  real(RP) :: EXP_u (EXP_klim+1) ! velocity u [m/s]
1252  real(RP) :: EXP_v (EXP_klim+1) ! velocity v [m/s]
1253 
1254  real(RP) :: fact1, fact2
1255  integer :: k, kref
1256  integer :: fid
1257  integer :: ierr
1258 
1259  namelist /param_mkinit_sounding/ &
1260  env_in_sounding_file
1261 
1262  !--- read namelist
1263  rewind(io_fid_conf)
1264  read(io_fid_conf,nml=param_mkinit_sounding,iostat=ierr)
1265 
1266  if( ierr < 0 ) then !--- missing
1267  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1268  elseif( ierr > 0 ) then !--- fatal error
1269  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_SOUNDING. Check!'
1270  call prc_mpistop
1271  endif
1272  if( io_nml ) write(io_fid_nml,nml=param_mkinit_sounding)
1273 
1274  !--- prepare sounding profile
1275  if( io_l ) write(io_fid_log,*) '+++ Input sounding file:', trim(env_in_sounding_file)
1276  fid = io_get_available_fid()
1277  open( fid, &
1278  file = trim(env_in_sounding_file), &
1279  form = 'formatted', &
1280  status = 'old', &
1281  iostat = ierr )
1282 
1283  if ( ierr /= 0 ) then
1284  write(*,*) 'xxx [mod_mkinit/read_sounding] Input file not found!'
1285  endif
1286 
1287  !--- read sounding file till end
1288  read(fid,*) sfc_pres, sfc_theta, sfc_qv
1289 
1290  if( io_l ) write(io_fid_log,*) '+++ Surface pressure [hPa]', sfc_pres
1291  if( io_l ) write(io_fid_log,*) '+++ Surface pot. temp [K]', sfc_theta
1292  if( io_l ) write(io_fid_log,*) '+++ Surface water vapor [g/kg]', sfc_qv
1293 
1294  do k = 2, exp_klim
1295  read(fid,*,iostat=ierr) exp_z(k), exp_pott(k), exp_qv(k), exp_u(k), exp_v(k)
1296  if ( ierr /= 0 ) exit
1297  enddo
1298 
1299  exp_kmax = k - 1
1300  close(fid)
1301 
1302  ! Boundary
1303  exp_z(1) = 0.0_rp
1304  exp_pott(1) = sfc_theta
1305  exp_qv(1) = sfc_qv
1306  exp_u(1) = exp_u(2)
1307  exp_v(1) = exp_v(2)
1308  exp_z(exp_kmax+1) = 100.e3_rp
1309  exp_pott(exp_kmax+1) = exp_pott(exp_kmax)
1310  exp_qv(exp_kmax+1) = exp_qv(exp_kmax)
1311  exp_u(exp_kmax+1) = exp_u(exp_kmax)
1312  exp_v(exp_kmax+1) = exp_v(exp_kmax)
1313 
1314  do k = 1, exp_kmax+1
1315  exp_qv(k) = exp_qv(k) * 1.e-3_rp ! [g/kg]->[kg/kg]
1316  enddo
1317 
1318  ! calc in dry condition
1319  pres_sfc = sfc_pres * 1.e2_rp ! [hPa]->[Pa]
1320  pott_sfc = sfc_theta
1321  if ( i_qv > 0 ) then
1322  qv_sfc = sfc_qv * 1.e-3_rp ! [g/kg]->[kg/kg]
1323  else
1324  qv_sfc = 0.0_rp
1325  end if
1326  qc_sfc = 0.0_rp
1327 
1328  !--- linear interpolate to model grid
1329  do k = ks, ke
1330  qc(k) = 0.0_rp
1331 
1332  do kref = 2, exp_kmax+1
1333  if ( grid_cz(k) > exp_z(kref-1) &
1334  .AND. grid_cz(k) <= exp_z(kref ) ) then
1335 
1336  fact1 = ( exp_z(kref) - grid_cz(k) ) / ( exp_z(kref)-exp_z(kref-1) )
1337  fact2 = ( grid_cz(k) - exp_z(kref-1) ) / ( exp_z(kref)-exp_z(kref-1) )
1338 
1339  pott(k) = exp_pott(kref-1) * fact1 &
1340  + exp_pott(kref ) * fact2
1341  qv(k) = exp_qv(kref-1) * fact1 &
1342  + exp_qv(kref ) * fact2
1343  velx(k) = exp_u(kref-1) * fact1 &
1344  + exp_u(kref ) * fact2
1345  vely(k) = exp_v(kref-1) * fact1 &
1346  + exp_v(kref ) * fact2
1347  endif
1348  enddo
1349  enddo
1350  if ( i_qv < 1 ) qv = 0.0_rp
1351 
1352  ! make density & pressure profile in moist condition
1353  call hydrostatic_buildrho( dens(:), & ! [OUT]
1354  temp(:), & ! [OUT]
1355  pres(:), & ! [OUT]
1356  pott(:), & ! [IN]
1357  qv(:), & ! [IN]
1358  qc(:), & ! [IN]
1359  temp_sfc(1,1,1), & ! [OUT]
1360  pres_sfc(1,1,1), & ! [IN]
1361  pott_sfc(1,1,1), & ! [IN]
1362  qv_sfc(1,1,1), & ! [IN]
1363  qc_sfc(1,1,1) ) ! [IN]
1364 
1365  return
1366  end subroutine read_sounding
1367 
1368  !-----------------------------------------------------------------------------
1370  subroutine mkinit_planestate
1371  use scale_atmos_hydrometeor, only: &
1372  i_qv
1373  implicit none
1374 
1375  ! Surface state
1376  real(RP) :: SFC_THETA ! surface potential temperature [K]
1377  real(RP) :: SFC_PRES ! surface pressure [Pa]
1378  real(RP) :: SFC_RH = 0.0_rp ! surface relative humidity [%]
1379  ! Environment state
1380  real(RP) :: ENV_THETA ! potential temperature of environment [K]
1381  real(RP) :: ENV_TLAPS = 0.0_rp ! Lapse rate of THETA [K/m]
1382  real(RP) :: ENV_U = 0.0_rp ! velocity u of environment [m/s]
1383  real(RP) :: ENV_V = 0.0_rp ! velocity v of environment [m/s]
1384  real(RP) :: ENV_RH = 0.0_rp ! relative humidity of environment [%]
1385  ! Disturbance
1386  real(RP) :: RANDOM_THETA = 0.0_rp ! amplitude of random disturbance theta
1387  real(RP) :: RANDOM_U = 0.0_rp ! amplitude of random disturbance u
1388  real(RP) :: RANDOM_V = 0.0_rp ! amplitude of random disturbance v
1389  real(RP) :: RANDOM_RH = 0.0_rp ! amplitude of random disturbance RH
1390 
1391  namelist / param_mkinit_planestate / &
1392  sfc_theta, &
1393  sfc_pres, &
1394  sfc_rh, &
1395  env_theta, &
1396  env_tlaps, &
1397  env_u, &
1398  env_v, &
1399  env_rh, &
1400  random_theta, &
1401  random_u, &
1402  random_v, &
1403  random_rh
1404 
1405  integer :: ierr
1406  integer :: k, i, j
1407  !---------------------------------------------------------------------------
1408 
1409  if( io_l ) write(io_fid_log,*)
1410  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit Horiz_UNIFORM] / Categ[preprocess] / Origin[SCALE-RM]'
1411 
1412  sfc_theta = thetastd
1413  sfc_pres = pstd
1414  env_theta = thetastd
1415 
1416  !--- read namelist
1417  rewind(io_fid_conf)
1418  read(io_fid_conf,nml=param_mkinit_planestate,iostat=ierr)
1419 
1420  if( ierr < 0 ) then !--- missing
1421  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1422  elseif( ierr > 0 ) then !--- fatal error
1423  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_PLANESTATE. Check!'
1424  call prc_mpistop
1425  endif
1426  if( io_nml ) write(io_fid_nml,nml=param_mkinit_planestate)
1427 
1428  ! calc in dry condition
1429  do j = js, je
1430  do i = is, ie
1431  pott_sfc(1,i,j) = sfc_theta
1432  pres_sfc(1,i,j) = sfc_pres
1433  qv_sfc(1,i,j) = 0.0_rp
1434  qc_sfc(1,i,j) = 0.0_rp
1435 
1436  do k = ks, ke
1437  qv(k,i,j) = 0.0_rp
1438  qc(k,i,j) = 0.0_rp
1439  enddo
1440  enddo
1441  enddo
1442 
1443  if ( env_theta < 0.0_rp ) then ! use isa profile
1444 
1445  call profile_isa( ka, ks, ke, & ! [IN]
1446  ia, is, ie, & ! [IN]
1447  ja, js, je, & ! [IN]
1448  pott_sfc(1,:,:), & ! [IN]
1449  pres_sfc(1,:,:), & ! [IN]
1450  real_cz(:,:,:), & ! [IN]
1451  pott(:,:,:) ) ! [OUT]
1452 
1453  else
1454 
1455  do j = js, je
1456  do i = is, ie
1457  do k = ks, ke
1458  pott(k,i,j) = env_theta + env_tlaps * real_cz(k,i,j)
1459  enddo
1460  enddo
1461  enddo
1462 
1463  endif
1464 
1465  ! make density & pressure profile in moist condition
1466  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
1467  temp(:,:,:), & ! [OUT]
1468  pres(:,:,:), & ! [OUT]
1469  pott(:,:,:), & ! [IN]
1470  qv(:,:,:), & ! [IN]
1471  qc(:,:,:), & ! [IN]
1472  temp_sfc(:,:,:), & ! [OUT]
1473  pres_sfc(:,:,:), & ! [IN]
1474  pott_sfc(:,:,:), & ! [IN]
1475  qv_sfc(:,:,:), & ! [IN]
1476  qc_sfc(:,:,:) ) ! [IN]
1477 
1478  if ( i_qv > 0 ) then
1479  ! calc QV from RH
1480  call saturation_pres2qsat_all( qsat_sfc(1,:,:), temp_sfc(1,:,:), pres_sfc(1,:,:) )
1481  call saturation_pres2qsat_all( qsat(:,:,:), temp(:,:,:), pres(:,:,:) )
1482 
1483  call random_get(rndm) ! make random
1484  do j = js, je
1485  do i = is, ie
1486  qv_sfc(1,i,j) = ( sfc_rh + rndm(ks-1,i,j) * random_rh ) * 1.e-2_rp * qsat_sfc(1,i,j)
1487 
1488  do k = ks, ke
1489  qv(k,i,j) = ( env_rh + rndm(k,i,j) * random_rh ) * 1.e-2_rp * qsat(k,i,j)
1490  enddo
1491  enddo
1492  enddo
1493  end if
1494 
1495  call random_get(rndm) ! make random
1496  do j = js, je
1497  do i = is, ie
1498  pott_sfc(1,i,j) = pott_sfc(1,i,j) + rndm(ks-1,i,j) * random_theta
1499 
1500  do k = ks, ke
1501  pott(k,i,j) = pott(k,i,j) + rndm(k,i,j) * random_theta
1502  enddo
1503  enddo
1504  enddo
1505 
1506  ! make density & pressure profile in moist condition
1507  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
1508  temp(:,:,:), & ! [OUT]
1509  pres(:,:,:), & ! [OUT]
1510  pott(:,:,:), & ! [IN]
1511  qv(:,:,:), & ! [IN]
1512  qc(:,:,:), & ! [IN]
1513  temp_sfc(:,:,:), & ! [OUT]
1514  pres_sfc(:,:,:), & ! [IN]
1515  pott_sfc(:,:,:), & ! [IN]
1516  qv_sfc(:,:,:), & ! [IN]
1517  qc_sfc(:,:,:) ) ! [IN]
1518 
1519  call comm_vars8( dens(:,:,:), 1 )
1520  call comm_wait ( dens(:,:,:), 1 )
1521 
1522  call random_get(rndm) ! make random
1523  do j = js, je
1524  do i = is, ie
1525  do k = ks, ke
1526  momx(k,i,j) = ( env_u + ( rndm(k,i,j) - 0.5_rp ) * 2.0_rp * random_u ) &
1527  * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
1528  enddo
1529  enddo
1530  enddo
1531 
1532  call random_get(rndm) ! make random
1533  do j = js, je
1534  do i = is, ie
1535  do k = ks, ke
1536  momy(k,i,j) = ( env_v + ( rndm(k,i,j) - 0.5_rp ) * 2.0_rp * random_v ) &
1537  * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
1538  enddo
1539  enddo
1540  enddo
1541 
1542  do j = js, je
1543  do i = is, ie
1544  do k = ks, ke
1545  momz(k,i,j) = 0.0_rp
1546  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
1547  enddo
1548  enddo
1549  enddo
1550 
1551  if ( i_qv > 0 ) then
1552  do j = js, je
1553  do i = is, ie
1554  do k = ks, ke
1555  qtrc(k,i,j,i_qv) = qv(k,i,j)
1556  enddo
1557  enddo
1558  enddo
1559  end if
1560 
1561  call flux_setup
1562 
1563  call ocean_setup
1564 
1565  return
1566  end subroutine mkinit_planestate
1567 
1568  !-----------------------------------------------------------------------------
1570  subroutine mkinit_tracerbubble
1571  use scale_atmos_hydrometeor, only: &
1572  i_nc
1573  use mod_atmos_admin, only: &
1574  atmos_phy_mp_type
1575  implicit none
1576 
1577 #ifndef DRY
1578  ! Surface state
1579  real(RP) :: SFC_THETA ! surface potential temperature [K]
1580  real(RP) :: SFC_PRES ! surface pressure [Pa]
1581  ! Environment state
1582  real(RP) :: ENV_THETA ! potential temperature of environment [K]
1583  real(RP) :: ENV_U = 0.0_rp ! velocity u of environment [m/s]
1584  real(RP) :: ENV_V = 0.0_rp ! velocity v of environment [m/s]
1585  ! Bubble
1586  character(len=H_SHORT) :: SHAPE_NC = 'BUBBLE' ! BUBBLE or RECT
1587  real(RP) :: BBL_NC = 1.0_rp ! extremum of NC in bubble [kg/kg]
1588 
1589  namelist / param_mkinit_tracerbubble / &
1590  sfc_theta, &
1591  sfc_pres, &
1592  env_theta, &
1593  env_u, &
1594  env_v, &
1595  shape_nc, &
1596  bbl_nc
1597 
1598  real(RP), pointer :: shapeFac(:,:,:) => null()
1599 
1600  integer :: k, i, j, iq
1601  integer :: ierr
1602  !---------------------------------------------------------------------------
1603 
1604  if( io_l ) write(io_fid_log,*)
1605  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit TRACERBUBBLE] / Categ[preprocess] / Origin[SCALE-RM]'
1606 
1607  sfc_theta = thetastd
1608  sfc_pres = pstd
1609  env_theta = thetastd
1610 
1611  !--- read namelist
1612  rewind(io_fid_conf)
1613  read(io_fid_conf,nml=param_mkinit_tracerbubble,iostat=ierr)
1614 
1615  if( ierr < 0 ) then !--- missing
1616  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1617  elseif( ierr > 0 ) then !--- fatal error
1618  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_TRACERBUBBLE. Check!'
1619  call prc_mpistop
1620  endif
1621  if( io_nml ) write(io_fid_nml,nml=param_mkinit_tracerbubble)
1622 
1623  ! calc in dry condition
1624  pres_sfc(1,1,1) = sfc_pres
1625  pott_sfc(1,1,1) = sfc_theta
1626  qv_sfc(1,1,1) = 0.0_rp
1627  qc_sfc(1,1,1) = 0.0_rp
1628 
1629  do k = ks, ke
1630  pott(k,1,1) = env_theta
1631  qv(k,1,1) = 0.0_rp
1632  qc(k,1,1) = 0.0_rp
1633  enddo
1634 
1635  ! make density & pressure profile in dry condition
1636  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
1637  temp(:,1,1), & ! [OUT]
1638  pres(:,1,1), & ! [OUT]
1639  pott(:,1,1), & ! [IN]
1640  qv(:,1,1), & ! [IN]
1641  qc(:,1,1), & ! [IN]
1642  temp_sfc(1,1,1), & ! [OUT]
1643  pres_sfc(1,1,1), & ! [IN]
1644  pott_sfc(1,1,1), & ! [IN]
1645  qv_sfc(1,1,1), & ! [IN]
1646  qc_sfc(1,1,1) ) ! [IN]
1647 
1648  do j = js, je
1649  do i = is, ie
1650  do k = ks, ke
1651  dens(k,i,j) = dens(k,1,1)
1652  momz(k,i,j) = 0.0_rp
1653  momx(k,i,j) = env_u * dens(k,1,1)
1654  momy(k,i,j) = env_v * dens(k,1,1)
1655  rhot(k,i,j) = pott(k,1,1) * dens(k,1,1)
1656 
1657  do iq = 1, qa
1658  qtrc(k,i,j,iq) = 0.0_rp
1659  enddo
1660  enddo
1661  enddo
1662  enddo
1663 
1664  ! make tracer bubble
1665  if ( i_nc > 0 ) then
1666 
1667  select case(shape_nc)
1668  case('BUBBLE')
1669  call bubble_setup
1670  shapefac => bubble
1671  case('RECT')
1672  call rect_setup
1673  shapefac => rect
1674  case default
1675  write(*,*) 'xxx SHAPE_NC=', trim(shape_nc), ' cannot be used on advect. Check!'
1676  call prc_mpistop
1677  end select
1678 
1679  do j = js, je
1680  do i = is, ie
1681  do k = ks, ke
1682  qtrc(k,i,j,i_nc) = bbl_nc * shapefac(k,i,j)
1683  enddo
1684  enddo
1685  enddo
1686  else
1687  write(*,*) 'xxx tracer I_NC is not defined. Check!'
1688  call prc_mpistop
1689  endif
1690 
1691  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
1692  write(*,*) 'xxx SBM cannot be used on tracerbubble. Check!'
1693  call prc_mpistop
1694  endif
1695 #endif
1696 
1697  return
1698  end subroutine mkinit_tracerbubble
1699 
1700  !-----------------------------------------------------------------------------
1710  subroutine mkinit_coldbubble
1711  use mod_atmos_admin, only: &
1712  atmos_phy_mp_type
1713  implicit none
1714 
1715  ! Surface state
1716  real(RP) :: SFC_THETA ! surface potential temperature [K]
1717  real(RP) :: SFC_PRES ! surface pressure [Pa]
1718  ! Environment state
1719  real(RP) :: ENV_THETA ! potential temperature of environment [K]
1720  ! Bubble
1721  real(RP) :: BBL_TEMP = -15.0_rp ! extremum of temperature in bubble [K]
1722 
1723  namelist / param_mkinit_coldbubble / &
1724  sfc_theta, &
1725  sfc_pres, &
1726  env_theta, &
1727  bbl_temp
1728 
1729  real(RP) :: RovCP
1730 
1731  integer :: ierr
1732  integer :: k, i, j, iq
1733  !---------------------------------------------------------------------------
1734 
1735  if( io_l ) write(io_fid_log,*)
1736  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit COLDBUBBLE] / Categ[preprocess] / Origin[SCALE-RM]'
1737 
1738  sfc_theta = thetastd
1739  sfc_pres = pstd
1740  env_theta = thetastd
1741 
1742  !--- read namelist
1743  rewind(io_fid_conf)
1744  read(io_fid_conf,nml=param_mkinit_coldbubble,iostat=ierr)
1745 
1746  if( ierr < 0 ) then !--- missing
1747  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1748  elseif( ierr > 0 ) then !--- fatal error
1749  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_COLDBUBBLE. Check!'
1750  call prc_mpistop
1751  endif
1752  if( io_nml ) write(io_fid_nml,nml=param_mkinit_coldbubble)
1753 
1754  rovcp = rdry / cpdry
1755 
1756  ! calc in dry condition
1757  pres_sfc(1,1,1) = sfc_pres
1758  pott_sfc(1,1,1) = sfc_theta
1759  qv_sfc(1,1,1) = 0.0_rp
1760  qc_sfc(1,1,1) = 0.0_rp
1761 
1762  do k = ks, ke
1763  pott(k,1,1) = env_theta
1764  qv(k,1,1) = 0.0_rp
1765  qc(k,1,1) = 0.0_rp
1766  enddo
1767 
1768  ! make density & pressure profile in dry condition
1769  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
1770  temp(:,1,1), & ! [OUT]
1771  pres(:,1,1), & ! [OUT]
1772  pott(:,1,1), & ! [IN]
1773  qv(:,1,1), & ! [IN]
1774  qc(:,1,1), & ! [IN]
1775  temp_sfc(1,1,1), & ! [OUT]
1776  pres_sfc(1,1,1), & ! [IN]
1777  pott_sfc(1,1,1), & ! [IN]
1778  qv_sfc(1,1,1), & ! [IN]
1779  qc_sfc(1,1,1) ) ! [IN]
1780 
1781  do j = 1, ja
1782  do i = 1, ia
1783  do k = ks, ke
1784  dens(k,i,j) = dens(k,1,1)
1785  momz(k,i,j) = 0.0_rp
1786  momx(k,i,j) = 0.0_rp
1787  momy(k,i,j) = 0.0_rp
1788 
1789  ! make cold bubble
1790  rhot(k,i,j) = dens(k,1,1) * ( pott(k,1,1) &
1791  + bbl_temp * ( p00/pres(k,1,1) )**rovcp * bubble(k,i,j) )
1792 
1793  do iq = 1, qa
1794  qtrc(k,i,j,iq) = 0.0_rp
1795  enddo
1796  enddo
1797  enddo
1798  enddo
1799 
1800  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
1801  write(*,*) 'xxx SBM cannot be used on coldbubble. Check!'
1802  call prc_mpistop
1803  endif
1804 
1805  return
1806  end subroutine mkinit_coldbubble
1807 
1808  !-----------------------------------------------------------------------------
1810  subroutine mkinit_lambwave
1811  use mod_atmos_admin, only: &
1812  atmos_phy_mp_type
1813  implicit none
1814 
1815  ! Surface state
1816  real(RP) :: SFC_PRES ! surface pressure [Pa]
1817  ! Environment state
1818  real(RP) :: ENV_U = 0.0_rp ! velocity u of environment [m/s]
1819  real(RP) :: ENV_V = 0.0_rp ! velocity v of environment [m/s]
1820  real(RP) :: ENV_TEMP = 300.0_rp ! temperature of environment [K]
1821  ! Bubble
1822  real(RP) :: BBL_PRES = 100._rp ! extremum of pressure in bubble [Pa]
1823 
1824  namelist / param_mkinit_lambwave / &
1825  sfc_pres, &
1826  env_u, &
1827  env_v, &
1828  env_temp, &
1829  bbl_pres
1830 
1831  real(RP) :: RovCP
1832 
1833  integer :: ierr
1834  integer :: k, i, j, iq
1835  !---------------------------------------------------------------------------
1836 
1837  if( io_l ) write(io_fid_log,*)
1838  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit LAMBWAVE] / Categ[preprocess] / Origin[SCALE-RM]'
1839 
1840  sfc_pres = pstd
1841 
1842  !--- read namelist
1843  rewind(io_fid_conf)
1844  read(io_fid_conf,nml=param_mkinit_lambwave,iostat=ierr)
1845 
1846  if( ierr < 0 ) then !--- missing
1847  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1848  elseif( ierr > 0 ) then !--- fatal error
1849  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_LAMBWAVE. Check!'
1850  call prc_mpistop
1851  endif
1852  if( io_nml ) write(io_fid_nml,nml=param_mkinit_lambwave)
1853 
1854  rovcp = rdry / cpdry
1855 
1856  do j = js, je
1857  do i = is, ie
1858  do k = ks, ke
1859  dens(k,i,j) = sfc_pres/(rdry*env_temp) * exp( - grav/(rdry*env_temp) * grid_cz(k) )
1860  momz(k,i,j) = 0.0_rp
1861  momx(k,i,j) = env_u * dens(k,i,j)
1862  momy(k,i,j) = env_v * dens(k,i,j)
1863 
1864  ! make pressure bubble
1865  pres(k,i,j) = dens(k,i,j) * env_temp * rdry + bbl_pres * bubble(k,i,j)
1866 
1867  rhot(k,i,j) = dens(k,i,j) * env_temp * ( p00/pres(k,i,j) )**rovcp
1868 
1869  do iq = 1, qa
1870  qtrc(k,i,j,iq) = 0.0_rp
1871  enddo
1872  enddo
1873  enddo
1874  enddo
1875 
1876  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
1877  write(*,*) 'xxx SBM cannot be used on lambwave. Check!'
1878  call prc_mpistop
1879  endif
1880 
1881  return
1882  end subroutine mkinit_lambwave
1883 
1884  !-----------------------------------------------------------------------------
1887  subroutine mkinit_gravitywave
1888  use mod_atmos_admin, only: &
1889  atmos_phy_mp_type
1890  implicit none
1891 
1892  ! Surface state
1893  real(RP) :: SFC_THETA ! surface potential temperature [K]
1894  real(RP) :: SFC_PRES ! surface pressure [Pa]
1895  ! Environment state
1896  real(RP) :: ENV_U = 20.0_rp ! velocity u of environment [m/s]
1897  real(RP) :: ENV_V = 0.0_rp ! velocity v of environment [m/s]
1898  real(RP) :: ENV_BVF = 0.01_rp ! Brunt Vaisala frequencies of environment [1/s]
1899  ! Bubble
1900  real(RP) :: BBL_THETA = 0.01_rp ! extremum of potential temperature in bubble [K]
1901 
1902  namelist / param_mkinit_gravitywave / &
1903  sfc_theta, &
1904  sfc_pres, &
1905  env_u, &
1906  env_v, &
1907  env_bvf, &
1908  bbl_theta
1909 
1910  integer :: ierr
1911  integer :: k, i, j, iq
1912  !---------------------------------------------------------------------------
1913 
1914  if( io_l ) write(io_fid_log,*)
1915  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit GRAVITYWAVE] / Categ[preprocess] / Origin[SCALE-RM]'
1916 
1917  sfc_theta = thetastd
1918  sfc_pres = pstd
1919 
1920  !--- read namelist
1921  rewind(io_fid_conf)
1922  read(io_fid_conf,nml=param_mkinit_gravitywave,iostat=ierr)
1923 
1924  if( ierr < 0 ) then !--- missing
1925  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
1926  elseif( ierr > 0 ) then !--- fatal error
1927  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_GRAVITYWAVE. Check!'
1928  call prc_mpistop
1929  endif
1930  if( io_nml ) write(io_fid_nml,nml=param_mkinit_gravitywave)
1931 
1932  ! calc in dry condition
1933  pres_sfc(1,1,1) = sfc_pres
1934  pott_sfc(1,1,1) = sfc_theta
1935  qv_sfc(1,1,1) = 0.0_rp
1936  qc_sfc(1,1,1) = 0.0_rp
1937 
1938  do k = ks, ke
1939  pott(k,1,1) = sfc_theta * exp( env_bvf*env_bvf / grav * grid_cz(k) )
1940  qv(k,1,1) = 0.0_rp
1941  qc(k,1,1) = 0.0_rp
1942  enddo
1943 
1944  ! make density & pressure profile in dry condition
1945  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
1946  temp(:,1,1), & ! [OUT]
1947  pres(:,1,1), & ! [OUT]
1948  pott(:,1,1), & ! [IN]
1949  qv(:,1,1), & ! [IN]
1950  qc(:,1,1), & ! [IN]
1951  temp_sfc(1,1,1), & ! [OUT]
1952  pres_sfc(1,1,1), & ! [IN]
1953  pott_sfc(1,1,1), & ! [IN]
1954  qv_sfc(1,1,1), & ! [IN]
1955  qc_sfc(1,1,1) ) ! [IN]
1956 
1957  do j = js, je
1958  do i = is, ie
1959  do k = ks, ke
1960  dens(k,i,j) = dens(k,1,1)
1961  momz(k,i,j) = 0.0_rp
1962  momx(k,i,j) = env_u * dens(k,1,1)
1963  momy(k,i,j) = env_v * dens(k,1,1)
1964 
1965  ! make warm bubble
1966  rhot(k,i,j) = dens(k,1,1) * ( pott(k,1,1) + bbl_theta * bubble(k,i,j) )
1967 
1968  do iq = 1, qa
1969  qtrc(k,i,j,iq) = 0.0_rp
1970  enddo
1971  enddo
1972  enddo
1973  enddo
1974 
1975  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
1976  write(*,*) 'xxx SBM cannot be used on gravitywave. Check!'
1977  call prc_mpistop
1978  endif
1979 
1980  return
1981  end subroutine mkinit_gravitywave
1982 
1983  !-----------------------------------------------------------------------------
1985  subroutine mkinit_khwave
1986  use mod_atmos_admin, only: &
1987  atmos_phy_mp_type
1988  implicit none
1989 
1990  ! Surface state
1991  real(RP) :: SFC_THETA ! surface potential temperature [K]
1992  real(RP) :: SFC_PRES ! surface pressure [Pa]
1993  ! Environment state
1994  real(RP) :: ENV_L1_ZTOP = 1900.0_rp ! top height of the layer1 (low THETA) [m]
1995  real(RP) :: ENV_L3_ZBOTTOM = 2100.0_rp ! bottom height of the layer3 (high THETA) [m]
1996  real(RP) :: ENV_L1_THETA = 300.0_rp ! THETA in the layer1 (low THETA) [K]
1997  real(RP) :: ENV_L3_THETA = 301.0_rp ! THETA in the layer3 (high THETA) [K]
1998  real(RP) :: ENV_L1_U = 0.0_rp ! velocity u in the layer1 (low THETA) [K]
1999  real(RP) :: ENV_L3_U = 20.0_rp ! velocity u in the layer3 (high THETA) [K]
2000  ! Disturbance
2001  real(RP) :: RANDOM_U = 0.0_rp ! amplitude of random disturbance u
2002 
2003  namelist / param_mkinit_khwave / &
2004  sfc_theta, &
2005  sfc_pres, &
2006  env_l1_ztop, &
2007  env_l3_zbottom, &
2008  env_l1_theta, &
2009  env_l3_theta, &
2010  env_l1_u, &
2011  env_l3_u, &
2012  random_u
2013 
2014  real(RP) :: fact
2015 
2016  integer :: ierr
2017  integer :: k, i, j, iq
2018  !---------------------------------------------------------------------------
2019 
2020  if( io_l ) write(io_fid_log,*)
2021  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit KHWAVE] / Categ[preprocess] / Origin[SCALE-RM]'
2022 
2023  sfc_theta = thetastd
2024  sfc_pres = pstd
2025 
2026  !--- read namelist
2027  rewind(io_fid_conf)
2028  read(io_fid_conf,nml=param_mkinit_khwave,iostat=ierr)
2029 
2030  if( ierr < 0 ) then !--- missing
2031  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
2032  elseif( ierr > 0 ) then !--- fatal error
2033  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_KHWAVE. Check!'
2034  call prc_mpistop
2035  endif
2036  if( io_nml ) write(io_fid_nml,nml=param_mkinit_khwave)
2037 
2038  ! calc in dry condition
2039  pres_sfc(1,1,1) = sfc_pres
2040  pott_sfc(1,1,1) = sfc_theta
2041  qv_sfc(1,1,1) = 0.0_rp
2042  qc_sfc(1,1,1) = 0.0_rp
2043 
2044  do k = ks, ke
2045  fact = ( grid_cz(k)-env_l1_ztop ) / ( env_l3_zbottom-env_l1_ztop )
2046  fact = max( min( fact, 1.0_rp ), 0.0_rp )
2047 
2048  pott(k,1,1) = env_l1_theta * ( 1.0_rp - fact ) &
2049  + env_l3_theta * ( fact )
2050 
2051  qv(k,1,1) = 0.0_rp
2052  qc(k,1,1) = 0.0_rp
2053  enddo
2054 
2055  ! make density & pressure profile in dry condition
2056  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
2057  temp(:,1,1), & ! [OUT]
2058  pres(:,1,1), & ! [OUT]
2059  pott(:,1,1), & ! [IN]
2060  qv(:,1,1), & ! [IN]
2061  qc(:,1,1), & ! [IN]
2062  temp_sfc(1,1,1), & ! [OUT]
2063  pres_sfc(1,1,1), & ! [IN]
2064  pott_sfc(1,1,1), & ! [IN]
2065  qv_sfc(1,1,1), & ! [IN]
2066  qc_sfc(1,1,1) ) ! [IN]
2067 
2068  do j = js, je
2069  do i = is, ie
2070  do k = ks, ke
2071  dens(k,i,j) = dens(k,1,1)
2072  momz(k,i,j) = 0.0_rp
2073  momy(k,i,j) = 0.0_rp
2074  rhot(k,i,j) = dens(k,1,1) * pott(k,1,1)
2075 
2076  do iq = 1, qa
2077  qtrc(k,i,j,iq) = 0.0_rp
2078  enddo
2079  enddo
2080  enddo
2081  enddo
2082 
2083  call random_get(rndm) ! make random
2084  do j = js, je
2085  do i = is, ie
2086  do k = ks, ke
2087  fact = ( grid_cz(k)-env_l1_ztop ) / ( env_l3_zbottom-env_l1_ztop )
2088  fact = max( min( fact, 1.0_rp ), 0.0_rp )
2089 
2090  momx(k,i,j) = ( env_l1_u * ( 1.0_rp - fact ) &
2091  + env_l3_u * ( fact ) &
2092  + ( rndm(k,i,j) - 0.5_rp ) * 2.0_rp * random_u &
2093  ) * dens(k,i,j)
2094  enddo
2095  enddo
2096  enddo
2097 
2098  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
2099  write(*,*) 'xxx SBM cannot be used on khwave. Check!'
2100  call prc_mpistop
2101  endif
2102 
2103  return
2104  end subroutine mkinit_khwave
2105 
2106  !-----------------------------------------------------------------------------
2108  subroutine mkinit_turbulence
2109  use scale_atmos_hydrometeor, only: &
2110  i_qv
2111  use mod_atmos_admin, only: &
2112  atmos_phy_mp_type
2113  implicit none
2114 
2115  ! Surface state
2116  real(RP) :: SFC_THETA ! surface potential temperature [K]
2117  real(RP) :: SFC_PRES ! surface pressure [Pa]
2118  real(RP) :: SFC_RH = 0.0_rp ! surface relative humidity [%]
2119  ! Environment state
2120  real(RP) :: ENV_THETA ! potential temperature of environment [K]
2121  real(RP) :: ENV_TLAPS = 4.e-3_rp ! Lapse rate of THETA [K/m]
2122  real(RP) :: ENV_U = 5.0_rp ! velocity u of environment [m/s]
2123  real(RP) :: ENV_V = 0.0_rp ! velocity v of environment [m/s]
2124  real(RP) :: ENV_RH = 0.0_rp ! relative humidity of environment [%]
2125  ! Disturbance
2126  real(RP) :: RANDOM_THETA = 1.0_rp ! amplitude of random disturbance theta
2127  real(RP) :: RANDOM_U = 0.0_rp ! amplitude of random disturbance u
2128  real(RP) :: RANDOM_V = 0.0_rp ! amplitude of random disturbance v
2129  real(RP) :: RANDOM_RH = 0.0_rp ! amplitude of random disturbance RH
2130 
2131  namelist / param_mkinit_turbulence / &
2132  sfc_theta, &
2133  sfc_pres, &
2134  sfc_rh, &
2135  env_theta, &
2136  env_tlaps, &
2137  env_u, &
2138  env_v, &
2139  env_rh, &
2140  random_theta, &
2141  random_u, &
2142  random_v, &
2143  random_rh
2144 
2145  integer :: ierr
2146  integer :: k, i, j, iq
2147  !---------------------------------------------------------------------------
2148 
2149  if( io_l ) write(io_fid_log,*)
2150  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit TURBULENCE] / Categ[preprocess] / Origin[SCALE-RM]'
2151 
2152  sfc_theta = thetastd
2153  sfc_pres = pstd
2154  env_theta = thetastd
2155 
2156  !--- read namelist
2157  rewind(io_fid_conf)
2158  read(io_fid_conf,nml=param_mkinit_turbulence,iostat=ierr)
2159 
2160  if( ierr < 0 ) then !--- missing
2161  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
2162  elseif( ierr > 0 ) then !--- fatal error
2163  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_TURBULENCE. Check!'
2164  call prc_mpistop
2165  endif
2166  if( io_nml ) write(io_fid_nml,nml=param_mkinit_turbulence)
2167 
2168  ! calc in dry condition
2169  pres_sfc(1,1,1) = sfc_pres
2170  pott_sfc(1,1,1) = sfc_theta
2171  qv_sfc(1,:,:) = 0.0_rp
2172  qc_sfc(1,:,:) = 0.0_rp
2173 
2174  do k = ks, ke
2175  pott(k,1,1) = env_theta + env_tlaps * grid_cz(k)
2176  qv(k,:,:) = 0.0_rp
2177  qc(k,:,:) = 0.0_rp
2178  enddo
2179 
2180  ! make density & pressure profile in dry condition
2181  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
2182  temp(:,1,1), & ! [OUT]
2183  pres(:,1,1), & ! [OUT]
2184  pott(:,1,1), & ! [IN]
2185  qv(:,1,1), & ! [IN]
2186  qc(:,1,1), & ! [IN]
2187  temp_sfc(1,1,1), & ! [OUT]
2188  pres_sfc(1,1,1), & ! [IN]
2189  pott_sfc(1,1,1), & ! [IN]
2190  qv_sfc(1,1,1), & ! [IN]
2191  qc_sfc(1,1,1) ) ! [IN]
2192 
2193  if ( i_qv > 0 ) then
2194  ! calc QV from RH
2195  call saturation_pres2qsat_all( qsat_sfc(1,1,1), temp_sfc(1,1,1), pres_sfc(1,1,1) )
2196  call saturation_pres2qsat_all( qsat(:,1,1), temp(:,1,1), pres(:,1,1) )
2197 
2198  call random_get(rndm) ! make random
2199  do j = js, je
2200  do i = is, ie
2201  qv_sfc(1,i,j) = ( sfc_rh + rndm(ks-1,i,j) * random_rh ) * 1.e-2_rp * qsat_sfc(1,1,1)
2202 
2203  do k = ks, ke
2204  qv(k,i,j) = ( env_rh + rndm(k,i,j) * random_rh ) * 1.e-2_rp * qsat(k,1,1)
2205  enddo
2206  enddo
2207  enddo
2208  end if
2209 
2210  call random_get(rndm) ! make random
2211  do j = js, je
2212  do i = is, ie
2213  pres_sfc(1,i,j) = sfc_pres
2214  pott_sfc(1,i,j) = sfc_theta + rndm(ks-1,i,j) * random_theta
2215 
2216  do k = ks, ke
2217  pott(k,i,j) = env_theta + env_tlaps * grid_cz(k) + rndm(k,i,j) * random_theta
2218  enddo
2219  enddo
2220  enddo
2221 
2222  ! make density & pressure profile in moist condition
2223  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
2224  temp(:,:,:), & ! [OUT]
2225  pres(:,:,:), & ! [OUT]
2226  pott(:,:,:), & ! [IN]
2227  qv(:,:,:), & ! [IN]
2228  qc(:,:,:), & ! [IN]
2229  temp_sfc(:,:,:), & ! [OUT]
2230  pres_sfc(:,:,:), & ! [IN]
2231  pott_sfc(:,:,:), & ! [IN]
2232  qv_sfc(:,:,:), & ! [IN]
2233  qc_sfc(:,:,:) ) ! [IN]
2234 
2235  call comm_vars8( dens(:,:,:), 1 )
2236  call comm_wait ( dens(:,:,:), 1 )
2237 
2238  call random_get(rndm) ! make random
2239  do j = js, je
2240  do i = is, ie
2241  do k = ks, ke
2242  momx(k,i,j) = ( env_u + ( rndm(k,i,j) - 0.5_rp ) * 2.0_rp * random_u ) &
2243  * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
2244  enddo
2245  enddo
2246  enddo
2247 
2248  call random_get(rndm) ! make random
2249  do j = js, je
2250  do i = is, ie
2251  do k = ks, ke
2252  momy(k,i,j) = ( env_v + ( rndm(k,i,j) - 0.5_rp ) * 2.0_rp * random_v ) &
2253  * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
2254  enddo
2255  enddo
2256  enddo
2257 
2258  do j = js, je
2259  do i = is, ie
2260  do k = ks, ke
2261  momz(k,i,j) = 0.0_rp
2262  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
2263  enddo
2264  enddo
2265  enddo
2266 
2267  if ( i_qv > 0.0_rp ) then
2268  do j = js, je
2269  do i = is, ie
2270  do k = ks, ke
2271  qtrc(k,i,j,i_qv) = qv(k,i,j)
2272  enddo
2273  enddo
2274  enddo
2275  end if
2276 
2277 #ifndef DRY
2278  do iq = 1, qa
2279  if ( iq == i_qv ) cycle
2280  do j = js, je
2281  do i = is, ie
2282  do k = ks, ke
2283  qtrc(k,i,j,iq) = 0.0_rp
2284  enddo
2285  enddo
2286  enddo
2287  enddo
2288 #endif
2289 
2290  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
2291  write(*,*) 'xxx SBM cannot be used on turbulence. Check!'
2292  call prc_mpistop
2293  endif
2294 
2295  return
2296  end subroutine mkinit_turbulence
2297 
2298  !-----------------------------------------------------------------------------
2300  subroutine mkinit_cavityflow
2301  implicit none
2302 
2303  ! Nondimenstional numbers for a cavity flow problem
2304  real(RP) :: REYNOLDS_NUM = 1.d03
2305  real(RP) :: MACH_NUM = 3.d-2
2306  real(RP) :: Ulid = 1.d01
2307  real(RP) :: PRES0 = 1.d05
2308 
2309  namelist / param_mkinit_cavityflow / &
2310  ulid , &
2311  pres0 , &
2312  reynolds_num, &
2313  mach_num
2314 
2315  real(RP) :: DENS0
2316  real(RP) :: TEMP
2317  real(RP) :: Gam
2318  real(RP) :: Cs2
2319 
2320  integer :: k, i, j
2321  integer :: ierr
2322  !---------------------------------------------------------------------------
2323 
2324  if( io_l ) write(io_fid_log,*)
2325  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit CAVITYFLOW] / Categ[preprocess] / Origin[SCALE-RM]'
2326 
2327  !--- read namelist
2328  rewind(io_fid_conf)
2329  read(io_fid_conf,nml=param_mkinit_cavityflow,iostat=ierr)
2330 
2331  if( ierr < 0 ) then !--- missing
2332  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
2333  elseif( ierr > 0 ) then !--- fatal error
2334  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_CAVITYFLOW. Check!'
2335  call prc_mpistop
2336  endif
2337  if( io_nml ) write(io_fid_nml,nml=param_mkinit_cavityflow)
2338 
2339  gam = cpdry / ( cpdry - rdry )
2340  cs2 = ( ulid / mach_num )**2
2341  temp = cs2 / ( gam * rdry )
2342  dens0 = pres0 / ( rdry * temp )
2343 
2344  if( io_l ) write(io_fid_log,*) "DENS = ", dens0
2345  if( io_l ) write(io_fid_log,*) "PRES = ", pres0
2346  if( io_l ) write(io_fid_log,*) "TEMP = ", rhot(10,10,4)/dens0, temp
2347  if( io_l ) write(io_fid_log,*) "Ulid = ", ulid
2348  if( io_l ) write(io_fid_log,*) "Cs = ", sqrt(cs2)
2349 
2350  do j = 1, ja
2351  do i = 1, ia
2352  do k = ks, ke
2353  dens(k,i,j) = dens0
2354  momz(k,i,j) = 0.0_rp
2355  momx(k,i,j) = 0.0_rp
2356  momy(k,i,j) = 0.0_rp
2357  pres(k,i,j) = pres0
2358  rhot(k,i,j) = p00/rdry * (p00/pres0)**((rdry - cpdry)/cpdry)
2359  qtrc(k,i,j,:) = 0.0_rp
2360  enddo
2361  enddo
2362  enddo
2363 
2364  momx(ke+1:ka,:,:) = dens0 * ulid
2365 
2366  return
2367  end subroutine mkinit_cavityflow
2368 
2369  !-----------------------------------------------------------------------------
2371  subroutine mkinit_mountainwave
2372  use scale_atmos_hydrometeor, only: &
2373  i_nc
2374  implicit none
2375 
2376  ! Surface state
2377  real(RP) :: SFC_THETA ! surface potential temperature [K]
2378  real(RP) :: SFC_PRES ! surface pressure [Pa]
2379  ! Environment state
2380  real(RP) :: ENV_U = 0.0_rp ! velocity u of environment [m/s]
2381  real(RP) :: ENV_V = 0.0_rp ! velocity v of environment [m/s]
2382 
2383  real(RP) :: SCORER = 2.e-3_rp ! Scorer parameter (~=N/U) [1/m]
2384  real(RP) :: BBL_NC = 0.0_rp ! extremum of NC in bubble [kg/kg]
2385 
2386  namelist / param_mkinit_mountainwave / &
2387  sfc_theta, &
2388  sfc_pres, &
2389  env_u, &
2390  env_v, &
2391  scorer, &
2392  bbl_nc
2393 
2394  real(RP) :: Ustar2, N2
2395 
2396  integer :: ierr
2397  integer :: k, i, j
2398  !---------------------------------------------------------------------------
2399 
2400  if( io_l ) write(io_fid_log,*)
2401  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit MOUNTAINWAVE] / Categ[preprocess] / Origin[SCALE-RM]'
2402 
2403  sfc_theta = thetastd
2404  sfc_pres = pstd
2405 
2406  !--- read namelist
2407  rewind(io_fid_conf)
2408  read(io_fid_conf,nml=param_mkinit_mountainwave,iostat=ierr)
2409 
2410  if( ierr < 0 ) then !--- missing
2411  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
2412  elseif( ierr > 0 ) then !--- fatal error
2413  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_MOUNTAINWAVE. Check!'
2414  call prc_mpistop
2415  endif
2416  if( io_nml ) write(io_fid_nml,nml=param_mkinit_mountainwave)
2417 
2418  ! calc in dry condition
2419  do j = js, je
2420  do i = is, ie
2421  pres_sfc(1,i,j) = sfc_pres
2422  pott_sfc(1,i,j) = sfc_theta
2423  qv_sfc(1,i,j) = 0.0_rp
2424  qc_sfc(1,i,j) = 0.0_rp
2425  enddo
2426  enddo
2427 
2428  do j = js, je
2429  do i = is, ie
2430  do k = ks, ke
2431  ustar2 = env_u * env_u + env_v * env_v
2432  n2 = ustar2 * (scorer*scorer)
2433 
2434  pott(k,i,j) = sfc_theta * exp( n2 / grav * real_cz(k,i,j) )
2435  qv(k,i,j) = 0.0_rp
2436  qc(k,i,j) = 0.0_rp
2437  enddo
2438  enddo
2439  enddo
2440 
2441  ! make density & pressure profile in dry condition
2442  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
2443  temp(:,:,:), & ! [OUT]
2444  pres(:,:,:), & ! [OUT]
2445  pott(:,:,:), & ! [IN]
2446  qv(:,:,:), & ! [IN]
2447  qc(:,:,:), & ! [IN]
2448  temp_sfc(:,:,:), & ! [OUT]
2449  pres_sfc(:,:,:), & ! [IN]
2450  pott_sfc(:,:,:), & ! [IN]
2451  qv_sfc(:,:,:), & ! [IN]
2452  qc_sfc(:,:,:) ) ! [IN]
2453 
2454  do j = js, je
2455  do i = is, ie
2456  do k = ks, ke
2457  dens(k,i,j) = dens(k,i,j)
2458  momz(k,i,j) = 0.0_rp
2459  momx(k,i,j) = env_u * dens(k,i,j)
2460  momy(k,i,j) = env_v * dens(k,i,j)
2461  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
2462 
2463  qtrc(k,i,j,:) = 0.0_rp
2464  enddo
2465  enddo
2466  enddo
2467 
2468  ! optional : add tracer bubble
2469  if ( bbl_nc > 0.0_rp ) then
2470  if ( i_nc > 0 ) then
2471  do j = js, je
2472  do i = is, ie
2473  do k = ks, ke
2474  qtrc(k,i,j,i_nc) = bbl_nc * bubble(k,i,j)
2475  enddo
2476  enddo
2477  enddo
2478  else
2479  write(*,*) 'xxx tracer I_NC is not defined. Check!'
2480  call prc_mpistop
2481  endif
2482  endif
2483 
2484  return
2485  end subroutine mkinit_mountainwave
2486 
2487  !-----------------------------------------------------------------------------
2492  subroutine mkinit_barocwave
2493  use scale_const, only: &
2494  ohm => const_ohm, &
2495  rplanet => const_radius, &
2496  grav => const_grav
2497  use scale_process
2498  use scale_grid, only: &
2499  y0 => grid_domain_center_y, &
2500  grid_fyg
2501  use scale_atmos_hydrometeor, only: &
2502  i_qv
2503 
2504  implicit none
2505 
2506  ! Parameters for global domain size
2507  real(RP) :: Ly ! The domain size in y-direction [m]
2508 
2509  ! Parameters for inital stratification
2510  real(RP) :: REF_TEMP = 288.e0_rp ! The reference temperature [K]
2511  real(RP) :: REF_PRES = 1.e5_rp ! The reference pressure [Pa]
2512  real(RP) :: LAPSE_RATE = 5.e-3_rp ! The lapse rate [K/m]
2513 
2514  ! Parameters associated with coriolis parameter on a beta-plane
2515  real(RP) :: Phi0Deg = 45.e0_rp ! The central latitude [degree_north]
2516 
2517  ! Parameters for background zonal jet
2518  real(RP) :: U0 = 35.e0_rp ! The parameter associated with zonal jet maximum amplitude [m/s]
2519  real(RP) :: b = 2.e0_rp ! The vertical half-width [1]
2520 
2521  ! Parameters for inital perturbation of zonal wind with a Gaussian profile
2522  !
2523  real(RP) :: Up = 1.e0_rp ! The maximum amplitude of zonal wind perturbation [m/s]
2524  real(RP) :: Lp = 600.e3_rp ! The width of Gaussian profile
2525  real(RP) :: Xc = 2000.e3_rp ! The center point (x) of inital perturbation
2526  real(RP) :: Yc = 2500.e3_rp ! The center point (y) of inital perturbation
2527 
2528  namelist / param_mkinit_barocwave / &
2529  ref_temp, ref_pres, lapse_rate, &
2530  phi0deg, &
2531  u0, b, &
2532  up, lp, xc, yc
2533 
2534  real(RP) :: f0, beta0
2535 
2536  real(RP) :: geopot(KA,IA,JA)
2537  real(RP) :: eta(KA,IA,JA)
2538  real(RP) :: temp(KA,IA,JA)
2539 
2540  real(RP) :: y
2541  real(RP) :: ln_eta
2542  real(RP) :: del_eta
2543  real(RP) :: yphase
2544  real(RP) :: temp_vfunc
2545  real(RP) :: geopot_hvari
2546 
2547  integer :: ierr
2548  integer :: k, i, j
2549 
2550  integer :: itr
2551 
2552  integer, parameter :: ITRMAX = 1000
2553  real(RP), parameter :: CONV_EPS = 1e-15_rp
2554  !---------------------------------------------------------------------------
2555 
2556  if( io_l ) write(io_fid_log,*)
2557  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit BAROCWAVE] / Categ[preprocess] / Origin[SCALE-RM]'
2558 
2559 
2560  !--- read namelist
2561  rewind(io_fid_conf)
2562  read(io_fid_conf,nml=param_mkinit_barocwave,iostat=ierr)
2563 
2564  if( ierr < 0 ) then !--- missing
2565  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
2566  elseif( ierr > 0 ) then !--- fatal error
2567  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_BAROCWAVE. Check!'
2568  call prc_mpistop
2569  endif
2570  if( io_nml ) write(io_fid_nml,nml=param_mkinit_barocwave)
2571 
2572  ly = grid_fyg(jag-jhalo) - grid_fyg(jhalo)
2573 
2574  ! Set coriolis parameters
2575  f0 = 2.0_rp*ohm*sin(phi0deg*pi/180.0_rp)
2576  beta0 = (2.0_rp*ohm/rplanet)*cos(phi0deg*pi/180.0_rp)
2577 
2578  ! Calculate eta(=p/p_s) level corresponding to z level of each (y,z) grid point
2579  ! using Newton's iteration method
2580 
2581  eta(:,:,:) = 1.0e-8_rp ! Set first guess of eta
2582 
2583  do j = js, je
2584  do i = is, is ! Note that initial fields are zonaly symmetric
2585 
2586  y = grid_cy(j)
2587  yphase = 2.0_rp*pi*y/ly
2588 
2589  ! Calc horizontal variation of geopotential height
2590  geopot_hvari = 0.5_rp*u0*( &
2591  (f0 - beta0*y0)*(y - 0.5_rp*ly*(1.0_rp + sin(yphase)/pi)) &
2592  + 0.5_rp*beta0*( y**2 - ly*y/pi*sin(yphase) - 0.5_rp*(ly/pi)**2*(cos(yphase) + 1.0_rp) &
2593  - ly**2/3.0_rp ) &
2594  )
2595 
2596  ! Set surface pressure and temperature
2597  pres_sfc(1,i,j) = ref_pres
2598  pott_sfc(1,i,j) = ref_temp - geopot_hvari/rdry
2599  ! Dry condition
2600  qv(:,i,j) = 0.0_rp
2601  qv_sfc(1,i,j) = 0.0_rp
2602  qc(:,i,j) = 0.0_rp
2603  qc_sfc(1,i,j) = 0.0_rp
2604 
2605  do k = ks, ke
2606  del_eta = 1.0_rp
2607 
2608  !-- The loop for iteration
2609  itr = 0
2610  do while( abs(del_eta) > conv_eps )
2611  ln_eta = log(eta(k,i,j))
2612 
2613  temp_vfunc = eta(k,i,j)**(rdry*lapse_rate/grav)
2614  temp(k,i,j) = &
2615  ref_temp*temp_vfunc &
2616  + geopot_hvari/rdry*(2.0_rp*(ln_eta/b)**2 - 1.0_rp)*exp(-(ln_eta/b)**2)
2617  geopot(k,i,j) = &
2618  ref_temp*grav/lapse_rate*(1.0_rp - temp_vfunc) &
2619  + geopot_hvari*ln_eta*exp(-(ln_eta/b)**2)
2620 
2621  del_eta = - ( - grav*grid_cz(k) + geopot(k,i,j) ) & ! <- F
2622  & *( - eta(k,i,j)/(rdry*temp(k,i,j)) ) ! <- (dF/deta)^-1
2623 
2624  eta(k,i,j) = eta(k,i,j) + del_eta
2625  itr = itr + 1
2626 
2627  if ( itr > itrmax ) then
2628  write(*,*) "* (X,Y,Z)=", grid_cx(i), grid_cy(j), grid_cz(k)
2629  write(*,*) "Fail the convergence of iteration. Check!"
2630  write(*,*) "itr=", itr, "del_eta=", del_eta, "eta=", eta(k,i,j), "temp=", temp(k,i,j)
2631  call prc_mpistop
2632  end if
2633  enddo !- End of loop for iteration ----------------------------
2634 
2635  pres(k,i,j) = eta(k,i,j)*ref_pres
2636  dens(k,i,j) = pres(k,i,j)/(rdry*temp(k,i,j))
2637  pott(k,i,j) = temp(k,i,j)*eta(k,i,j)**(-rdry/cpdry)
2638 
2639  enddo
2640 
2641  ! Make density & pressure profile in dry condition using the profile of
2642  ! potential temperature calculated above.
2643  call hydrostatic_buildrho( dens(:,i,j), & ! [OUT]
2644  temp(:,i,j), & ! [OUT]
2645  pres(:,i,j), & ! [OUT]
2646  pott(:,i,j), & ! [IN]
2647  qv(:,i,j), & ! [IN]
2648  qc(:,i,j), & ! [IN]
2649  temp_sfc(1,i,j), & ! [OUT]
2650  pres_sfc(1,i,j), & ! [IN]
2651  pott_sfc(1,i,j), & ! [IN]
2652  qv_sfc(1,i,j), & ! [IN]
2653  qc_sfc(1,i,j) ) ! [IN]
2654  enddo
2655  enddo
2656 
2657  !-----------------------------------------------------------------------------------
2658 
2659  do j = js, je
2660  do k = ks, ke
2661 
2662  eta(k,is,j) = pres(k,is,j)/ref_pres
2663  ln_eta = log(eta(k,is,j))
2664  yphase = 2.0_rp*pi*grid_cy(j)/ly
2665 !!$ PRES(k,IS:IE,j) = eta(k,IS,j)*REF_PRES
2666 !!$ DENS(k,IS:IE,j) = PRES(k,IS,j)/(Rdry*temp(k,IS,j))
2667  dens(k,is:ie,j) = dens(k,is,j)
2668  pres(k,is:ie,j) = pres(k,is,j)
2669  momx(k,is-1:ie,j) = dens(k,is,j)*(-u0*sin(0.5_rp*yphase)**2*ln_eta*exp(-(ln_eta/b)**2))
2670  rhot(k,is:ie,j) = dens(k,is,j)*pott(k,is,j) !temp(k,IS,j)*eta(k,IS,j)**(-Rdry/CPdry)
2671  if ( i_qv > 0 ) qtrc(k,is:ie,j,i_qv) = 0.0_rp
2672  enddo
2673  enddo
2674  momy(:,:,:) = 0.0_rp
2675  momz(:,:,:) = 0.0_rp
2676 
2677  !---------------------------------------------------------------------------------------
2678 
2679  ! Add the inital perturbation for zonal velocity
2680  do j = js, je
2681  do i = is-1, ie
2682  momx(ks:ke,i,j) = momx(ks:ke,i,j) &
2683  + dens(ks:ke,i,j)* up*exp( - ((grid_fx(i) - xc)**2 + (grid_cy(j) - yc)**2)/lp**2 )
2684  enddo
2685  enddo
2686 
2687  return
2688  end subroutine mkinit_barocwave
2689 
2690  !-----------------------------------------------------------------------------
2692  subroutine mkinit_warmbubble
2693  use scale_atmos_hydrometeor, only: &
2694  i_qv
2695  implicit none
2696 
2697  ! Surface state
2698  real(RP) :: SFC_THETA ! surface potential temperature [K]
2699  real(RP) :: SFC_PRES ! surface pressure [Pa]
2700  real(RP) :: SFC_RH = 80.0_rp ! surface relative humidity [%]
2701  ! Environment state
2702  real(RP) :: ENV_U = 0.0_rp ! velocity u of environment [m/s]
2703  real(RP) :: ENV_V = 0.0_rp ! velocity v of environment [m/s]
2704  real(RP) :: ENV_RH = 80.0_rp ! Relative Humidity of environment [%]
2705  real(RP) :: ENV_L1_ZTOP = 1.e3_rp ! top height of the layer1 (constant THETA) [m]
2706  real(RP) :: ENV_L2_ZTOP = 14.e3_rp ! top height of the layer2 (small THETA gradient) [m]
2707  real(RP) :: ENV_L2_TLAPS = 4.e-3_rp ! Lapse rate of THETA in the layer2 (small THETA gradient) [K/m]
2708  real(RP) :: ENV_L3_TLAPS = 3.e-2_rp ! Lapse rate of THETA in the layer3 (large THETA gradient) [K/m]
2709  ! Bubble
2710  real(RP) :: BBL_THETA = 1.0_rp ! extremum of temperature in bubble [K]
2711 
2712  namelist / param_mkinit_warmbubble / &
2713  sfc_theta, &
2714  sfc_pres, &
2715  env_u, &
2716  env_v, &
2717  env_rh, &
2718  env_l1_ztop, &
2719  env_l2_ztop, &
2720  env_l2_tlaps, &
2721  env_l3_tlaps, &
2722  bbl_theta
2723 
2724  integer :: ierr
2725  integer :: k, i, j
2726  !---------------------------------------------------------------------------
2727 
2728  if( io_l ) write(io_fid_log,*)
2729  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit WARMBUBBLE] / Categ[preprocess] / Origin[SCALE-RM]'
2730 
2731  if ( i_qv < 1 ) then
2732  write(*,*) 'xxx QV is not registered'
2733  call prc_mpistop
2734  end if
2735 
2736  sfc_theta = thetastd
2737  sfc_pres = pstd
2738 
2739  !--- read namelist
2740  rewind(io_fid_conf)
2741  read(io_fid_conf,nml=param_mkinit_warmbubble,iostat=ierr)
2742 
2743  if( ierr < 0 ) then !--- missing
2744  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
2745  elseif( ierr > 0 ) then !--- fatal error
2746  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_WARMBUBBLE. Check!'
2747  call prc_mpistop
2748  endif
2749  if( io_nml ) write(io_fid_nml,nml=param_mkinit_warmbubble)
2750 
2751  ! calc in dry condition
2752  pres_sfc(1,1,1) = sfc_pres
2753  pott_sfc(1,1,1) = sfc_theta
2754  qv_sfc(1,:,:) = 0.0_rp
2755  qc_sfc(1,:,:) = 0.0_rp
2756 
2757  do k = ks, ke
2758  if ( grid_cz(k) <= env_l1_ztop ) then ! Layer 1
2759  pott(k,1,1) = sfc_theta
2760  elseif( grid_cz(k) < env_l2_ztop ) then ! Layer 2
2761  pott(k,1,1) = pott(k-1,1,1) + env_l2_tlaps * ( grid_cz(k)-grid_cz(k-1) )
2762  else ! Layer 3
2763  pott(k,1,1) = pott(k-1,1,1) + env_l3_tlaps * ( grid_cz(k)-grid_cz(k-1) )
2764  endif
2765  qv(k,1,1) = 0.0_rp
2766  qc(k,1,1) = 0.0_rp
2767  enddo
2768 
2769  ! make density & pressure profile in dry condition
2770  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
2771  temp(:,1,1), & ! [OUT]
2772  pres(:,1,1), & ! [OUT]
2773  pott(:,1,1), & ! [IN]
2774  qv(:,1,1), & ! [IN]
2775  qc(:,1,1), & ! [IN]
2776  temp_sfc(1,1,1), & ! [OUT]
2777  pres_sfc(1,1,1), & ! [IN]
2778  pott_sfc(1,1,1), & ! [IN]
2779  qv_sfc(1,1,1), & ! [IN]
2780  qc_sfc(1,1,1) ) ! [IN]
2781 
2782  ! calc QV from RH
2783  call saturation_pres2qsat_all( qsat_sfc(1,1,1), temp_sfc(1,1,1), pres_sfc(1,1,1) )
2784  call saturation_pres2qsat_all( qsat(:,1,1), temp(:,1,1), pres(:,1,1) )
2785 
2786  qv_sfc(1,1,1) = sfc_rh * 1.e-2_rp * qsat_sfc(1,1,1)
2787  do k = ks, ke
2788  if ( grid_cz(k) <= env_l1_ztop ) then ! Layer 1
2789  qv(k,1,1) = env_rh * 1.e-2_rp * qsat(k,1,1)
2790  elseif( grid_cz(k) <= env_l2_ztop ) then ! Layer 2
2791  qv(k,1,1) = env_rh * 1.e-2_rp * qsat(k,1,1)
2792  else ! Layer 3
2793  qv(k,1,1) = 0.0_rp
2794  endif
2795  enddo
2796 
2797  ! make density & pressure profile in moist condition
2798  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
2799  temp(:,1,1), & ! [OUT]
2800  pres(:,1,1), & ! [OUT]
2801  pott(:,1,1), & ! [IN]
2802  qv(:,1,1), & ! [IN]
2803  qc(:,1,1), & ! [IN]
2804  temp_sfc(1,1,1), & ! [OUT]
2805  pres_sfc(1,1,1), & ! [IN]
2806  pott_sfc(1,1,1), & ! [IN]
2807  qv_sfc(1,1,1), & ! [IN]
2808  qc_sfc(1,1,1) ) ! [IN]
2809 
2810  do j = js, je
2811  do i = is, ie
2812  do k = ks, ke
2813  dens(k,i,j) = dens(k,1,1)
2814  momz(k,i,j) = 0.0_rp
2815  momx(k,i,j) = env_u * dens(k,i,j)
2816  momy(k,i,j) = env_v * dens(k,i,j)
2817 
2818  ! make warm bubble
2819  rhot(k,i,j) = dens(k,1,1) * ( pott(k,1,1) + bbl_theta * bubble(k,i,j) )
2820 
2821  qtrc(k,i,j,i_qv) = qv(k,1,1)
2822  enddo
2823  enddo
2824  enddo
2825 
2826  call flux_setup
2827 
2828  return
2829  end subroutine mkinit_warmbubble
2830 
2831  !-----------------------------------------------------------------------------
2833  subroutine mkinit_supercell
2834  use scale_atmos_hydrometeor, only: &
2835  i_qv
2836  implicit none
2837 
2838  real(RP) :: RHO(KA)
2839  real(RP) :: VELX(KA)
2840  real(RP) :: VELY(KA)
2841  real(RP) :: POTT(KA)
2842  real(RP) :: QV(KA)
2843 
2844  ! Bubble
2845  real(RP) :: BBL_THETA = 3.d0 ! extremum of temperature in bubble [K]
2846 
2847  namelist / param_mkinit_supercell / &
2848  bbl_theta
2849 
2850  integer :: ierr
2851  integer :: k, i, j
2852  !---------------------------------------------------------------------------
2853 
2854  if( io_l ) write(io_fid_log,*)
2855  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit SUPERCELL] / Categ[preprocess] / Origin[SCALE-RM]'
2856 
2857  if ( i_qv < 1 ) then
2858  write(*,*) 'xxx QV is not registered'
2859  call prc_mpistop
2860  end if
2861 
2862  !--- read namelist
2863  rewind(io_fid_conf)
2864  read(io_fid_conf,nml=param_mkinit_supercell,iostat=ierr)
2865 
2866  if( ierr < 0 ) then !--- missing
2867  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
2868  elseif( ierr > 0 ) then !--- fatal error
2869  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_SUPERCELL. Check!'
2870  call prc_mpistop
2871  endif
2872  if( io_nml ) write(io_fid_nml,nml=param_mkinit_supercell)
2873 
2874  call read_sounding( rho, velx, vely, pott, qv ) ! (out)
2875 
2876  do j = js, je
2877  do i = is, ie
2878  do k = ks, ke
2879  dens(k,i,j) = rho(k)
2880  momz(k,i,j) = 0.0_rp
2881  momx(k,i,j) = rho(k) * velx(k)
2882  momy(k,i,j) = rho(k) * vely(k)
2883 
2884  ! make warm bubble
2885  rhot(k,i,j) = rho(k) * ( pott(k) + bbl_theta * bubble(k,i,j) )
2886 
2887  qtrc(k,i,j,i_qv) = qv(k)
2888  enddo
2889  enddo
2890  enddo
2891 
2892  call flux_setup
2893 
2894  return
2895  end subroutine mkinit_supercell
2896 
2897  !-----------------------------------------------------------------------------
2899  subroutine mkinit_squallline
2900  use scale_atmos_hydrometeor, only: &
2901  i_qv
2902  implicit none
2903 
2904  real(RP) :: RHO(KA)
2905  real(RP) :: VELX(KA)
2906  real(RP) :: VELY(KA)
2907  real(RP) :: POTT(KA)
2908  real(RP) :: QV(KA)
2909 
2910  real(RP) :: RANDOM_THETA = 0.01_rp
2911  real(RP) :: OFFSET_velx = 12.0_rp
2912  real(RP) :: OFFSET_vely = -2.0_rp
2913 
2914  namelist / param_mkinit_squallline / &
2915  random_theta, &
2916  offset_velx, &
2917  offset_vely
2918 
2919  integer :: ierr
2920  integer :: k, i, j
2921  !---------------------------------------------------------------------------
2922 
2923  if( io_l ) write(io_fid_log,*)
2924  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit SQUALLLINE] / Categ[preprocess] / Origin[SCALE-RM]'
2925 
2926  if ( i_qv < 1 ) then
2927  write(*,*) 'xxx QV is not registered'
2928  call prc_mpistop
2929  end if
2930 
2931  !--- read namelist
2932  rewind(io_fid_conf)
2933  read(io_fid_conf,nml=param_mkinit_squallline,iostat=ierr)
2934 
2935  if( ierr < 0 ) then !--- missing
2936  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
2937  elseif( ierr > 0 ) then !--- fatal error
2938  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_SQUALLLINE. Check!'
2939  call prc_mpistop
2940  endif
2941  if( io_nml ) write(io_fid_nml,nml=param_mkinit_squallline)
2942 
2943  call read_sounding( rho, velx, vely, pott, qv ) ! (out)
2944 
2945  call random_get(rndm) ! make random
2946  do j = js, je
2947  do i = is, ie
2948  do k = ks, ke
2949  dens(k,i,j) = rho(k)
2950  momz(k,i,j) = 0.0_rp
2951  momx(k,i,j) = ( velx(k) - offset_velx ) * rho(k)
2952  momy(k,i,j) = ( vely(k) - offset_vely ) * rho(k)
2953  rhot(k,i,j) = rho(k) * ( pott(k) + rndm(k,i,j) * random_theta )
2954 
2955  qtrc(k,i,j,i_qv) = qv(k)
2956  enddo
2957  enddo
2958  enddo
2959 
2960  call flux_setup
2961 
2962  call ocean_setup
2963 
2964  return
2965  end subroutine mkinit_squallline
2966 
2967  !-----------------------------------------------------------------------------
2969  subroutine mkinit_wk1982
2970  use scale_atmos_hydrometeor, only: &
2971  i_qv
2972  implicit none
2973 
2974  ! Surface state
2975  real(RP) :: SFC_THETA = 300.0_rp ! surface pot. temperature [K]
2976  real(RP) :: SFC_PRES ! surface pressure [Pa]
2977  ! Parameter in Weisman and Klemp (1982)
2978  real(RP) :: TR_Z = 12000.0_rp ! height of tropopause [m]
2979  real(RP) :: TR_THETA = 343.0_rp ! pot. temperature at tropopause [K]
2980  real(RP) :: TR_TEMP = 213.0_rp ! temperature at tropopause [K]
2981  real(RP) :: SHEAR_Z = 3000.0_rp ! center height of shear layer [m]
2982  real(RP) :: SHEAR_U = 15.0_rp ! velocity u over the shear layer [m/s]
2983  ! Bubble
2984  real(RP) :: BBL_THETA = 3.d0 ! extremum of temperature in bubble [K]
2985 
2986  namelist / param_mkinit_wk1982 / &
2987  sfc_theta, &
2988  sfc_pres, &
2989  tr_z, &
2990  tr_theta, &
2991  tr_temp, &
2992  shear_z, &
2993  shear_u, &
2994  bbl_theta
2995 
2996  real(RP) :: rh (KA,IA,JA)
2997  real(RP) :: rh_sfc(1 ,IA,JA)
2998 
2999  integer :: ierr
3000  integer :: k, i, j, iq
3001  !---------------------------------------------------------------------------
3002 
3003  if( io_l ) write(io_fid_log,*)
3004  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit WK1982] / Categ[preprocess] / Origin[SCALE-RM]'
3005 
3006  if ( i_qv < 1 ) then
3007  write(*,*) 'xxx QV is not registered'
3008  call prc_mpistop
3009  end if
3010 
3011  sfc_pres = pstd
3012 
3013  rewind(io_fid_conf)
3014  read(io_fid_conf,nml=param_mkinit_wk1982,iostat=ierr)
3015  if( ierr < 0 ) then !--- missing
3016  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
3017  elseif( ierr > 0 ) then !--- fatal error
3018  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_WK1982. Check!'
3019  call prc_mpistop
3020  endif
3021  if( io_nml ) write(io_fid_nml,nml=param_mkinit_wk1982)
3022 
3023  ! calc in dry condition
3024  do j = js, je
3025  do i = is, ie
3026  pres_sfc(1,i,j) = sfc_pres
3027  pott_sfc(1,i,j) = sfc_theta
3028  qv_sfc(1,i,j) = 0.0_rp
3029  qc_sfc(1,i,j) = 0.0_rp
3030 
3031  do k = ks, ke
3032  if ( real_cz(k,i,j) <= tr_z ) then ! below initial cloud top
3033  pott(k,i,j) = pott_sfc(1,i,j) &
3034  + ( tr_theta - pott_sfc(1,i,j) ) * ( real_cz(k,i,j) / tr_z )**1.25_rp
3035  else
3036  pott(k,i,j) = tr_theta * exp( grav * ( real_cz(k,i,j) - tr_z ) / cpdry / tr_temp )
3037  endif
3038 
3039  qv(k,i,j) = 0.0_rp
3040  qc(k,i,j) = 0.0_rp
3041  enddo
3042  enddo
3043  enddo
3044 
3045  ! make density & pressure profile in dry condition
3046  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
3047  temp(:,:,:), & ! [OUT]
3048  pres(:,:,:), & ! [OUT]
3049  pott(:,:,:), & ! [IN]
3050  qv(:,:,:), & ! [IN]
3051  qc(:,:,:), & ! [IN]
3052  temp_sfc(:,:,:), & ! [OUT]
3053  pres_sfc(:,:,:), & ! [IN]
3054  pott_sfc(:,:,:), & ! [IN]
3055  qv_sfc(:,:,:), & ! [IN]
3056  qc_sfc(:,:,:) ) ! [IN]
3057 
3058  ! calc QV from RH
3059  do j = js, je
3060  do i = is, ie
3061  rh_sfc(1,i,j) = 1.0_rp - 0.75_rp * ( real_fz(ks-1,i,j) / tr_z )**1.25_rp
3062 
3063  do k = ks, ke
3064  if ( real_cz(k,i,j) <= tr_z ) then ! below initial cloud top
3065  rh(k,i,j) = 1.0_rp - 0.75_rp * ( real_cz(k,i,j) / tr_z )**1.25_rp
3066  else
3067  rh(k,i,j) = 0.25_rp
3068  endif
3069  enddo
3070  enddo
3071  enddo
3072 
3073  call saturation_pres2qsat_all( qsat_sfc(1,:,:), temp_sfc(1,:,:), pres_sfc(1,:,:) )
3074  call saturation_pres2qsat_all( qsat(:,:,:), temp(:,:,:), pres(:,:,:) )
3075 
3076  do j = js, je
3077  do i = is, ie
3078  qv_sfc(1,i,j) = rh_sfc(1,i,j) * qsat_sfc(1,i,j)
3079  do k = ks, ke
3080  qv(k,i,j) = rh(k,i,j) * qsat(k,i,j)
3081  enddo
3082  enddo
3083  enddo
3084 
3085  ! make density & pressure profile in moist condition
3086  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
3087  temp(:,:,:), & ! [OUT]
3088  pres(:,:,:), & ! [OUT]
3089  pott(:,:,:), & ! [IN]
3090  qv(:,:,:), & ! [IN]
3091  qc(:,:,:), & ! [IN]
3092  temp_sfc(:,:,:), & ! [OUT]
3093  pres_sfc(:,:,:), & ! [IN]
3094  pott_sfc(:,:,:), & ! [IN]
3095  qv_sfc(:,:,:), & ! [IN]
3096  qc_sfc(:,:,:) ) ! [IN]
3097 
3098  do k = ks, ke
3099  if( io_l ) write(io_fid_log,*) k, real_cz(k,is,js), pres(k,is,js), pott(k,is,js), rh(k,is,js), qv(k,is,js)*1000
3100  enddo
3101 
3102  call comm_vars8( dens(:,:,:), 1 )
3103  call comm_wait ( dens(:,:,:), 1 )
3104 
3105  do j = js, je
3106  do i = is, ie
3107  do k = ks, ke
3108  momx(k,i,j) = shear_u * tanh( real_cz(k,i,j) / shear_z ) &
3109  * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
3110  enddo
3111  enddo
3112  enddo
3113 
3114  do j = js, je
3115  do i = is, ie
3116  do k = ks, ke
3117  momy(k,i,j) = 0.0_rp
3118  momz(k,i,j) = 0.0_rp
3119  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
3120 
3121  ! make warm bubble
3122  rhot(k,i,j) = dens(k,i,j) * ( pott(k,i,j) + bbl_theta * bubble(k,i,j) )
3123 
3124  qtrc(k,i,j,i_qv) = qv(k,i,j)
3125  enddo
3126  enddo
3127  enddo
3128 
3129 #ifndef DRY
3130  if ( qa >= 2 ) then
3131  do iq = 2, qa
3132  do j = js, je
3133  do i = is, ie
3134  do k = ks, ke
3135  qtrc(k,i,j,iq) = 0.0_rp
3136  enddo
3137  enddo
3138  enddo
3139  enddo
3140  endif
3141 #endif
3142 
3143  call flux_setup
3144 
3145  return
3146  end subroutine mkinit_wk1982
3147 
3148  !-----------------------------------------------------------------------------
3150  subroutine mkinit_dycoms2_rf01
3151  use scale_atmos_hydrometeor, only: &
3152  i_qv, &
3153  i_qc, &
3154  i_nc, &
3155  qhe
3156  use scale_atmos_phy_mp_suzuki10, only: &
3157  nccn
3158  use mod_atmos_admin, only: &
3159  atmos_phy_mp_type, &
3160  atmos_phy_ae_type
3161  implicit none
3162 
3163 #ifndef DRY
3164  real(RP) :: PERTURB_AMP = 0.0_rp
3165  integer :: RANDOM_LIMIT = 5
3166  integer :: RANDOM_FLAG = 0 ! 0 -> no perturbation
3167  ! 1 -> petrurbation for pt
3168  ! 2 -> perturbation for u, v, w
3169  logical :: USE_LWSET = .false. ! use liq. water. static energy temp.?
3170 
3171  namelist / param_mkinit_rf01 / &
3172  perturb_amp, &
3173  random_limit, &
3174  random_flag, &
3175  use_lwset
3176 
3177  real(RP) :: potl(KA,IA,JA) ! liquid potential temperature
3178  real(RP) :: LHV (KA,IA,JA) ! latent heat of vaporization [J/kg]
3179 
3180  real(RP) :: qall ! QV+QC
3181  real(RP) :: fact
3182  real(RP) :: pi2
3183  real(RP) :: sint
3184  real(RP) :: GEOP_sw ! switch for geopotential energy correction
3185 
3186  integer :: ierr
3187  integer :: k, i, j, iq
3188  !---------------------------------------------------------------------------
3189 
3190  pi2 = atan(1.0_rp) * 2.0_rp ! pi/2
3191 
3192  if( io_l ) write(io_fid_log,*)
3193  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit DYCOMS2RF01] / Categ[preprocess] / Origin[SCALE-RM]'
3194 
3195  rewind(io_fid_conf)
3196 
3197  if ( i_qv < 1 ) then
3198  write(*,*) 'xxx QV is not registered'
3199  call prc_mpistop
3200  end if
3201 
3202  read(io_fid_conf,nml=param_mkinit_rf01,iostat=ierr)
3203  if( ierr < 0 ) then !--- missing
3204  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
3205  elseif( ierr > 0 ) then !--- fatal error
3206  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_RF01. Check!'
3207  call prc_mpistop
3208  endif
3209  if( io_nml ) write(io_fid_nml,nml=param_mkinit_rf01)
3210 
3211  if ( use_lwset ) then
3212  geop_sw = 1.0_rp
3213  else
3214  geop_sw = 0.0_rp
3215  endif
3216 
3217  ! calc in dry condition
3218  do j = js, je
3219  do i = is, ie
3220 
3221  pres_sfc(1,i,j) = 1017.8e2_rp ! [Pa]
3222  pott_sfc(1,i,j) = 289.0_rp ! [K]
3223  qv_sfc(1,i,j) = 0.0_rp
3224  qc_sfc(1,i,j) = 0.0_rp
3225 
3226  do k = ks, ke
3227  velx(k,i,j) = 7.0_rp
3228  vely(k,i,j) = -5.5_rp
3229  if ( grid_cz(k) < 820.0_rp ) then ! below initial cloud top
3230  potl(k,i,j) = 289.0_rp - grav / cpdry * grid_cz(k) * geop_sw
3231  elseif( grid_cz(k) <= 860.0_rp ) then
3232  sint = sin( pi2 * ( grid_cz(k)-840.0_rp ) / 20.0_rp ) * 0.5_rp
3233  potl(k,i,j) = ( 289.0_rp - grav / cpdry * grid_cz(k) * geop_sw ) * (0.5_rp-sint) &
3234  + ( 297.5_rp+sign(abs(grid_cz(k)-840.0_rp)**(1.0_rp/3.0_rp),grid_cz(k)-840.0_rp) &
3235  - grav / cpdry * grid_cz(k) * geop_sw ) * (0.5_rp+sint)
3236  else
3237  potl(k,i,j) = 297.5_rp + ( grid_cz(k)-840.0_rp )**(1.0_rp/3.0_rp) &
3238  - grav / cpdry * grid_cz(k) * geop_sw
3239  endif
3240 
3241  qv(k,i,j) = 0.0_rp
3242  qc(k,i,j) = 0.0_rp
3243  enddo
3244 
3245  enddo
3246  enddo
3247 
3248  ! make density & pressure profile in dry condition
3249  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
3250  temp(:,:,:), & ! [OUT]
3251  pres(:,:,:), & ! [OUT]
3252  potl(:,:,:), & ! [IN]
3253  qv(:,:,:), & ! [IN]
3254  qc(:,:,:), & ! [IN]
3255  temp_sfc(:,:,:), & ! [OUT]
3256  pres_sfc(:,:,:), & ! [IN]
3257  pott_sfc(:,:,:), & ! [IN]
3258  qv_sfc(:,:,:), & ! [IN]
3259  qc_sfc(:,:,:) ) ! [IN]
3260 
3261  ! calc in moist condition
3262  do j = js, je
3263  do i = is, ie
3264  qv_sfc(1,i,j) = 9.0e-3_rp ! [kg/kg]
3265  qc_sfc(1,i,j) = 0.0_rp
3266 
3267  do k = ks, ke
3268  if ( grid_cz(k) < 820.0_rp ) then ! below initial cloud top
3269  qall = 9.0e-3_rp
3270  elseif( grid_cz(k) <= 860.0_rp ) then ! boundary
3271  sint = sin( pi2 * ( grid_cz(k)-840.0_rp ) / 20.0_rp ) * 0.5_rp
3272  qall = 9.0e-3_rp * (0.5_rp-sint) &
3273  + 1.5e-3_rp * (0.5_rp+sint)
3274  elseif( grid_cz(k) <= 5000.0_rp ) then
3275  qall = 1.5e-3_rp
3276  else
3277  qall = 0.0_rp
3278  endif
3279 
3280  if ( grid_cz(k) <= 600.0_rp ) then
3281  qc(k,i,j) = 0.0_rp
3282  elseif( grid_cz(k) < 820.0_rp ) then ! in the cloud
3283  fact = ( grid_cz(k)-600.0_rp ) / ( 840.0_rp-600.0_rp )
3284  qc(k,i,j) = 0.45e-3_rp * fact
3285  elseif( grid_cz(k) <= 860.0_rp ) then ! boundary
3286  sint = sin( pi2 * ( grid_cz(k)-840.0_rp ) / 20.0_rp ) * 0.5_rp
3287  fact = ( grid_cz(k)-600.0_rp ) / ( 840.0_rp-600.0_rp )
3288  qc(k,i,j) = 0.45e-3_rp * fact * (0.5_rp-sint)
3289  else
3290  qc(k,i,j) = 0.0_rp
3291  endif
3292 
3293  qv(k,i,j) = qall - qc(k,i,j)
3294  enddo
3295 
3296  enddo
3297  enddo
3298 
3299  call hydrometeor_lhv( lhv(:,:,:), temp(:,:,:) )
3300 
3301  do j = js, je
3302  do i = is, ie
3303  do k = ks, ke
3304  temp(k,i,j) = temp(k,i,j) + lhv(k,i,j) / cpdry * qc(k,i,j)
3305  enddo
3306  enddo
3307  enddo
3308 
3309  ! make density & pressure profile in moist condition
3310  call hydrostatic_buildrho_bytemp( dens(:,:,:), & ! [OUT]
3311  pott(:,:,:), & ! [OUT]
3312  pres(:,:,:), & ! [OUT]
3313  temp(:,:,:), & ! [IN]
3314  qv(:,:,:), & ! [IN]
3315  qc(:,:,:), & ! [IN]
3316  pott_sfc(:,:,:), & ! [OUT]
3317  pres_sfc(:,:,:), & ! [IN]
3318  temp_sfc(:,:,:), & ! [IN]
3319  qv_sfc(:,:,:), & ! [IN]
3320  qc_sfc(:,:,:) ) ! [IN]
3321 
3322  do j = js, je
3323  do i = is, ie
3324  dens( 1:ks-1,i,j) = dens(ks,i,j)
3325  dens(ke+1:ka, i,j) = dens(ke,i,j)
3326  enddo
3327  enddo
3328 
3329  call comm_vars8( dens(:,:,:), 1 )
3330  call comm_wait ( dens(:,:,:), 1 )
3331 
3332  call random_get(rndm) ! make random
3333  do j = js, je
3334  do i = is, ie
3335  do k = ks, ke
3336  if ( random_flag == 2 .and. k <= random_limit ) then ! below initial cloud top
3337  momz(k,i,j) = ( 2.0_rp * ( rndm(k,i,j)-0.5_rp ) * perturb_amp ) &
3338  * 0.5_rp * ( dens(k+1,i,j) + dens(k,i,j) )
3339  else
3340  momz(k,i,j) = 0.0_rp
3341  endif
3342  enddo
3343  enddo
3344  enddo
3345 
3346  call random_get(rndm) ! make random
3347  do j = js, je
3348  do i = is, ie
3349  do k = ks, ke
3350  if ( random_flag == 2 .AND. k <= random_limit ) then ! below initial cloud top
3351  momx(k,i,j) = ( velx(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.5_rp ) * perturb_amp ) &
3352  * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
3353  else
3354  momx(k,i,j) = velx(k,i,j) * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
3355  endif
3356  enddo
3357  enddo
3358  enddo
3359 
3360  call random_get(rndm) ! make random
3361  do j = js, je
3362  do i = is, ie
3363  do k = ks, ke
3364  if ( random_flag == 2 .AND. k <= random_limit ) then ! below initial cloud top
3365  momy(k,i,j) = ( vely(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.5_rp ) * perturb_amp ) &
3366  * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
3367  else
3368  momy(k,i,j) = vely(k,i,j) * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
3369  endif
3370  enddo
3371  enddo
3372  enddo
3373 
3374  call random_get(rndm) ! make random
3375  do j = js, je
3376  do i = is, ie
3377  do k = ks, ke
3378  if ( random_flag == 1 .and. k <= random_limit ) then ! below initial cloud top
3379  rhot(k,i,j) = ( pott(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.5_rp ) * perturb_amp ) &
3380  * dens(k,i,j)
3381  else
3382  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
3383  endif
3384  enddo
3385  enddo
3386  enddo
3387 
3388  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
3389  do j = js, je
3390  do i = is, ie
3391  do k = ks, ke
3392  qtrc(k,i,j,i_qv) = qv(k,i,j) + qc(k,i,j) !--- Super saturated air at initial
3393  do iq = qhe+1, qhe+nccn
3394  qtrc(k,i,j,iq) = qtrc(k,i,j,iq) / dens(k,i,j)
3395  enddo
3396  enddo
3397  enddo
3398  enddo
3399  else
3400  do j = js, je
3401  do i = is, ie
3402  do k = ks, ke
3403  qtrc(k,i,j,i_qv) = qv(k,i,j)
3404  qtrc(k,i,j,i_qc) = qc(k,i,j)
3405  enddo
3406  enddo
3407  enddo
3408 
3409  if ( i_nc > 0 ) then
3410  do j = js, je
3411  do i = is, ie
3412  do k = ks, ke
3413  if ( qc(k,i,j) > 0.0_rp ) then
3414  qtrc(k,i,j,i_nc) = 120.e6_rp / dens(k,i,j) ! [number/m3] / [kg/m3]
3415  endif
3416  enddo
3417  enddo
3418  enddo
3419  endif
3420  endif
3421 
3422 #endif
3423  if ( atmos_phy_ae_type == 'KAJINO13' ) then
3424  call aerosol_setup
3425  endif
3426  return
3427  end subroutine mkinit_dycoms2_rf01
3428 
3429  !-----------------------------------------------------------------------------
3431  subroutine mkinit_dycoms2_rf02
3432  use scale_atmos_hydrometeor, only: &
3433  i_qv, &
3434  i_qc, &
3435  i_nc, &
3436  qhe
3437  use scale_atmos_phy_mp_suzuki10, only: &
3438  nccn
3439  use mod_atmos_admin, only: &
3440  atmos_phy_mp_type, &
3441  atmos_phy_ae_type
3442  implicit none
3443 
3444 #ifndef DRY
3445  real(RP) :: PERTURB_AMP = 0.0_rp
3446  integer :: RANDOM_LIMIT = 5
3447  integer :: RANDOM_FLAG = 0 ! 0 -> no perturbation
3448  ! 1 -> perturbation for PT
3449  ! 2 -> perturbation for u,v,w
3450 
3451  namelist / param_mkinit_rf02 / &
3452  perturb_amp, &
3453  random_limit, &
3454  random_flag
3455 
3456  real(RP) :: potl(KA,IA,JA) ! liquid potential temperature
3457  real(RP) :: LHV (KA,IA,JA) ! latent heat of vaporization [J/kg]
3458 
3459  real(RP) :: qall ! QV+QC
3460  real(RP) :: fact
3461  real(RP) :: pi2
3462  real(RP) :: sint
3463 
3464  integer :: ierr
3465  integer :: k, i, j, iq
3466  !---------------------------------------------------------------------------
3467 
3468  pi2 = atan(1.0_rp) * 2.0_rp ! pi/2
3469  if( io_l ) write(io_fid_log,*)
3470  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit DYCOMS2RF02] / Categ[preprocess] / Origin[SCALE-RM]'
3471 
3472  if ( i_qv < 1 ) then
3473  write(*,*) 'xxx QV is not registered'
3474  call prc_mpistop
3475  end if
3476 
3477  rewind(io_fid_conf)
3478  read(io_fid_conf,nml=param_mkinit_rf02,iostat=ierr)
3479  if( ierr < 0 ) then !--- missing
3480  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
3481  elseif( ierr > 0 ) then !--- fatal error
3482  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_RF02. Check!'
3483  call prc_mpistop
3484  endif
3485  if( io_nml ) write(io_fid_nml,nml=param_mkinit_rf02)
3486 
3487  ! calc in dry condition
3488  call random_get(rndm) ! make random
3489  do j = js, je
3490  do i = is, ie
3491 
3492  pres_sfc(1,i,j) = 1017.8e2_rp ! [Pa]
3493  pott_sfc(1,i,j) = 288.3_rp ! [K]
3494  qv_sfc(1,i,j) = 0.0_rp
3495  qc_sfc(1,i,j) = 0.0_rp
3496 
3497  do k = ks, ke
3498  velx(k,i,j) = 3.0_rp + 4.3 * grid_cz(k)*1.e-3_rp
3499  vely(k,i,j) = -9.0_rp + 5.6 * grid_cz(k)*1.e-3_rp
3500 
3501  if ( grid_cz(k) < 775.0_rp ) then ! below initial cloud top
3502  potl(k,i,j) = 288.3_rp ! [K]
3503  else if ( grid_cz(k) <= 815.0_rp ) then
3504  sint = sin( pi2 * (grid_cz(k) - 795.0_rp)/20.0_rp )
3505  potl(k,i,j) = 288.3_rp * (1.0_rp-sint)*0.5_rp &
3506  + ( 295.0_rp+sign(abs(grid_cz(k)-795.0_rp)**(1.0_rp/3.0_rp),grid_cz(k)-795.0_rp) ) &
3507  * (1.0_rp+sint)*0.5_rp
3508  else
3509  potl(k,i,j) = 295.0_rp + ( grid_cz(k)-795.0_rp )**(1.0_rp/3.0_rp)
3510  endif
3511 
3512  qv(k,i,j) = 0.0_rp
3513  qc(k,i,j) = 0.0_rp
3514  enddo
3515  enddo
3516  enddo
3517 
3518  ! make density & pressure profile in dry condition
3519  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
3520  temp(:,:,:), & ! [OUT]
3521  pres(:,:,:), & ! [OUT]
3522  potl(:,:,:), & ! [IN]
3523  qv(:,:,:), & ! [IN]
3524  qc(:,:,:), & ! [IN]
3525  temp_sfc(:,:,:), & ! [OUT]
3526  pres_sfc(:,:,:), & ! [IN]
3527  pott_sfc(:,:,:), & ! [IN]
3528  qv_sfc(:,:,:), & ! [IN]
3529  qc_sfc(:,:,:) ) ! [IN]
3530 
3531  ! calc in moist condition
3532  do j = js, je
3533  do i = is, ie
3534  qv_sfc(1,i,j) = 9.45e-3_rp
3535  qc_sfc(1,i,j) = 0.0_rp
3536 
3537  do k = ks, ke
3538  if ( grid_cz(k) < 775.0_rp ) then ! below initial cloud top
3539  qall = 9.45e-3_rp ! [kg/kg]
3540  else if ( grid_cz(k) <= 815.0_rp ) then
3541  sint = sin( pi2 * (grid_cz(k) - 795.0_rp)/20.0_rp )
3542  qall = 9.45e-3_rp * (1.0_rp-sint)*0.5_rp + &
3543  ( 5.e-3_rp - 3.e-3_rp * ( 1.0_rp - exp( (795.0_rp-grid_cz(k))/500.0_rp ) ) ) * (1.0_rp+sint)*0.5_rp
3544  else
3545  qall = 5.e-3_rp - 3.e-3_rp * ( 1.0_rp - exp( (795.0_rp-grid_cz(k))/500.0_rp ) ) ! [kg/kg]
3546  endif
3547 
3548  if( grid_cz(k) < 400.0_rp ) then
3549  qc(k,i,j) = 0.0_rp
3550  elseif( grid_cz(k) < 775.0_rp ) then
3551  fact = ( grid_cz(k)-400.0_rp ) / ( 795.0_rp-400.0_rp )
3552  qc(k,i,j) = 0.65e-3_rp * fact
3553  elseif( grid_cz(k) <= 815.0_rp ) then
3554  sint = sin( pi2 * ( grid_cz(k)-795.0_rp )/20.0_rp )
3555  fact = ( grid_cz(k)-400.0_rp ) / ( 795.0_rp-400.0_rp )
3556  qc(k,i,j) = 0.65e-3_rp * fact * (1.0_rp-sint) * 0.5_rp
3557  else
3558  qc(k,i,j) = 0.0_rp
3559  endif
3560  qv(k,i,j) = qall - qc(k,i,j)
3561  enddo
3562 
3563  enddo
3564  enddo
3565 
3566  call hydrometeor_lhv( lhv(:,:,:), temp(:,:,:) )
3567 
3568  do j = js, je
3569  do i = is, ie
3570  do k = ks, ke
3571  temp(k,i,j) = temp(k,i,j) + lhv(k,i,j) / cpdry * qc(k,i,j)
3572  enddo
3573  enddo
3574  enddo
3575 
3576  ! make density & pressure profile in moist condition
3577  call hydrostatic_buildrho_bytemp( dens(:,:,:), & ! [OUT]
3578  pott(:,:,:), & ! [OUT]
3579  pres(:,:,:), & ! [OUT]
3580  temp(:,:,:), & ! [IN]
3581  qv(:,:,:), & ! [IN]
3582  qc(:,:,:), & ! [IN]
3583  pott_sfc(:,:,:), & ! [OUT]
3584  pres_sfc(:,:,:), & ! [IN]
3585  temp_sfc(:,:,:), & ! [IN]
3586  qv_sfc(:,:,:), & ! [IN]
3587  qc_sfc(:,:,:) ) ! [IN]
3588 
3589  do j = js, je
3590  do i = is, ie
3591  dens( 1:ks-1,i,j) = dens(ks,i,j)
3592  dens(ke+1:ka, i,j) = dens(ke,i,j)
3593  enddo
3594  enddo
3595 
3596  call comm_vars8( dens(:,:,:), 1 )
3597  call comm_wait ( dens(:,:,:), 1 )
3598 
3599  call random_get(rndm) ! make random
3600  do j = js, je
3601  do i = is, ie
3602  do k = ks, ke
3603  if( random_flag == 2 .and. k <= random_limit ) then
3604  momz(k,i,j) = ( 0.0_rp + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp ) &
3605  * 0.5_rp * ( dens(k+1,i,j) + dens(k,i,j) )
3606  else
3607  momz(k,i,j) = 0.0_rp
3608  endif
3609  enddo
3610  enddo
3611  enddo
3612 
3613  call random_get(rndm) ! make random
3614  do j = js, je
3615  do i = is, ie
3616  do k = ks, ke
3617  if( random_flag == 2 .and. k <= random_limit ) then
3618  momx(k,i,j) = ( velx(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp ) &
3619  * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
3620  else
3621  momx(k,i,j) = ( velx(k,i,j) ) * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
3622  endif
3623  enddo
3624  enddo
3625  enddo
3626 
3627  call random_get(rndm) ! make random
3628  do j = js, je
3629  do i = is, ie
3630  do k = ks, ke
3631  if( random_flag == 2 .and. k <= random_limit ) then
3632  momy(k,i,j) = ( vely(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp ) &
3633  * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
3634  else
3635  momy(k,i,j) = vely(k,i,j) * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
3636  endif
3637  enddo
3638  enddo
3639  enddo
3640 
3641  call random_get(rndm) ! make random
3642  do j = js, je
3643  do i = is, ie
3644  do k = ks, ke
3645  if( random_flag == 1 .and. k <= random_limit ) then
3646  rhot(k,i,j) = ( pott(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp ) &
3647  * dens(k,i,j)
3648  else
3649  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
3650  endif
3651  enddo
3652  enddo
3653  enddo
3654 
3655  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
3656  do j = js, je
3657  do i = is, ie
3658  do k = ks, ke
3659  !--- Super saturated air at initial
3660  qtrc(k,i,j,i_qv) = qv(k,i,j) + qc(k,i,j)
3661  do iq = qhe+1, qhe+nccn
3662  qtrc(k,i,j,iq) = qtrc(k,i,j,iq) / dens(k,i,j)
3663  enddo
3664  enddo
3665  enddo
3666  enddo
3667  else
3668  do j = js, je
3669  do i = is, ie
3670  do k = ks, ke
3671  qtrc(k,i,j,i_qv) = qv(k,i,j)
3672  qtrc(k,i,j,i_qc) = qc(k,i,j)
3673  enddo
3674  enddo
3675  enddo
3676 
3677  if ( i_nc > 0 ) then
3678  do j = js, je
3679  do i = is, ie
3680  do k = ks, ke
3681  if ( qc(k,i,j) > 0.0_rp ) then
3682  qtrc(k,i,j,i_nc) = 55.0e6_rp / dens(k,i,j) ! [number/m3] / [kg/m3]
3683  endif
3684  enddo
3685  enddo
3686  enddo
3687  endif
3688  endif
3689 
3690 #endif
3691  if ( atmos_phy_ae_type == 'KAJINO13' ) then
3692  call aerosol_setup
3693  endif
3694  return
3695  end subroutine mkinit_dycoms2_rf02
3696 
3697  !-----------------------------------------------------------------------------
3699  subroutine mkinit_dycoms2_rf02_dns
3700  use scale_atmos_hydrometeor, only: &
3701  i_qv, &
3702  i_qc, &
3703  i_nc, &
3704  qhe
3705  use scale_atmos_phy_mp_suzuki10, only: &
3706  nccn
3707  use mod_atmos_admin, only: &
3708  atmos_phy_mp_type
3709  implicit none
3710 
3711 #ifndef DRY
3712  real(RP) :: ZB = 750.0_rp ! domain bottom
3713 ! real(RP) :: ZT = 900.0_RP ! domain top
3714  real(RP) :: CONST_U = 0.0_rp
3715  real(RP) :: CONST_V = 0.0_rp
3716  real(RP) :: PRES_ZB = 93060.0_rp
3717  real(RP) :: PERTURB_AMP = 0.0_rp
3718  integer :: RANDOM_LIMIT = 5
3719  integer :: RANDOM_FLAG = 0 ! 0 -> no perturbation
3720  ! 1 -> perturbation for PT
3721  ! 2 -> perturbation for u,v,w
3722 
3723  namelist / param_mkinit_rf02_dns / &
3724  zb, const_u, const_v,pres_zb,&
3725  perturb_amp, &
3726  random_limit, &
3727  random_flag
3728 
3729  real(RP) :: potl(KA,IA,JA) ! liquid potential temperature
3730  real(RP) :: LHV (KA,IA,JA) ! latent heat of vaporization [J/kg]
3731 
3732  real(RP) :: qall ! QV+QC
3733  real(RP) :: fact
3734  real(RP) :: pi2
3735  real(RP) :: RovCP
3736 
3737  integer :: ierr
3738  integer :: k, i, j, iq
3739  !---------------------------------------------------------------------------
3740 
3741  pi2 = atan(1.0_rp) * 2.0_rp ! pi/2
3742 
3743  if( io_l ) write(io_fid_log,*)
3744  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit DYCOMS2RF02_DNS] / Categ[preprocess] / Origin[SCALE-RM]'
3745 
3746  if ( i_qv < 1 ) then
3747  write(*,*) 'xxx QV is not registered'
3748  call prc_mpistop
3749  end if
3750 
3751  rewind(io_fid_conf)
3752  read(io_fid_conf,nml=param_mkinit_rf02_dns,iostat=ierr)
3753  if( ierr < 0 ) then !--- missing
3754  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
3755  elseif( ierr > 0 ) then !--- fatal error
3756  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_RF02_DNS. Check!'
3757  call prc_mpistop
3758  endif
3759  if( io_nml ) write(io_fid_nml,nml=param_mkinit_rf02_dns)
3760 
3761  ! calc in dry condition
3762  call random_get(rndm) ! make random
3763  do j = js, je
3764  do i = is, ie
3765 
3766  pres_sfc(1,i,j) = pres_zb
3767 ! pott_sfc(1,i,j) = 288.3_RP ! [K]
3768 ! qv_sfc (1,i,j) = 9.45E-3_RP
3769 ! qc_sfc (1,i,j) = 0.0_RP
3770 
3771  do k = ks, ke
3772 
3773  velx(k,i,j) = const_u
3774  vely(k,i,j) = const_v
3775 
3776 ! if ( ZB+GRID_CZ(k) < 775.0_RP ) then ! below initial cloud top
3777  if ( zb+grid_cz(k) <= 795.0_rp ) then ! below initial cloud top
3778  potl(k,i,j) = 288.3_rp ! [K]
3779  qall = 9.45e-3_rp ! [kg/kg]
3780 ! necessary?
3781 ! else if ( GRID_CZ(k) <= 815.0_RP ) then
3782 ! sint = sin( pi2 * (GRID_CZ(k) - 795.0_RP)/20.0_RP )
3783 ! potl(k,i,j) = 288.3_RP * (1.0_RP-sint)*0.5_RP + &
3784 ! ( 295.0_RP+sign(abs(GRID_CZ(k)-795.0_RP)**(1.0_RP/3.0_RP),GRID_CZ(k)-795.0_RP) ) * (1.0_RP+sint)*0.5_RP
3785 ! qall = 9.45E-3_RP * (1.0_RP-sint)*0.5_RP + &
3786 ! ( 5.E-3_RP - 3.E-3_RP * ( 1.0_RP - exp( (795.0_RP-GRID_CZ(k))/500.0_RP ) ) ) * (1.0_RP+sint)*0.5_RP
3787  else
3788  potl(k,i,j) = 295.0_rp + ( zb+grid_cz(k)-795.0_rp )**(1.0_rp/3.0_rp)
3789  qall = 5.e-3_rp - 3.e-3_rp * ( 1.0_rp - exp( (795.0_rp-(zb+grid_cz(k)))/500.0_rp ) ) ! [kg/kg]
3790  endif
3791 
3792  if( zb+grid_cz(k) < 400.0_rp ) then
3793  qc(k,i,j) = 0.0_rp
3794  elseif( zb+grid_cz(k) <= 795.0_rp ) then
3795  fact = ( (zb+grid_cz(k))-400.0_rp ) / ( 795.0_rp-400.0_rp )
3796  qc(k,i,j) = 0.8e-3_rp * fact
3797  else
3798  qc(k,i,j) = 0.0_rp
3799  endif
3800  qv(k,i,j) = qall - qc(k,i,j)
3801 
3802 !if(i==is.and.j==js)write(*,*)'chkk',k,cz(k)+zb,qc(k,i,j),qv(k,i,j)
3803  enddo
3804  enddo
3805  enddo
3806 
3807 !write(*,*)'chk3',ks,ke
3808  ! extrapolation (temtative)
3809  pott_sfc(1,:,:) = potl(ks,:,:)-0.5*(potl(ks+1,:,:)-potl(ks,:,:))
3810  qv_sfc(1,:,:) = qv(ks,:,:)-0.5*(qv(ks+1,:,:)-qv(ks,:,:))
3811  qc_sfc(1,:,:) = qc(ks,:,:)-0.5*(qc(ks+1,:,:)-qc(ks,:,:))
3812 
3813  ! make density & pressure profile in moist condition
3814  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
3815  temp(:,:,:), & ! [OUT]
3816  pres(:,:,:), & ! [OUT]
3817  potl(:,:,:), & ! [IN]
3818  qv(:,:,:), & ! [IN]
3819  qc(:,:,:), & ! [IN]
3820  temp_sfc(:,:,:), & ! [OUT]
3821  pres_sfc(:,:,:), & ! [IN]
3822  pott_sfc(:,:,:), & ! [IN]
3823  qv_sfc(:,:,:), & ! [IN]
3824  qc_sfc(:,:,:) ) ! [IN]
3825 
3826  call hydrometeor_lhv( lhv(:,:,:), temp(:,:,:) )
3827 
3828  rovcp = rdry / cpdry
3829  do j = js, je
3830  do i = is, ie
3831  do k = ks, ke
3832  pott(k,i,j) = potl(k,i,j) + lhv(k,i,j) / cpdry * qc(k,i,j) * ( p00/pres(k,i,j) )**rovcp
3833  enddo
3834  enddo
3835  enddo
3836 
3837  ! make density & pressure profile in moist condition
3838  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
3839  temp(:,:,:), & ! [OUT]
3840  pres(:,:,:), & ! [OUT]
3841  pott(:,:,:), & ! [IN]
3842  qv(:,:,:), & ! [IN]
3843  qc(:,:,:), & ! [IN]
3844  temp_sfc(:,:,:), & ! [OUT]
3845  pres_sfc(:,:,:), & ! [IN]
3846  pott_sfc(:,:,:), & ! [IN]
3847  qv_sfc(:,:,:), & ! [IN]
3848  qc_sfc(:,:,:) ) ! [IN]
3849 
3850  do j = js, je
3851  do i = is, ie
3852  dens( 1:ks-1,i,j) = dens(ks,i,j)
3853  dens(ke+1:ka, i,j) = dens(ke,i,j)
3854  enddo
3855  enddo
3856 
3857  call comm_vars8( dens(:,:,:), 1 )
3858  call comm_wait ( dens(:,:,:), 1 )
3859 
3860  call random_get(rndm) ! make random
3861  do j = js, je
3862  do i = is, ie
3863  do k = ks, ke
3864  if( random_flag == 2 .and. k <= random_limit ) then
3865  momz(k,i,j) = ( 0.0_rp + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp ) &
3866  * 0.5_rp * ( dens(k+1,i,j) + dens(k,i,j) )
3867  else
3868  momz(k,i,j) = 0.0_rp
3869  endif
3870  enddo
3871  enddo
3872  enddo
3873 
3874 !write(*,*)'chk8'
3875  call random_get(rndm) ! make random
3876  do j = js, je
3877  do i = is, ie
3878  do k = ks, ke
3879  if( random_flag == 2 .and. k <= random_limit ) then
3880  momx(k,i,j) = ( velx(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp ) &
3881  * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
3882  else
3883  momx(k,i,j) = ( velx(k,i,j) ) * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
3884  endif
3885  enddo
3886  enddo
3887  enddo
3888 !write(*,*)'chk9'
3889 
3890  call random_get(rndm) ! make random
3891  do j = js, je
3892  do i = is, ie
3893  do k = ks, ke
3894  if( random_flag == 2 .and. k <= random_limit ) then
3895  momy(k,i,j) = ( vely(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp ) &
3896  * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
3897  else
3898  momy(k,i,j) = vely(k,i,j) * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
3899  endif
3900  enddo
3901  enddo
3902  enddo
3903 !write(*,*)'chk10'
3904 
3905  call random_get(rndm) ! make random
3906  do j = js, je
3907  do i = is, ie
3908  do k = ks, ke
3909  if( random_flag == 1 .and. k <= random_limit ) then
3910  rhot(k,i,j) = ( pott(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp ) &
3911  * dens(k,i,j)
3912  else
3913  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
3914  endif
3915  enddo
3916  enddo
3917  enddo
3918 
3919 !write(*,*)'chk11'
3920  do iq = 1, qa
3921  do j = js, je
3922  do i = is, ie
3923  do k = ks, ke
3924  qtrc(k,i,j,iq) = 0.0_rp
3925  enddo
3926  enddo
3927  enddo
3928  enddo
3929 
3930 !write(*,*)'chk12'
3931  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
3932  do j = js, je
3933  do i = is, ie
3934  do k = ks, ke
3935  !--- Super saturated air at initial
3936  qtrc(k,i,j,i_qv) = qv(k,i,j) + qc(k,i,j)
3937 
3938  !--- for aerosol
3939  do iq = 1, nccn
3940  qtrc(k,i,j,qhe+iq) = gan(iq) / dens(k,i,j)
3941  enddo
3942  enddo
3943  enddo
3944  enddo
3945  else
3946  do j = js, je
3947  do i = is, ie
3948  do k = ks, ke
3949  qtrc(k,i,j,i_qv) = qv(k,i,j)
3950  qtrc(k,i,j,i_qc) = qc(k,i,j)
3951  enddo
3952  enddo
3953  enddo
3954 
3955  if ( i_nc > 0 ) then
3956  do j = js, je
3957  do i = is, ie
3958  do k = ks, ke
3959  if ( qc(k,i,j) > 0.0_rp ) then
3960  qtrc(k,i,j,i_nc) = 55.0e6_rp / dens(k,i,j) ! [number/m3] / [kg/m3]
3961  endif
3962  enddo
3963  enddo
3964  enddo
3965  endif
3966  endif
3967 
3968 #endif
3969  return
3970  end subroutine mkinit_dycoms2_rf02_dns
3971 
3972  !-----------------------------------------------------------------------------
3974  subroutine mkinit_rico
3975  use scale_atmos_hydrometeor, only: &
3976  i_qv, &
3977  i_qc, &
3978  i_nc, &
3979  qhe
3980  use scale_atmos_phy_mp_suzuki10, only: &
3981  nccn
3982  use mod_atmos_admin, only: &
3983  atmos_phy_mp_type, &
3984  atmos_phy_ae_type
3985  implicit none
3986 
3987 #ifndef DRY
3988  real(RP):: PERTURB_AMP_PT = 0.1_rp
3989  real(RP):: PERTURB_AMP_QV = 2.5e-5_rp
3990 
3991  namelist / param_mkinit_rico / &
3992  perturb_amp_pt, &
3993  perturb_amp_qv
3994 
3995  real(RP) :: LHV (KA,IA,JA) ! latent heat of vaporization [J/kg]
3996  real(RP) :: potl(KA,IA,JA) ! liquid potential temperature
3997  real(RP) :: qall ! QV+QC
3998  real(RP) :: fact
3999 
4000  integer :: ierr
4001  integer :: k, i, j, iq
4002  !---------------------------------------------------------------------------
4003 
4004  if( io_l ) write(io_fid_log,*)
4005  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit RICO] / Categ[preprocess] / Origin[SCALE-RM]'
4006 
4007  if ( i_qv < 1 ) then
4008  write(*,*) 'xxx QV is not registered'
4009  call prc_mpistop
4010  end if
4011 
4012  rewind(io_fid_conf)
4013  read(io_fid_conf,nml=param_mkinit_rico,iostat=ierr)
4014  if( ierr < 0 ) then !--- missing
4015  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
4016  elseif( ierr > 0 ) then !--- fatal error
4017  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_RICO. Check!'
4018  call prc_mpistop
4019  endif
4020  if( io_nml ) write(io_fid_nml,nml=param_mkinit_rico)
4021 
4022  ! calc in moist condition
4023  do j = js, je
4024  do i = is, ie
4025 
4026  pres_sfc(1,i,j) = 1015.4e2_rp ! [Pa]
4027  pott_sfc(1,i,j) = 297.9_rp
4028  qv_sfc(1,i,j) = 0.0_rp
4029  qc_sfc(1,i,j) = 0.0_rp
4030 
4031  do k = ks, ke
4032  !--- potential temperature
4033  if ( grid_cz(k) < 740.0_rp ) then ! below initial cloud top
4034  potl(k,i,j) = 297.9_rp
4035  else
4036  fact = ( grid_cz(k)-740.0_rp ) * ( 317.0_rp-297.9_rp ) / ( 4000.0_rp-740.0_rp )
4037  potl(k,i,j) = 297.9_rp + fact
4038  endif
4039 
4040  !--- horizontal wind velocity
4041  if ( grid_cz(k) <= 4000.0_rp ) then ! below initial cloud top
4042  fact = ( grid_cz(k)-0.0_rp ) * ( -1.9_rp+9.9_rp ) / ( 4000.0_rp-0.0_rp )
4043  velx(k,i,j) = -9.9_rp + fact
4044  vely(k,i,j) = -3.8_rp
4045  else
4046  velx(k,i,j) = -1.9_rp
4047  vely(k,i,j) = -3.8_rp
4048  endif
4049 
4050  qv(k,i,j) = 0.0_rp
4051  qc(k,i,j) = 0.0_rp
4052 
4053  enddo
4054 
4055  enddo
4056  enddo
4057 
4058  ! make density & pressure profile in moist condition
4059  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
4060  temp(:,:,:), & ! [OUT]
4061  pres(:,:,:), & ! [OUT]
4062  potl(:,:,:), & ! [IN]
4063  qv(:,:,:), & ! [IN]
4064  qc(:,:,:), & ! [IN]
4065  temp_sfc(:,:,:), & ! [OUT]
4066  pres_sfc(:,:,:), & ! [IN]
4067  pott_sfc(:,:,:), & ! [IN]
4068  qv_sfc(:,:,:), & ! [IN]
4069  qc_sfc(:,:,:) ) ! [IN]
4070 
4071 
4072  do j = js, je
4073  do i = is, ie
4074  qv_sfc(1,i,j) = 16.0e-3_rp ! [kg/kg]
4075  qc_sfc(1,i,j) = 0.0_rp
4076 
4077  do k = ks, ke
4078  !--- mixing ratio of vapor
4079  if ( grid_cz(k) <= 740.0_rp ) then ! below initial cloud top
4080  fact = ( grid_cz(k)-0.0_rp ) * ( 13.8e-3_rp-16.0e-3_rp ) / ( 740.0_rp-0.0_rp )
4081  qall = 16.0e-3_rp + fact
4082  elseif ( grid_cz(k) <= 3260.0_rp ) then ! boundary
4083  fact = ( grid_cz(k)-740.0_rp ) * ( 2.4e-3_rp-13.8e-3_rp ) / ( 3260.0_rp-740.0_rp )
4084  qall = 13.8e-3_rp + fact
4085  elseif( grid_cz(k) <= 4000.0_rp ) then
4086  fact = ( grid_cz(k)-3260.0_rp ) * ( 1.8e-3_rp-2.4e-3_rp ) / ( 4000.0_rp-3260.0_rp )
4087  qall = 2.4e-3_rp + fact
4088  else
4089  qall = 0.0_rp
4090  endif
4091 
4092  qc(k,i,j) = 0.0_rp
4093  qv(k,i,j) = qall - qc(k,i,j)
4094  enddo
4095 
4096  enddo
4097  enddo
4098 
4099  call hydrometeor_lhv( lhv(:,:,:), temp(:,:,:) )
4100 
4101  do j = js, je
4102  do i = is, ie
4103  do k = ks, ke
4104  temp(k,i,j) = temp(k,i,j) + lhv(k,i,j) / cpdry * qc(k,i,j)
4105  enddo
4106  enddo
4107  enddo
4108 
4109  ! make density & pressure profile in moist condition
4110  call hydrostatic_buildrho_bytemp( dens(:,:,:), & ! [OUT]
4111  pott(:,:,:), & ! [OUT]
4112  pres(:,:,:), & ! [OUT]
4113  temp(:,:,:), & ! [IN]
4114  qv(:,:,:), & ! [IN]
4115  qc(:,:,:), & ! [IN]
4116  pott_sfc(:,:,:), & ! [OUT]
4117  pres_sfc(:,:,:), & ! [IN]
4118  temp_sfc(:,:,:), & ! [IN]
4119  qv_sfc(:,:,:), & ! [IN]
4120  qc_sfc(:,:,:) ) ! [IN]
4121 
4122 
4123  do j = js, je
4124  do i = is, ie
4125  dens( 1:ks-1,i,j) = dens(ks,i,j)
4126  dens(ke+1:ka ,i,j) = dens(ke,i,j)
4127  enddo
4128  enddo
4129 
4130  call comm_vars8( dens(:,:,:), 1 )
4131  call comm_wait ( dens(:,:,:), 1 )
4132 
4133  do j = js, je
4134  do i = is, ie
4135  do k = ks, ke
4136  momz(k,i,j) = 0.0_rp
4137  enddo
4138  enddo
4139  enddo
4140 
4141  do j = js, je
4142  do i = is, ie
4143  do k = ks, ke
4144  momx(k,i,j) = velx(k,i,j) * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
4145  enddo
4146  enddo
4147  enddo
4148 
4149  do j = js, je
4150  do i = is, ie
4151  do k = ks, ke
4152  momy(k,i,j) = vely(k,i,j) * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
4153  enddo
4154  enddo
4155  enddo
4156 
4157  call random_get(rndm) ! make random
4158  do j = js, je
4159  do i = is, ie
4160  do k = ks, ke
4161  rhot(k,i,j) = ( pott(k,i,j)+2.0_rp*( rndm(k,i,j)-0.5_rp )*perturb_amp_pt ) * dens(k,i,j)
4162  enddo
4163  enddo
4164  enddo
4165 
4166  call random_get(rndm) ! make random
4167 
4168  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
4169  do j = js, je
4170  do i = is, ie
4171  do k = ks, ke
4172  !--- Super saturated air at initial
4173  qtrc(k,i,j,i_qv) = qv(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp_qv &
4174  + qc(k,i,j)
4175  do iq = qhe+1, qhe+nccn
4176  qtrc(k,i,j,iq) = qtrc(k,i,j,iq) / dens(k,i,j)
4177  enddo
4178  enddo
4179  enddo
4180  enddo
4181  else
4182  do j = js, je
4183  do i = is, ie
4184  do k = ks, ke
4185  qtrc(k,i,j,i_qv) = qv(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp_qv
4186  qtrc(k,i,j,i_qc) = qc(k,i,j)
4187  enddo
4188  enddo
4189  enddo
4190 
4191  if ( i_nc > 0 ) then
4192  do j = js, je
4193  do i = is, ie
4194  do k = ks, ke
4195  if ( qc(k,i,j) > 0.0_rp ) then
4196  qtrc(k,i,j,i_nc) = 70.e6_rp / dens(k,i,j) ! [number/m3] / [kg/m3]
4197  endif
4198  enddo
4199  enddo
4200  enddo
4201  endif
4202  endif
4203 
4204 #endif
4205  if ( atmos_phy_ae_type == 'KAJINO13' ) then
4206  call aerosol_setup
4207  endif
4208  return
4209  end subroutine mkinit_rico
4210 
4211  !-----------------------------------------------------------------------------
4213  subroutine mkinit_bomex
4214  use scale_atmos_hydrometeor, only: &
4215  i_qv, &
4216  i_qc, &
4217  i_nc, &
4218  qhe
4219  use scale_atmos_phy_mp_suzuki10, only: &
4220  nccn
4221  use mod_atmos_admin, only: &
4222  atmos_phy_mp_type, &
4223  atmos_phy_ae_type
4224  implicit none
4225 
4226 #ifndef DRY
4227  real(RP):: PERTURB_AMP_PT = 0.1_rp
4228  real(RP):: PERTURB_AMP_QV = 2.5e-5_rp
4229 
4230  namelist / param_mkinit_bomex / &
4231  perturb_amp_pt, &
4232  perturb_amp_qv
4233 
4234  real(RP) :: LHV (KA,IA,JA) ! latent heat of vaporization [J/kg]
4235  real(RP) :: potl(KA,IA,JA) ! liquid potential temperature
4236  real(RP) :: qall ! QV+QC
4237  real(RP) :: fact
4238 
4239  integer :: ierr
4240  integer :: k, i, j, iq
4241  !---------------------------------------------------------------------------
4242 
4243  if( io_l ) write(io_fid_log,*)
4244  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit BOMEX] / Categ[preprocess] / Origin[SCALE-RM]'
4245 
4246  if ( i_qv < 1 ) then
4247  write(*,*) 'xxx QV is not registered'
4248  call prc_mpistop
4249  end if
4250 
4251  rewind(io_fid_conf)
4252  read(io_fid_conf,nml=param_mkinit_bomex,iostat=ierr)
4253  if( ierr < 0 ) then !--- missing
4254  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
4255  elseif( ierr > 0 ) then !--- fatal error
4256  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_BOMEX. Check!'
4257  call prc_mpistop
4258  endif
4259  if( io_nml ) write(io_fid_nml,nml=param_mkinit_bomex)
4260 
4261  ! calc in moist condition
4262  do j = js, je
4263  do i = is, ie
4264 
4265  pres_sfc(1,i,j) = 1015.e2_rp ! [Pa]
4266  pott_sfc(1,i,j) = 299.1_rp
4267  qv_sfc(1,i,j) = 0.0_rp
4268  qc_sfc(1,i,j) = 0.0_rp
4269 
4270  do k = ks, ke
4271  !--- potential temperature
4272  if ( grid_cz(k) < 520.0_rp ) then ! below initial cloud top
4273  potl(k,i,j) = 298.7_rp
4274  elseif( grid_cz(k) < 1480.0_rp ) then
4275  fact = ( grid_cz(k)-520.0_rp ) * ( 302.4_rp-298.7_rp ) / ( 1480.0_rp-520.0_rp )
4276  potl(k,i,j) = 298.7_rp + fact
4277  elseif( grid_cz(k) < 2000.0_rp ) then
4278  fact = ( grid_cz(k)-1480.0_rp ) * ( 308.2_rp-302.4_rp ) / ( 2000.0_rp-1480.0_rp )
4279  potl(k,i,j) = 302.4_rp + fact
4280  else
4281  fact = ( grid_cz(k)-2000.0_rp ) * 3.65e-3_rp
4282  potl(k,i,j) = 308.2_rp + fact
4283  endif
4284 
4285  !--- horizontal wind velocity
4286  if ( grid_cz(k) <= 700.0_rp ) then ! below initial cloud top
4287  velx(k,i,j) = -8.75_rp
4288  vely(k,i,j) = 0.0_rp
4289  else
4290  fact = 1.8e-3_rp * ( grid_cz(k)-700.0_rp )
4291  velx(k,i,j) = -8.75_rp + fact
4292  vely(k,i,j) = 0.0_rp
4293  endif
4294 
4295  qv(k,i,j) = 0.0_rp
4296  qc(k,i,j) = 0.0_rp
4297 
4298  enddo
4299 
4300  enddo
4301  enddo
4302 
4303  ! make density & pressure profile in moist condition
4304  call hydrostatic_buildrho( dens(:,:,:), & ! [OUT]
4305  temp(:,:,:), & ! [OUT]
4306  pres(:,:,:), & ! [OUT]
4307  potl(:,:,:), & ! [IN]
4308  qv(:,:,:), & ! [IN]
4309  qc(:,:,:), & ! [IN]
4310  temp_sfc(:,:,:), & ! [OUT]
4311  pres_sfc(:,:,:), & ! [IN]
4312  pott_sfc(:,:,:), & ! [IN]
4313  qv_sfc(:,:,:), & ! [IN]
4314  qc_sfc(:,:,:) ) ! [IN]
4315 
4316 
4317  do j = js, je
4318  do i = is, ie
4319  qv_sfc(1,i,j) = 22.45e-3_rp ! [kg/kg]
4320  qc_sfc(1,i,j) = 0.0_rp
4321 
4322  do k = ks, ke
4323  !--- mixing ratio of vapor
4324  if ( grid_cz(k) <= 520.0_rp ) then ! below initial cloud top
4325  fact = ( grid_cz(k)-0.0_rp ) * ( 16.3e-3_rp-17.0e-3_rp ) / ( 520.0_rp-0.0_rp )
4326  qall = 17.0e-3_rp + fact
4327  elseif ( grid_cz(k) <= 1480.0_rp ) then ! boundary
4328  fact = ( grid_cz(k)-520.0_rp ) * ( 10.7e-3_rp-16.3e-3_rp ) / ( 1480.0_rp-520.0_rp )
4329  qall = 16.3e-3_rp + fact
4330  elseif( grid_cz(k) <= 2000.0_rp ) then
4331  fact = ( grid_cz(k)-1480.0_rp ) * ( 4.2e-3_rp-10.7e-3_rp ) / ( 2000.0_rp-1480.0_rp )
4332  qall = 10.7e-3_rp + fact
4333  else
4334  fact = ( grid_cz(k)-2000.0_rp ) * ( -1.2e-6_rp )
4335  qall = 4.2e-3_rp + fact
4336  endif
4337 
4338  qc(k,i,j) = 0.0_rp
4339  qv(k,i,j) = qall - qc(k,i,j)
4340  enddo
4341 
4342  enddo
4343  enddo
4344 
4345  call hydrometeor_lhv( lhv(:,:,:), temp(:,:,:) )
4346 
4347  do j = js, je
4348  do i = is, ie
4349  do k = ks, ke
4350  temp(k,i,j) = temp(k,i,j) + lhv(k,i,j) / cpdry * qc(k,i,j)
4351  enddo
4352  enddo
4353  enddo
4354 
4355  ! make density & pressure profile in moist condition
4356  call hydrostatic_buildrho_bytemp( dens(:,:,:), & ! [OUT]
4357  pott(:,:,:), & ! [OUT]
4358  pres(:,:,:), & ! [OUT]
4359  temp(:,:,:), & ! [IN]
4360  qv(:,:,:), & ! [IN]
4361  qc(:,:,:), & ! [IN]
4362  pott_sfc(:,:,:), & ! [OUT]
4363  pres_sfc(:,:,:), & ! [IN]
4364  temp_sfc(:,:,:), & ! [IN]
4365  qv_sfc(:,:,:), & ! [IN]
4366  qc_sfc(:,:,:) ) ! [IN]
4367 
4368 
4369  do j = js, je
4370  do i = is, ie
4371  dens( 1:ks-1,i,j) = dens(ks,i,j)
4372  dens(ke+1:ka ,i,j) = dens(ke,i,j)
4373  enddo
4374  enddo
4375 
4376  call comm_vars8( dens(:,:,:), 1 )
4377  call comm_wait ( dens(:,:,:), 1 )
4378 
4379  do j = js, je
4380  do i = is, ie
4381  do k = ks, ke
4382  momz(k,i,j) = 0.0_rp
4383  enddo
4384  enddo
4385  enddo
4386 
4387  do j = js, je
4388  do i = is, ie
4389  do k = ks, ke
4390  momx(k,i,j) = velx(k,i,j) * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
4391  enddo
4392  enddo
4393  enddo
4394 
4395  do j = js, je
4396  do i = is, ie
4397  do k = ks, ke
4398  momy(k,i,j) = vely(k,i,j) * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
4399  enddo
4400  enddo
4401  enddo
4402 
4403  call random_get(rndm) ! make random
4404  do j = js, je
4405  do i = is, ie
4406  do k = ks, ke
4407  if( grid_cz(k) <= 1600.0_rp ) then !--- lowest 40 model layer when dz=40m
4408  rhot(k,i,j) = ( pott(k,i,j)+2.0_rp*( rndm(k,i,j)-0.5_rp )*perturb_amp_pt ) * dens(k,i,j)
4409  else
4410  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
4411  endif
4412  enddo
4413  enddo
4414  enddo
4415 
4416  call random_get(rndm) ! make random
4417 
4418  if ( atmos_phy_mp_type == 'SUZUKI10' ) then
4419  do j = js, je
4420  do i = is, ie
4421  do k = ks, ke
4422  !--- Super saturated air at initial
4423  if( grid_cz(k) <= 1600.0_rp ) then !--- lowest 40 model layer when dz=40m
4424  qtrc(k,i,j,i_qv) = qv(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp_qv &
4425  + qc(k,i,j)
4426  else
4427  qtrc(k,i,j,i_qv) = qv(k,i,j) + qc(k,i,j)
4428  endif
4429  do iq = qhe+1, qhe+nccn
4430  qtrc(k,i,j,iq) = qtrc(k,i,j,iq) / dens(k,i,j)
4431  enddo
4432  enddo
4433  enddo
4434  enddo
4435  else
4436  do j = js, je
4437  do i = is, ie
4438  do k = ks, ke
4439  if( grid_cz(k) <= 1600.0_rp ) then !--- lowest 40 model layer when dz=40m
4440  qtrc(k,i,j,i_qv) = qv(k,i,j) + 2.0_rp * ( rndm(k,i,j)-0.50_rp ) * perturb_amp_qv
4441  else
4442  qtrc(k,i,j,i_qv) = qv(k,i,j)
4443  endif
4444  qtrc(k,i,j,i_qc) = qc(k,i,j)
4445  enddo
4446  enddo
4447  enddo
4448 
4449  if ( i_nc > 0 ) then
4450  do j = js, je
4451  do i = is, ie
4452  do k = ks, ke
4453  if ( qc(k,i,j) > 0.0_rp ) then
4454  qtrc(k,i,j,i_nc) = 70.e6_rp / dens(k,i,j) ! [number/m3] / [kg/m3]
4455  endif
4456  enddo
4457  enddo
4458  enddo
4459  endif
4460  endif
4461 
4462 #endif
4463  if ( atmos_phy_ae_type == 'KAJINO13' ) then
4464  call aerosol_setup
4465  endif
4466  return
4467  end subroutine mkinit_bomex
4468 
4469  !-----------------------------------------------------------------------------
4470  subroutine mkinit_interporation
4471  use gtool_file, only: &
4472  filegetshape, &
4473  fileread
4474  use scale_atmos_hydrometeor, only: &
4475  i_qv, &
4476  i_qc
4477  implicit none
4478 
4479  real(RP) :: dz(KA,IA,JA)
4480 
4481  real(RP) :: W(KA,IA,JA)
4482  real(RP) :: U(KA,IA,JA)
4483  real(RP) :: V(KA,IA,JA)
4484 
4485  real(RP) :: fact_cz0(KA)
4486  real(RP) :: fact_cz1(KA)
4487  real(RP) :: fact_fz0(KA)
4488  real(RP) :: fact_fz1(KA)
4489  real(RP) :: fact_cx0(IA)
4490  real(RP) :: fact_cx1(IA)
4491  real(RP) :: fact_fx0(IA)
4492  real(RP) :: fact_fx1(IA)
4493  real(RP) :: fact_cy0(JA)
4494  real(RP) :: fact_cy1(JA)
4495  real(RP) :: fact_fy0(JA)
4496  real(RP) :: fact_fy1(JA)
4497 
4498  integer :: idx_cz0(KA)
4499  integer :: idx_cz1(KA)
4500  integer :: idx_fz0(KA)
4501  integer :: idx_fz1(KA)
4502  integer :: idx_cx0(IA)
4503  integer :: idx_cx1(IA)
4504  integer :: idx_fx0(IA)
4505  integer :: idx_fx1(IA)
4506  integer :: idx_cy0(JA)
4507  integer :: idx_cy1(JA)
4508  integer :: idx_fy0(JA)
4509  integer :: idx_fy1(JA)
4510 
4511  real(RP), allocatable :: DENS_ORG(:,:,:)
4512  real(RP), allocatable :: MOMZ_ORG(:,:,:)
4513  real(RP), allocatable :: MOMX_ORG(:,:,:)
4514  real(RP), allocatable :: MOMY_ORG(:,:,:)
4515  real(RP), allocatable :: RHOT_ORG(:,:,:)
4516  real(RP), allocatable :: QTRC_ORG(:,:,:,:)
4517 
4518  real(RP), allocatable :: W_ORG(:,:,:)
4519  real(RP), allocatable :: U_ORG(:,:,:)
4520  real(RP), allocatable :: V_ORG(:,:,:)
4521  real(RP), allocatable :: POTT_ORG(:,:,:)
4522 
4523  real(RP), allocatable :: CZ_ORG(:)
4524  real(RP), allocatable :: FZ_ORG(:)
4525  real(RP), allocatable :: CX_ORG(:)
4526  real(RP), allocatable :: FX_ORG(:)
4527  real(RP), allocatable :: CY_ORG(:)
4528  real(RP), allocatable :: FY_ORG(:)
4529 
4530  integer :: dims(3)
4531 
4532  character(len=H_LONG) :: BASENAME_ORG = ''
4533 
4534  namelist / param_mkinit_interporation / &
4535  basename_org
4536 
4537  integer :: ierr
4538  integer :: k, i, j, iq
4539  !---------------------------------------------------------------------------
4540 
4541  if( io_l ) write(io_fid_log,*)
4542  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit INTERPORATION] / Categ[preprocess] / Origin[SCALE-RM]'
4543 
4544  !--- read namelist
4545  rewind(io_fid_conf)
4546  read(io_fid_conf,nml=param_mkinit_interporation,iostat=ierr)
4547 
4548  if( ierr < 0 ) then !--- missing
4549  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
4550  elseif( ierr > 0 ) then !--- fatal error
4551  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_INTERPORATION. Check!'
4552  call prc_mpistop
4553  endif
4554  if( io_nml ) write(io_fid_nml,nml=param_mkinit_interporation)
4555 
4556  call filegetshape( dims(:), &
4557  basename_org, "DENS", 1, single=.true. )
4558 
4559  allocate( dens_org(dims(1),dims(2),dims(3)) )
4560  allocate( momz_org(dims(1),dims(2),dims(3)) )
4561  allocate( momx_org(dims(1),dims(2),dims(3)) )
4562  allocate( momy_org(dims(1),dims(2),dims(3)) )
4563  allocate( rhot_org(dims(1),dims(2),dims(3)) )
4564  allocate( qtrc_org(dims(1),dims(2),dims(3),qa) )
4565 
4566  allocate( w_org(dims(1),dims(2),dims(3)) )
4567  allocate( u_org(dims(1),dims(2),dims(3)) )
4568  allocate( v_org(dims(1),dims(2),dims(3)) )
4569  allocate( pott_org(dims(1),dims(2),dims(3)) )
4570 
4571  allocate( cz_org(dims(1)) )
4572  allocate( fz_org(dims(1)) )
4573  allocate( cx_org(dims(2)) )
4574  allocate( fx_org(dims(2)) )
4575  allocate( cy_org(dims(3)) )
4576  allocate( fy_org(dims(3)) )
4577 
4578  call fileread( dens_org(:,:,:), &
4579  basename_org, "DENS", 1, 1, single=.true. )
4580  call fileread( momz_org(:,:,:), &
4581  basename_org, "MOMZ", 1, 1, single=.true. )
4582  call fileread( momx_org(:,:,:), &
4583  basename_org, "MOMX", 1, 1, single=.true. )
4584  call fileread( momy_org(:,:,:), &
4585  basename_org, "MOMY", 1, 1, single=.true. )
4586  call fileread( rhot_org(:,:,:), &
4587  basename_org, "RHOT", 1, 1, single=.true. )
4588  do iq = 1, qa
4589  call fileread( qtrc_org(:,:,:,iq), &
4590  basename_org, tracer_name(iq), 1, 1, single=.true. )
4591  end do
4592 
4593  call fileread( cz_org(:), &
4594  basename_org, "z" , 1, 1, single=.true. )
4595  call fileread( cx_org(:), &
4596  basename_org, "x" , 1, 1, single=.true. )
4597  call fileread( cy_org(:), &
4598  basename_org, "y" , 1, 1, single=.true. )
4599  call fileread( fx_org(:), &
4600  basename_org, "xh", 1, 1, single=.true. )
4601  call fileread( fy_org(:), &
4602  basename_org, "yh", 1, 1, single=.true. )
4603 
4604  do k = ks, ke
4605  call interporation_fact( fact_cz0(k), fact_cz1(k), & ! (OUT)
4606  idx_cz0( k), idx_cz1(k), & ! (OUT)
4607  grid_cz(k), cz_org, dims(1), & ! (IN)
4608  .false. ) ! (IN)
4609  call interporation_fact( fact_fz0(k), fact_fz1(k), & ! (OUT)
4610  idx_fz0(k), idx_fz1(k), & ! (OUT)
4611  grid_fz(k), fz_org, dims(1), & ! (IN)
4612  .false. ) ! (IN)
4613  enddo
4614  do i = is, ie
4615  call interporation_fact( fact_cx0(i), fact_cx1(i), & ! (OUT)
4616  idx_cx0(i), idx_cx1(i), & ! (OUT)
4617  grid_cx(i), cx_org, dims(2), & ! (IN)
4618  .true. ) ! (IN)
4619  call interporation_fact( fact_fx0(i), fact_fx1(i), & ! (OUT)
4620  idx_fx0(i), idx_fx1(i), & ! (OUT)
4621  grid_fx(i), fx_org, dims(2), & ! (IN)
4622  .true. ) ! (IN)
4623  enddo
4624  do j = js, je
4625  call interporation_fact( fact_cy0(j), fact_cy1(j), & ! (OUT)
4626  idx_cy0(j), idx_cy1(j), & ! (OUT)
4627  grid_cy(j), cy_org, dims(3), & ! (IN)
4628  .true. ) ! (IN)
4629  call interporation_fact( fact_fy0(j), fact_fy1(j), & ! (OUT)
4630  idx_fy0(j), idx_fy1(j), & ! (OUT)
4631  grid_fy(j), fy_org, dims(3), & ! (IN)
4632  .true. ) ! (IN)
4633  enddo
4634 
4635 
4636  do j = 1, dims(3)
4637  do i = 1, dims(2)
4638  do k = 1, dims(1)-1
4639  w_org(k,i,j) = 2.0_rp * momz_org(k,i,j) / ( dens_org(k+1,i,j) + dens_org(k,i,j) )
4640  end do
4641  end do
4642  end do
4643  do j = 1, dims(3)
4644  do i = 1, dims(2)
4645  w_org(dims(1),i,j) = 0.0_rp
4646  end do
4647  end do
4648 
4649  do j = 1, dims(3)
4650  do i = 1, dims(2)-1
4651  do k = 1, dims(1)
4652  u_org(k,i,j) = 2.0_rp * momx_org(k,i,j) / ( dens_org(k,i+1,j) + dens_org(k,i,j) )
4653  end do
4654  end do
4655  end do
4656  do j = 1, dims(3)
4657  do k = 1, dims(1)
4658  u_org(k,dims(2),j) = 2.0_rp * momx_org(k,dims(2),j) / ( dens_org(k,1,j) + dens_org(k,dims(2),j) )
4659  end do
4660  end do
4661 
4662  do j = 1, dims(3)-1
4663  do i = 1, dims(2)
4664  do k = 1, dims(1)
4665  v_org(k,i,j) = 2.0_rp * momy_org(k,i,j) / ( dens_org(k,i,j+1) + dens_org(k,i,j) )
4666  end do
4667  end do
4668  end do
4669  do i = 1, dims(2)
4670  do k = 1, dims(1)
4671  v_org(k,i,dims(3)) = 2.0_rp * momy_org(k,i,dims(3)) / ( dens_org(k,i,1) + dens_org(k,i,dims(3)) )
4672  end do
4673  end do
4674 
4675  do j = 1, dims(3)
4676  do i = 1, dims(2)
4677  do k = 1, dims(1)
4678  pott_org(k,i,j) = rhot_org(k,i,j) / dens_org(k,i,j)
4679  end do
4680  end do
4681  end do
4682 
4683  do j = 1, dims(3)
4684  do i = 1, dims(2)
4685  do k = 1, dims(1)
4686  dens_org(k,i,j) = log( dens_org(k,i,j) )
4687  end do
4688  end do
4689  end do
4690 
4691  do j = js, je
4692  do i = is, ie
4693  do k = ks, ie
4694  dens(k,i,j) = exp( &
4695  fact_cz0(k)*fact_cx0(i)*fact_cy0(j)*dens_org(idx_cz0(k),idx_cx0(i),idx_cy0(j)) &
4696  + fact_cz1(k)*fact_cx0(i)*fact_cy0(j)*dens_org(idx_cz1(k),idx_cx0(i),idx_cy0(j)) &
4697  + fact_cz0(k)*fact_cx1(i)*fact_cy0(j)*dens_org(idx_cz0(k),idx_cx1(i),idx_cy0(j)) &
4698  + fact_cz1(k)*fact_cx1(i)*fact_cy0(j)*dens_org(idx_cz1(k),idx_cx1(i),idx_cy0(j)) &
4699  + fact_cz0(k)*fact_cx0(i)*fact_cy1(j)*dens_org(idx_cz0(k),idx_cx0(i),idx_cy1(j)) &
4700  + fact_cz1(k)*fact_cx0(i)*fact_cy1(j)*dens_org(idx_cz1(k),idx_cx0(i),idx_cy1(j)) &
4701  + fact_cz0(k)*fact_cx1(i)*fact_cy1(j)*dens_org(idx_cz0(k),idx_cx1(i),idx_cy1(j)) &
4702  + fact_cz1(k)*fact_cx1(i)*fact_cy1(j)*dens_org(idx_cz1(k),idx_cx1(i),idx_cy1(j)) &
4703  )
4704 
4705  w(k,i,j) = fact_fz0(k)*fact_cx0(i)*fact_cy0(j)*w_org(idx_fz0(k),idx_cx0(i),idx_cy0(j)) &
4706  + fact_fz1(k)*fact_cx0(i)*fact_cy0(j)*w_org(idx_fz1(k),idx_cx0(i),idx_cy0(j)) &
4707  + fact_fz0(k)*fact_cx1(i)*fact_cy0(j)*w_org(idx_fz0(k),idx_cx1(i),idx_cy0(j)) &
4708  + fact_fz1(k)*fact_cx1(i)*fact_cy0(j)*w_org(idx_fz1(k),idx_cx1(i),idx_cy0(j)) &
4709  + fact_fz0(k)*fact_cx0(i)*fact_cy1(j)*w_org(idx_fz0(k),idx_cx0(i),idx_cy1(j)) &
4710  + fact_fz1(k)*fact_cx0(i)*fact_cy1(j)*w_org(idx_fz1(k),idx_cx0(i),idx_cy1(j)) &
4711  + fact_fz0(k)*fact_cx1(i)*fact_cy1(j)*w_org(idx_fz0(k),idx_cx1(i),idx_cy1(j)) &
4712  + fact_fz1(k)*fact_cx1(i)*fact_cy1(j)*w_org(idx_fz1(k),idx_cx1(i),idx_cy1(j))
4713 
4714  u(k,i,j) = fact_cz0(k)*fact_fx0(i)*fact_cy0(j)*u_org(idx_cz0(k),idx_fx0(i),idx_cy0(j)) &
4715  + fact_cz1(k)*fact_fx0(i)*fact_cy0(j)*u_org(idx_cz1(k),idx_fx0(i),idx_cy0(j)) &
4716  + fact_cz0(k)*fact_fx1(i)*fact_cy0(j)*u_org(idx_cz0(k),idx_fx1(i),idx_cy0(j)) &
4717  + fact_cz1(k)*fact_fx1(i)*fact_cy0(j)*u_org(idx_cz1(k),idx_fx1(i),idx_cy0(j)) &
4718  + fact_cz0(k)*fact_fx0(i)*fact_cy1(j)*u_org(idx_cz0(k),idx_fx0(i),idx_cy1(j)) &
4719  + fact_cz1(k)*fact_fx0(i)*fact_cy1(j)*u_org(idx_cz1(k),idx_fx0(i),idx_cy1(j)) &
4720  + fact_cz0(k)*fact_fx1(i)*fact_cy1(j)*u_org(idx_cz0(k),idx_fx1(i),idx_cy1(j)) &
4721  + fact_cz1(k)*fact_fx1(i)*fact_cy1(j)*u_org(idx_cz1(k),idx_fx1(i),idx_cy1(j))
4722 
4723  v(k,i,j) = fact_cz0(k)*fact_cx0(i)*fact_fy0(j)*v_org(idx_cz0(k),idx_cx0(i),idx_fy0(j)) &
4724  + fact_cz1(k)*fact_cx0(i)*fact_fy0(j)*v_org(idx_cz1(k),idx_cx0(i),idx_fy0(j)) &
4725  + fact_cz0(k)*fact_cx1(i)*fact_fy0(j)*v_org(idx_cz0(k),idx_cx1(i),idx_fy0(j)) &
4726  + fact_cz1(k)*fact_cx1(i)*fact_fy0(j)*v_org(idx_cz1(k),idx_cx1(i),idx_fy0(j)) &
4727  + fact_cz0(k)*fact_cx0(i)*fact_fy1(j)*v_org(idx_cz0(k),idx_cx0(i),idx_fy1(j)) &
4728  + fact_cz1(k)*fact_cx0(i)*fact_fy1(j)*v_org(idx_cz1(k),idx_cx0(i),idx_fy1(j)) &
4729  + fact_cz0(k)*fact_cx1(i)*fact_fy1(j)*v_org(idx_cz0(k),idx_cx1(i),idx_fy1(j)) &
4730  + fact_cz1(k)*fact_cx1(i)*fact_fy1(j)*v_org(idx_cz1(k),idx_cx1(i),idx_fy1(j))
4731 
4732  pott(k,i,j) = fact_cz0(k)*fact_cx0(i)*fact_cy0(j)*pott_org(idx_cz0(k),idx_cx0(i),idx_cy0(j)) &
4733  + fact_cz1(k)*fact_cx0(i)*fact_cy0(j)*pott_org(idx_cz1(k),idx_cx0(i),idx_cy0(j)) &
4734  + fact_cz0(k)*fact_cx1(i)*fact_cy0(j)*pott_org(idx_cz0(k),idx_cx1(i),idx_cy0(j)) &
4735  + fact_cz1(k)*fact_cx1(i)*fact_cy0(j)*pott_org(idx_cz1(k),idx_cx1(i),idx_cy0(j)) &
4736  + fact_cz0(k)*fact_cx0(i)*fact_cy1(j)*pott_org(idx_cz0(k),idx_cx0(i),idx_cy1(j)) &
4737  + fact_cz1(k)*fact_cx0(i)*fact_cy1(j)*pott_org(idx_cz1(k),idx_cx0(i),idx_cy1(j)) &
4738  + fact_cz0(k)*fact_cx1(i)*fact_cy1(j)*pott_org(idx_cz0(k),idx_cx1(i),idx_cy1(j)) &
4739  + fact_cz1(k)*fact_cx1(i)*fact_cy1(j)*pott_org(idx_cz1(k),idx_cx1(i),idx_cy1(j))
4740 
4741  do iq = 1, qa
4742  qtrc(k,i,j,iq) = fact_cz0(k)*fact_cx0(i)*fact_cy0(j)*qtrc_org(idx_cz0(k),idx_cx0(i),idx_cy0(j),iq) &
4743  + fact_cz1(k)*fact_cx0(i)*fact_cy0(j)*qtrc_org(idx_cz1(k),idx_cx0(i),idx_cy0(j),iq) &
4744  + fact_cz0(k)*fact_cx1(i)*fact_cy0(j)*qtrc_org(idx_cz0(k),idx_cx1(i),idx_cy0(j),iq) &
4745  + fact_cz1(k)*fact_cx1(i)*fact_cy0(j)*qtrc_org(idx_cz1(k),idx_cx1(i),idx_cy0(j),iq) &
4746  + fact_cz0(k)*fact_cx0(i)*fact_cy1(j)*qtrc_org(idx_cz0(k),idx_cx0(i),idx_cy1(j),iq) &
4747  + fact_cz1(k)*fact_cx0(i)*fact_cy1(j)*qtrc_org(idx_cz1(k),idx_cx0(i),idx_cy1(j),iq) &
4748  + fact_cz0(k)*fact_cx1(i)*fact_cy1(j)*qtrc_org(idx_cz0(k),idx_cx1(i),idx_cy1(j),iq) &
4749  + fact_cz1(k)*fact_cx1(i)*fact_cy1(j)*qtrc_org(idx_cz1(k),idx_cx1(i),idx_cy1(j),iq)
4750  enddo
4751  enddo
4752  enddo
4753  enddo
4754 
4755  deallocate( dens_org )
4756  deallocate( momz_org )
4757  deallocate( momx_org )
4758  deallocate( momy_org )
4759  deallocate( rhot_org )
4760  deallocate( qtrc_org )
4761 
4762  deallocate( w_org )
4763  deallocate( u_org )
4764  deallocate( v_org )
4765  deallocate( pott_org )
4766 
4767  deallocate( cz_org )
4768  deallocate( fz_org )
4769  deallocate( cx_org )
4770  deallocate( fx_org )
4771  deallocate( cy_org )
4772  deallocate( fy_org )
4773 
4774  if ( i_qv > 0 ) then
4775  do j = js, je
4776  do i = is, ie
4777  do k = ks, ke
4778  qv(k,i,j) = qtrc(k,i,j,i_qv)
4779  end do
4780  end do
4781  end do
4782  else
4783  do j = js, je
4784  do i = is, ie
4785  do k = ks, ke
4786  qv(k,i,j) = 0.0_rp
4787  end do
4788  end do
4789  end do
4790  end if
4791 
4792  if ( i_qc > 0 ) then
4793  do j = js, je
4794  do i = is, ie
4795  do k = ks, ke
4796  qc(k,i,j) = qtrc(k,i,j,i_qc)
4797  end do
4798  end do
4799  end do
4800  else
4801  do j = js, je
4802  do i = is, ie
4803  do k = ks, ke
4804  qc(k,i,j) = 0.0_rp
4805  end do
4806  end do
4807  end do
4808  end if
4809 
4810  do j = js, je
4811  do i = is, ie
4812  do k = ks+1, ke
4813  dz(k,i,j) = real_cz(k,i,j) - real_cz(k-1,i,j) ! distance from cell center to cell center
4814  enddo
4815  enddo
4816  enddo
4817 
4818  ! make density & pressure profile in moist condition
4819  call hydrostatic_buildrho_atmos( dens(:,:,:), & ! [INOUT]
4820  temp(:,:,:), & ! [OUT]
4821  pres(:,:,:), & ! [OUT]
4822  pott(:,:,:), & ! [IN]
4823  qv(:,:,:), & ! [IN]
4824  qc(:,:,:), & ! [IN]
4825  dz(:,:,:) ) ! [IN]
4826 
4827  do j = js, je
4828  do i = is, ie
4829  do k = ks, ke
4830  momz(k,i,j) = 0.5_rp * w(k,i,j) * ( dens(k,i,j) + dens(k+1,i,j) )
4831  enddo
4832  enddo
4833  enddo
4834 
4835  do j = js, je
4836  do i = is, ie
4837  do k = ks, ke
4838  momx(k,i,j) = 0.5_rp * u(k,i,j) * ( dens(k,i+1,j) + dens(k,i,j) )
4839  enddo
4840  enddo
4841  enddo
4842 
4843  do j = js, je
4844  do i = is, ie
4845  do k = ks, ke
4846  momy(k,i,j) = 0.5_rp * v(k,i,j) * ( dens(k,i,j+1) + dens(k,i,j) )
4847  enddo
4848  enddo
4849  enddo
4850 
4851  do j = js, je
4852  do i = is, ie
4853  do k = ks, ke
4854  rhot(k,i,j) = pott(k,i,j) * dens(k,i,j)
4855  enddo
4856  enddo
4857  enddo
4858 
4859  return
4860  end subroutine mkinit_interporation
4861 
4862  !-----------------------------------------------------------------------------
4863  subroutine interporation_fact( &
4864  fact0, fact1, &
4865  idx0, idx1, &
4866  x, x_org, nx, &
4867  loop )
4868  implicit none
4869 
4870  real(RP), intent(out) :: fact0
4871  real(RP), intent(out) :: fact1
4872  integer, intent(out) :: idx0
4873  integer, intent(out) :: idx1
4874  real(RP), intent(in) :: x
4875  integer, intent(in) :: nx
4876  real(RP), intent(in) :: x_org(nx)
4877  logical, intent(in) :: loop
4878 
4879  real(RP) :: xwork
4880  integer :: i
4881 
4882  if ( x < x_org(1) ) then
4883  if ( loop ) then
4884  xwork = x_org(1) - ( x_org(2) - x_org(1) )**2 / ( x_org(3) - x_org(2) )
4885  fact0 = ( x_org(1) - x ) / ( x_org(1) - xwork )
4886  fact1 = ( x - xwork ) / ( x_org(1) - xwork )
4887  idx0 = nx
4888  idx1 = 1
4889  else
4890  fact0 = ( x_org(2) - x ) / ( x_org(2) - x_org(1) )
4891  fact1 = ( x - x_org(1) ) / ( x_org(2) - x_org(1) )
4892  idx0 = 1
4893  idx1 = 2
4894  end if
4895  else if ( x > x_org(nx) ) then
4896  if ( loop ) then
4897  xwork = x_org(nx) + ( x_org(nx) - x_org(nx-1) )**2 / ( x_org(nx-1) - x_org(nx-2) )
4898  fact0 = ( xwork - x ) / ( xwork - x_org(nx) )
4899  fact1 = ( x - x_org(nx) ) / ( xwork - x_org(nx) )
4900  idx0 = nx
4901  idx1 = 1
4902  else
4903  fact0 = ( x_org(nx) - x ) / ( x_org(nx) - x_org(nx-1) )
4904  fact1 = ( x - x_org(nx-1) ) / ( x_org(nx) - x_org(nx-1) )
4905  idx0 = nx-1
4906  idx1 = nx
4907  end if
4908  else
4909  do i = 2, nx
4910  if ( x <= x_org(i) ) then
4911  fact0 = ( x_org(i) - x ) / ( x_org(i) - x_org(i-1) )
4912  fact1 = ( x - x_org(i-1) ) / ( x_org(i) - x_org(i-1) )
4913  idx0 = i-1
4914  idx1 = i
4915  exit
4916  end if
4917  end do
4918  end if
4919 
4920  return
4921  end subroutine interporation_fact
4922  !-----------------------------------------------------------------------------
4924  subroutine mkinit_oceancouple
4925  implicit none
4926 
4927  if( io_l ) write(io_fid_log,*)
4928  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit OceanCouple] / Categ[preprocess] / Origin[SCALE-RM]'
4929 
4930  call flux_setup
4931 
4932  call ocean_setup
4933 
4934  return
4935  end subroutine mkinit_oceancouple
4936 
4937  !-----------------------------------------------------------------------------
4939  subroutine mkinit_landcouple
4940  implicit none
4941 
4942  if( io_l ) write(io_fid_log,*)
4943  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit LandCouple] / Categ[preprocess] / Origin[SCALE-RM]'
4944 
4945  call flux_setup
4946 
4947  call land_setup
4948 
4949  return
4950  end subroutine mkinit_landcouple
4951 
4952  !-----------------------------------------------------------------------------
4954  subroutine mkinit_urbancouple
4955  implicit none
4956 
4957  if( io_l ) write(io_fid_log,*)
4958  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit UrbanCouple] / Categ[preprocess] / Origin[SCALE-RM]'
4959 
4960  call flux_setup
4961 
4962  call urban_setup
4963 
4964  return
4965  end subroutine mkinit_urbancouple
4966 
4967  !-----------------------------------------------------------------------------
4969  subroutine mkinit_seabreeze
4970  use scale_rm_process, only: &
4971  prc_num_x
4972  use scale_landuse, only: &
4973  landuse_calc_fact, &
4974  landuse_frac_land
4975  use scale_grid, only: &
4976  grid_domain_center_x
4977  implicit none
4978 
4979  real(RP) :: LAND_SIZE
4980 
4981  namelist / param_mkinit_seabreeze / &
4982  land_size
4983 
4984  integer :: ierr
4985  integer :: i, j
4986  !---------------------------------------------------------------------------
4987 
4988  if( io_l ) write(io_fid_log,*)
4989  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit SEABREEZE] / Categ[preprocess] / Origin[SCALE-RM]'
4990 
4991  land_size = 0.0_rp
4992 
4993  !--- read namelist
4994  rewind(io_fid_conf)
4995  read(io_fid_conf,nml=param_mkinit_seabreeze,iostat=ierr)
4996 
4997  if( ierr < 0 ) then !--- missing
4998  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
4999  elseif( ierr > 0 ) then !--- fatal error
5000  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_SEABREEZE. Check!'
5001  call prc_mpistop
5002  endif
5003  if( io_nml ) write(io_fid_nml,nml=param_mkinit_seabreeze)
5004 
5005  call flux_setup
5006 
5007  call land_setup
5008 
5009  call ocean_setup
5010 
5011  ! make landuse conditions
5012  do j = js, je
5013  do i = is, ie
5014  if ( abs( grid_cx(i) - grid_domain_center_x ) < land_size ) then
5015  landuse_frac_land(i,j) = 1.0_rp
5016  else
5017  landuse_frac_land(i,j) = 0.0_rp
5018  endif
5019  enddo
5020  enddo
5021 
5022  call landuse_calc_fact
5023 
5024  return
5025  end subroutine mkinit_seabreeze
5026 
5027  !-----------------------------------------------------------------------------
5029  subroutine mkinit_heatisland
5030  use scale_rm_process, only: &
5031  prc_num_x
5032  use scale_landuse, only: &
5033  landuse_calc_fact, &
5034  landuse_frac_land, &
5035  landuse_frac_urban
5036  implicit none
5037 
5038  real(RP) :: dist
5039 
5040  integer :: i, j
5041  !---------------------------------------------------------------------------
5042 
5043  if( io_l ) write(io_fid_log,*)
5044  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit HEATISLAND] / Categ[preprocess] / Origin[SCALE-RM]'
5045 
5046  call flux_setup
5047 
5048  call land_setup
5049 
5050  call urban_setup
5051 
5052  ! 1/9 size of domain
5053  dist = ( grid_cxg(imax*prc_num_x) - grid_cxg(1) ) / 9.0_rp
5054 
5055  ! make landuse conditions
5056  do j = js, je
5057  do i = is, ie
5058  if ( grid_cx(i) >= dist * 4.0_rp &
5059  .AND. grid_cx(i) < dist * 5.0_rp ) then
5060  landuse_frac_land(i,j) = 1.0_rp
5061  landuse_frac_urban(i,j) = 1.0_rp
5062  else
5063  landuse_frac_land(i,j) = 1.0_rp
5064  landuse_frac_urban(i,j) = 0.0_rp
5065  endif
5066  enddo
5067  enddo
5068 
5069  call landuse_calc_fact
5070 
5071  return
5072  end subroutine mkinit_heatisland
5073 
5074  !-----------------------------------------------------------------------------
5076  subroutine mkinit_grayzone
5077  use scale_atmos_hydrometeor, only: &
5078  i_qv
5079  implicit none
5080 
5081  real(RP) :: RHO(KA)
5082  real(RP) :: VELX(KA)
5083  real(RP) :: VELY(KA)
5084  real(RP) :: POTT(KA)
5085  real(RP) :: QV(KA)
5086 
5087  real(RP) :: PERTURB_AMP = 0.0_rp
5088  integer :: RANDOM_LIMIT = 0
5089  integer :: RANDOM_FLAG = 0 ! 0 -> no perturbation
5090  ! 1 -> petrurbation for pt
5091  ! 2 -> perturbation for u, v, w
5092 
5093  namelist / param_mkinit_grayzone / &
5094  perturb_amp, &
5095  random_limit, &
5096  random_flag
5097 
5098  integer :: ierr
5099  integer :: k, i, j
5100  !---------------------------------------------------------------------------
5101 
5102  if( io_l ) write(io_fid_log,*)
5103  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit GRAYZONE] / Categ[preprocess] / Origin[SCALE-RM]'
5104 
5105  if ( i_qv < 1 ) then
5106  write(*,*) 'xxx QV is not registered'
5107  call prc_mpistop
5108  end if
5109 
5110  !--- read namelist
5111  rewind(io_fid_conf)
5112  read(io_fid_conf,nml=param_mkinit_grayzone,iostat=ierr)
5113 
5114  if( ierr < 0 ) then !--- missing
5115  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
5116  elseif( ierr > 0 ) then !--- fatal error
5117  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_GRAYZONE. Check!'
5118  call prc_mpistop
5119  endif
5120  if( io_nml ) write(io_fid_nml,nml=param_mkinit_grayzone)
5121 
5122  call read_sounding( rho, velx, vely, pott, qv ) ! (out)
5123 
5124 ! do j = JS, JE
5125 ! do i = IS, IE
5126  do j = 1, ja
5127  do i = 1, ia
5128  do k = ks, ke
5129  dens(k,i,j) = rho(k)
5130 ! MOMZ(k,i,j) = 0.0_RP
5131 ! MOMX(k,i,j) = RHO(k) * VELX(k)
5132 ! MOMY(k,i,j) = RHO(k) * VELY(k)
5133 
5134 ! RHOT(k,i,j) = RHO(k) * POTT(k)
5135 
5136  qtrc(k,i,j,i_qv) = qv(k)
5137  enddo
5138  enddo
5139  enddo
5140 
5141  do j = js, je
5142  do i = is, ie
5143  dens( 1:ks-1,i,j) = dens(ks,i,j)
5144  dens(ke+1:ka, i,j) = dens(ke,i,j)
5145  enddo
5146  enddo
5147 
5148  call random_get(rndm) ! make random
5149  do j = js, je
5150  do i = is, ie
5151  do k = ks, ke
5152  if ( random_flag == 2 .and. k <= random_limit ) then ! below initial cloud top
5153  momz(k,i,j) = ( 2.0_rp * ( rndm(k,i,j)-0.5_rp ) * perturb_amp ) &
5154  * 0.5_rp * ( dens(k+1,i,j) + dens(k,i,j) )
5155  else
5156  momz(k,i,j) = 0.0_rp
5157  endif
5158  enddo
5159  enddo
5160  enddo
5161 
5162  call random_get(rndm) ! make random
5163  do j = js, je
5164  do i = is, ie
5165  do k = ks, ke
5166  if ( random_flag == 2 .AND. k <= random_limit ) then ! below initial cloud top
5167  momx(k,i,j) = ( velx(k) + 2.0_rp * ( rndm(k,i,j)-0.5_rp ) * perturb_amp ) &
5168  * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
5169  else
5170  momx(k,i,j) = velx(k) * 0.5_rp * ( dens(k,i+1,j) + dens(k,i,j) )
5171  endif
5172  enddo
5173  enddo
5174  enddo
5175 
5176  call random_get(rndm) ! make random
5177  do j = js, je
5178  do i = is, ie
5179  do k = ks, ke
5180  if ( random_flag == 2 .AND. k <= random_limit ) then ! below initial cloud top
5181  momy(k,i,j) = ( vely(k) + 2.0_rp * ( rndm(k,i,j)-0.5_rp ) * perturb_amp ) &
5182  * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
5183  else
5184  momy(k,i,j) = vely(k) * 0.5_rp * ( dens(k,i,j+1) + dens(k,i,j) )
5185  endif
5186  enddo
5187  enddo
5188  enddo
5189 
5190  call random_get(rndm) ! make random
5191  do j = js, je
5192  do i = is, ie
5193  do k = ks, ke
5194  if ( random_flag == 1 .and. k <= random_limit ) then ! below initial cloud top
5195  rhot(k,i,j) = ( pott(k) + 2.0_rp * ( rndm(k,i,j)-0.5_rp ) * perturb_amp ) &
5196  * dens(k,i,j)
5197  else
5198  rhot(k,i,j) = pott(k) * dens(k,i,j)
5199  endif
5200  enddo
5201  enddo
5202  enddo
5203 
5204  return
5205  end subroutine mkinit_grayzone
5206 
5207  !-----------------------------------------------------------------------------
5209  subroutine mkinit_boxaero
5210  use scale_atmos_hydrometeor, only: &
5211  i_qv
5212  use mod_atmos_admin, only: &
5213  atmos_phy_ae_type
5214  implicit none
5215 
5216  real(RP) :: init_dens = 1.12_rp ![kg/m3]
5217  real(RP) :: init_temp = 298.18_rp ![K]
5218  real(RP) :: init_pres = 1.e+5_rp ![Pa]
5219  real(RP) :: init_ssliq = 0.01_rp ![%]
5220 
5221  namelist / param_mkinit_boxaero / &
5222  init_dens, &
5223  init_temp, &
5224  init_pres, &
5225  init_ssliq
5226 
5227  real(RP) :: qsat
5228  integer :: i, j, k, ierr
5229  !---------------------------------------------------------------------------
5230 
5231  if ( atmos_phy_ae_type /= 'KAJINO13' ) then
5232  if( io_l ) write(io_fid_log,*) '+++ For [Box model of aerosol],'
5233  if( io_l ) write(io_fid_log,*) '+++ ATMOS_PHY_AE_TYPE should be KAJINO13. Stop!'
5234  call prc_mpistop
5235  endif
5236 
5237  if ( i_qv < 1 ) then
5238  write(*,*) 'xxx QV is not registered'
5239  call prc_mpistop
5240  end if
5241 
5242  if( io_l ) write(io_fid_log,*)
5243  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit Box model of aerosol] / Categ[preprocess] / Origin[SCALE-RM]'
5244 
5245  !--- read namelist
5246  rewind(io_fid_conf)
5247  read(io_fid_conf,nml=param_mkinit_boxaero,iostat=ierr)
5248  if( ierr < 0 ) then !--- missing
5249  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
5250  elseif( ierr > 0 ) then !--- fatal error
5251  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_BOXAERO. Check!'
5252  call prc_mpistop
5253  endif
5254  if( io_nml ) write(io_fid_nml,nml=param_mkinit_boxaero)
5255 
5256  qtrc(:,:,:,:) = 0.0_rp
5257  call saturation_pres2qsat_all( qsat, init_temp, init_pres )
5258 
5259  do j = 1, ja
5260  do i = 1, ia
5261  do k = 1, ka
5262  dens(k,i,j) = init_dens
5263  momx(k,i,j) = 0.0_rp
5264  momy(k,i,j) = 0.0_rp
5265  momz(k,i,j) = 0.0_rp
5266  pott(k,i,j) = init_temp * ( p00/init_pres )**(rdry/cpdry)
5267  rhot(k,i,j) = init_dens * pott(k,i,j)
5268 
5269  qtrc(k,i,j,i_qv) = ( init_ssliq + 1.0_rp ) * qsat
5270  enddo
5271  enddo
5272  enddo
5273 
5274  if ( atmos_phy_ae_type == 'KAJINO13' ) then
5275  call aerosol_setup
5276  endif
5277 
5278  return
5279  end subroutine mkinit_boxaero
5280 
5281  !-----------------------------------------------------------------------------
5283  subroutine mkinit_warmbubbleaero
5284  use scale_atmos_hydrometeor, only: &
5285  i_qv
5286  use mod_atmos_admin, only: &
5287  atmos_phy_ae_type
5288  implicit none
5289 
5290  ! Surface state
5291  real(RP) :: SFC_THETA ! surface potential temperature [K]
5292  real(RP) :: SFC_PRES ! surface pressure [Pa]
5293  real(RP) :: SFC_RH = 80.0_rp ! surface relative humidity [%]
5294  ! Environment state
5295  real(RP) :: ENV_U = 0.0_rp ! velocity u of environment [m/s]
5296  real(RP) :: ENV_V = 0.0_rp ! velocity v of environment [m/s]
5297  real(RP) :: ENV_RH = 80.0_rp ! Relative Humidity of environment [%]
5298  real(RP) :: ENV_L1_ZTOP = 1.e3_rp ! top height of the layer1 (constant THETA) [m]
5299  real(RP) :: ENV_L2_ZTOP = 14.e3_rp ! top height of the layer2 (small THETA gradient) [m]
5300  real(RP) :: ENV_L2_TLAPS = 4.e-3_rp ! Lapse rate of THETA in the layer2 (small THETA gradient) [K/m]
5301  real(RP) :: ENV_L3_TLAPS = 3.e-2_rp ! Lapse rate of THETA in the layer3 (large THETA gradient) [K/m]
5302  ! Bubble
5303  real(RP) :: BBL_THETA = 1.0_rp ! extremum of temperature in bubble [K]
5304 
5305  namelist / param_mkinit_warmbubble / &
5306  sfc_theta, &
5307  sfc_pres, &
5308  env_u, &
5309  env_v, &
5310  env_rh, &
5311  env_l1_ztop, &
5312  env_l2_ztop, &
5313  env_l2_tlaps, &
5314  env_l3_tlaps, &
5315  bbl_theta
5316 
5317  integer :: ierr
5318  integer :: k, i, j
5319  !---------------------------------------------------------------------------
5320 
5321  if( io_l ) write(io_fid_log,*)
5322  if( io_l ) write(io_fid_log,*) '++++++ Module[mkinit WARMBUBBLEAERO] / Categ[preprocess] / Origin[SCALE-RM]'
5323 
5324  if ( i_qv < 1 ) then
5325  write(*,*) 'xxx QV is not registerd'
5326  call prc_mpistop
5327  end if
5328 
5329 
5330  sfc_theta = thetastd
5331  sfc_pres = pstd
5332 
5333  !--- read namelist
5334  rewind(io_fid_conf)
5335  read(io_fid_conf,nml=param_mkinit_warmbubble,iostat=ierr)
5336 
5337  if( ierr < 0 ) then !--- missing
5338  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
5339  elseif( ierr > 0 ) then !--- fatal error
5340  write(*,*) 'xxx Not appropriate names in namelist PARAM_MKINIT_WARMBUBBLE. Check!'
5341  call prc_mpistop
5342  endif
5343  if( io_nml ) write(io_fid_nml,nml=param_mkinit_warmbubble)
5344 
5345  ! calc in dry condition
5346  pres_sfc(1,1,1) = sfc_pres
5347  pott_sfc(1,1,1) = sfc_theta
5348  qv_sfc(1,1,1) = 0.0_rp
5349  qc_sfc(1,1,1) = 0.0_rp
5350 
5351  do k = ks, ke
5352  if ( grid_cz(k) <= env_l1_ztop ) then ! Layer 1
5353  pott(k,1,1) = sfc_theta
5354  elseif( grid_cz(k) < env_l2_ztop ) then ! Layer 2
5355  pott(k,1,1) = pott(k-1,1,1) + env_l2_tlaps * ( grid_cz(k)-grid_cz(k-1) )
5356  else ! Layer 3
5357  pott(k,1,1) = pott(k-1,1,1) + env_l3_tlaps * ( grid_cz(k)-grid_cz(k-1) )
5358  endif
5359  qv(k,1,1) = 0.0_rp
5360  qc(k,1,1) = 0.0_rp
5361  enddo
5362 
5363  ! make density & pressure profile in dry condition
5364  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
5365  temp(:,1,1), & ! [OUT]
5366  pres(:,1,1), & ! [OUT]
5367  pott(:,1,1), & ! [IN]
5368  qv(:,1,1), & ! [IN]
5369  qc(:,1,1), & ! [IN]
5370  temp_sfc(1,1,1), & ! [OUT]
5371  pres_sfc(1,1,1), & ! [IN]
5372  pott_sfc(1,1,1), & ! [IN]
5373  qv_sfc(1,1,1), & ! [IN]
5374  qc_sfc(1,1,1) ) ! [IN]
5375 
5376  ! calc QV from RH
5377  call saturation_pres2qsat_all( qsat_sfc(1,1,1), temp_sfc(1,1,1), pres_sfc(1,1,1) )
5378  call saturation_pres2qsat_all( qsat(:,1,1), temp(:,1,1), pres(:,1,1) )
5379 
5380  qv_sfc(1,1,1) = sfc_rh * 1.e-2_rp * qsat_sfc(1,1,1)
5381  do k = ks, ke
5382  if ( grid_cz(k) <= env_l1_ztop ) then ! Layer 1
5383  qv(k,1,1) = env_rh * 1.e-2_rp * qsat(k,1,1)
5384  elseif( grid_cz(k) <= env_l2_ztop ) then ! Layer 2
5385  qv(k,1,1) = env_rh * 1.e-2_rp * qsat(k,1,1)
5386  else ! Layer 3
5387  qv(k,1,1) = 0.0_rp
5388  endif
5389  enddo
5390 
5391  ! make density & pressure profile in moist condition
5392  call hydrostatic_buildrho( dens(:,1,1), & ! [OUT]
5393  temp(:,1,1), & ! [OUT]
5394  pres(:,1,1), & ! [OUT]
5395  pott(:,1,1), & ! [IN]
5396  qv(:,1,1), & ! [IN]
5397  qc(:,1,1), & ! [IN]
5398  temp_sfc(1,1,1), & ! [OUT]
5399  pres_sfc(1,1,1), & ! [IN]
5400  pott_sfc(1,1,1), & ! [IN]
5401  qv_sfc(1,1,1), & ! [IN]
5402  qc_sfc(1,1,1) ) ! [IN]
5403 
5404  do j = js, je
5405  do i = is, ie
5406  do k = ks, ke
5407  dens(k,i,j) = dens(k,1,1)
5408  momz(k,i,j) = 0.0_rp
5409  momx(k,i,j) = env_u * dens(k,i,j)
5410  momy(k,i,j) = env_v * dens(k,i,j)
5411 
5412  ! make warm bubble
5413  rhot(k,i,j) = dens(k,1,1) * ( pott(k,1,1) + bbl_theta * bubble(k,i,j) )
5414 
5415  qtrc(k,i,j,i_qv) = qv(k,1,1)
5416  enddo
5417  enddo
5418  enddo
5419 
5420  call flux_setup
5421 
5422  if( atmos_phy_ae_type == 'KAJINO13' ) then
5423  call aerosol_setup
5424  endif
5425 
5426  return
5427  end subroutine mkinit_warmbubbleaero
5428 
5429  !-----------------------------------------------------------------------------
5431  subroutine mkinit_real
5432  use mod_realinput, only: &
5433  realinput_atmos, &
5434  realinput_surface
5435  implicit none
5436 
5437  call realinput_atmos( flg_intrp )
5438 
5439  call realinput_surface
5440 
5441  call flux_setup
5442 
5443  return
5444  end subroutine mkinit_real
5445 
5446 end module mod_mkinit
mod_atmos_admin
module ATMOS admin
Definition: mod_atmos_admin.f90:10
scale_grid_index::imax
integer, public imax
of computational cells: x, local
Definition: scale_grid_index.F90:44
scale_rm_process::prc_num_x
integer, public prc_num_x
x length of 2D processor topology
Definition: scale_rm_process.F90:43
scale_grid_index::is
integer, public is
start point of inner domain: x, local
Definition: scale_grid_index.F90:60
gtool_file
module GTOOL_FILE
Definition: gtool_file.f90:17
mod_atmos_vars::momz
real(rp), dimension(:,:,:), allocatable, target, public momz
Definition: mod_atmos_vars.f90:73
scale_grid::grid_domain_center_x
real(rp), public grid_domain_center_x
center position of global domain [m]: x
Definition: scale_grid_cartesian.F90:53
scale_grid_index::je
integer, public je
end point of inner domain: y, local
Definition: scale_grid_index.F90:63
scale_const::const_cpdry
real(rp), public const_cpdry
specific heat (dry air,constant pressure) [J/kg/K]
Definition: scale_const.F90:58
scale_atmos_saturation
module ATMOSPHERE / Saturation adjustment
Definition: scale_atmos_sub_saturation.F90:16
scale_atmos_phy_tb::i_tke
integer, public i_tke
Definition: scale_atmos_phy_tb.F90:106
mod_ocean_vars::ocean_sfc_albedo
real(rp), dimension(:,:,:), allocatable, public ocean_sfc_albedo
ocean surface albedo (0-1)
Definition: mod_ocean_vars.f90:61
scale_const::const_i_lw
integer, public const_i_lw
long-wave radiation index
Definition: scale_const.F90:95
scale_process::prc_mpistop
subroutine, public prc_mpistop
Abort MPI.
Definition: scale_process.F90:267
mod_urban_vars::urban_qc
real(rp), dimension(:,:), allocatable, public urban_qc
Definition: mod_urban_vars.f90:64
scale_atmos_phy_mp_suzuki10::nbin
integer, public nbin
Definition: scale_atmos_phy_mp_suzuki10.F90:109
mod_mkinit::i_oceancouple
integer, parameter, public i_oceancouple
Definition: mod_mkinit.f90:129
mod_atmos_vars::rhot
real(rp), dimension(:,:,:), allocatable, target, public rhot
Definition: mod_atmos_vars.f90:76
scale_atmos_hydrometeor::i_qc
integer, public i_qc
Definition: scale_atmos_sub_hydrometeor.F90:79
scale_grid::grid_cxg
real(rp), dimension(:), allocatable, public grid_cxg
center coordinate [m]: x, global
Definition: scale_grid_cartesian.F90:85
mod_atmos_phy_mp_vars
module Atmosphere / Physics Cloud Microphysics
Definition: mod_atmos_phy_mp_vars.f90:15
mod_urban_driver
module URBAN driver
Definition: mod_urban_driver.f90:10
scale_stdio::io_l
logical, public io_l
output log or not? (this process)
Definition: scale_stdio.F90:61
mod_ocean_vars::ocean_sfc_z0e
real(rp), dimension(:,:), allocatable, public ocean_sfc_z0e
ocean surface roughness length for vapor [m]
Definition: mod_ocean_vars.f90:64
scale_grid::grid_cz
real(rp), dimension(:), allocatable, public grid_cz
center coordinate [m]: z, local=global
Definition: scale_grid_cartesian.F90:56
mod_ocean_vars::ocean_temp
real(rp), dimension(:,:), allocatable, public ocean_temp
temperature at uppermost ocean layer [K]
Definition: mod_ocean_vars.f90:59
mod_mkinit::i_boxaero
integer, parameter, public i_boxaero
Definition: mod_mkinit.f90:142
scale_grid::grid_fxg
real(rp), dimension(:), allocatable, public grid_fxg
face coordinate [m]: x, global
Definition: scale_grid_cartesian.F90:83
mod_mkinit::land_setup
subroutine land_setup
Land setup.
Definition: mod_mkinit.f90:957
mod_atmos_vars
module ATMOSPHERIC Variables
Definition: mod_atmos_vars.f90:18
scale_const::const_radius
real(rp), public const_radius
radius of the planet [m]
Definition: scale_const.F90:46
mod_atmos_vars::momx
real(rp), dimension(:,:,:), allocatable, target, public momx
Definition: mod_atmos_vars.f90:74
scale_stdio
module STDIO
Definition: scale_stdio.F90:12
mod_mkinit::ocean_setup
subroutine ocean_setup
Ocean setup.
Definition: mod_mkinit.f90:1011
scale_landuse::landuse_calc_fact
subroutine, public landuse_calc_fact
Definition: scale_landuse.F90:156
mod_mkinit::i_dycoms2_rf02_dns
integer, parameter, public i_dycoms2_rf02_dns
Definition: mod_mkinit.f90:137
scale_grid_index::ke
integer, public ke
end point of inner domain: z, local
Definition: scale_grid_index.F90:59
mod_admin_time::time_doatmos_restart
logical, public time_doatmos_restart
execute atmosphere restart output in this step?
Definition: mod_admin_time.f90:54
mod_admin_time::time_doocean_restart
logical, public time_doocean_restart
execute ocean restart output in this step?
Definition: mod_admin_time.f90:56
scale_tracer::qa
integer, public qa
Definition: scale_tracer.F90:35
scale_grid_real::real_fz
real(rp), dimension(:,:,:), allocatable, public real_fz
geopotential height [m] (cell face )
Definition: scale_grid_real.F90:41
mod_mkinit::i_turbulence
integer, parameter, public i_turbulence
Definition: mod_mkinit.f90:115
mod_mkinit::mkinit
subroutine, public mkinit
Driver.
Definition: mod_mkinit.f90:358
mod_mkinit::i_urbancouple
integer, parameter, public i_urbancouple
Definition: mod_mkinit.f90:130
scale_grid_real::real_cz
real(rp), dimension(:,:,:), allocatable, public real_cz
geopotential height [m] (cell center)
Definition: scale_grid_real.F90:40
mod_mkinit::i_gravitywave
integer, parameter, public i_gravitywave
Definition: mod_mkinit.f90:113
mod_atmos_admin::atmos_phy_ae_type
character(len=h_short), public atmos_phy_ae_type
Definition: mod_atmos_admin.f90:36
mod_mkinit::faero
real(rp) function faero(f0, r0, x, alpha, rhoa)
Definition: mod_mkinit.f90:868
mod_urban_vars::urban_tb
real(rp), dimension(:,:), allocatable, public urban_tb
Definition: mod_urban_vars.f90:61
mod_mkinit::i_grayzone
integer, parameter, public i_grayzone
Definition: mod_mkinit.f90:141
mod_atmos_vars::dens
real(rp), dimension(:,:,:), allocatable, target, public dens
Definition: mod_atmos_vars.f90:72
mod_mkinit::interporation_fact
subroutine interporation_fact(fact0, fact1, idx0, idx1, x, x_org, nx, loop)
Definition: mod_mkinit.f90:4868
scale_landuse::landuse_frac_urban
real(rp), dimension(:,:), allocatable, public landuse_frac_urban
urban fraction
Definition: scale_landuse.F90:41
scale_atmos_hydrometeor::qhs
integer, public qhs
Definition: scale_atmos_sub_hydrometeor.F90:94
scale_atmos_hydrometeor::i_nc
integer, public i_nc
Definition: scale_atmos_sub_hydrometeor.F90:86
scale_atmos_hydrometeor::qhe
integer, public qhe
Definition: scale_atmos_sub_hydrometeor.F90:95
mod_urban_vars
module URBAN Variables
Definition: mod_urban_vars.f90:11
mod_urban_vars::urban_raing
real(rp), dimension(:,:), allocatable, public urban_raing
Definition: mod_urban_vars.f90:71
mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_lw_up
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_lw_up
Definition: mod_atmos_phy_rd_vars.f90:58
scale_tracer::tracer_name
character(len=h_short), dimension(qa_max), public tracer_name
Definition: scale_tracer.F90:44
mod_mkinit::i_landcouple
integer, parameter, public i_landcouple
Definition: mod_mkinit.f90:128
mod_atmos_phy_mp_vars::atmos_phy_mp_sflx_rain
real(rp), dimension(:,:), allocatable, public atmos_phy_mp_sflx_rain
Definition: mod_atmos_phy_mp_vars.f90:64
mod_urban_vars::urban_uc
real(rp), dimension(:,:), allocatable, public urban_uc
Definition: mod_urban_vars.f90:65
scale_grid::grid_fx
real(rp), dimension(:), allocatable, public grid_fx
face coordinate [m]: x, local
Definition: scale_grid_cartesian.F90:67
scale_const::const_rdry
real(rp), public const_rdry
specific gas constant (dry air) [J/kg/K]
Definition: scale_const.F90:57
mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_sw_up
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_sw_up
Definition: mod_atmos_phy_rd_vars.f90:65
mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_lw_dn
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_lw_dn
Definition: mod_atmos_phy_rd_vars.f90:59
mod_urban_driver::urban_surface_set
subroutine, public urban_surface_set(countup)
Set surface boundary to other model.
Definition: mod_urban_driver.f90:290
mod_atmos_phy_rd_vars
module Atmosphere / Physics Radiation
Definition: mod_atmos_phy_rd_vars.f90:14
scale_grid_index
module grid index
Definition: scale_grid_index.F90:14
mod_land_vars::land_sfc_temp
real(rp), dimension(:,:), allocatable, public land_sfc_temp
land surface skin temperature [K]
Definition: mod_land_vars.f90:64
scale_atmos_phy_ae_kajino13::atmos_phy_ae_kajino13_mkinit
subroutine, public atmos_phy_ae_kajino13_mkinit(QTRC, CCN, DENS, RHOT, m0_init, dg_init, sg_init, d_min_inp, d_max_inp, k_min_inp, k_max_inp, n_kap_inp)
Definition: scale_atmos_phy_ae_kajino13.F90:2274
scale_const::const_ohm
real(rp), public const_ohm
angular velocity of the planet [1/s]
Definition: scale_const.F90:47
scale_stdio::io_nml
logical, public io_nml
output log or not? (for namelist, this process)
Definition: scale_stdio.F90:62
mod_atmos_driver::atmos_surface_get
subroutine, public atmos_surface_get
Get surface boundary condition.
Definition: mod_atmos_driver.f90:471
scale_tracer
module TRACER
Definition: scale_tracer.F90:15
mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_sw_up
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_sw_up
Definition: mod_atmos_phy_rd_vars.f90:60
mod_mkinit::tke_setup
subroutine tke_setup
TKE setup.
Definition: mod_mkinit.f90:1180
mod_urban_vars::urban_tr
real(rp), dimension(:,:), allocatable, public urban_tr
Definition: mod_urban_vars.f90:60
mod_realinput
module REAL input
Definition: mod_realinput.f90:13
scale_atmos_profile
module ATMOSPHERE / Typical vertical profile
Definition: scale_atmos_sub_profile.F90:14
mod_mkinit::i_seabreeze
integer, parameter, public i_seabreeze
Definition: mod_mkinit.f90:134
scale_grid_index::ia
integer, public ia
of whole cells: x, local, with HALO
Definition: scale_grid_index.F90:55
mod_urban_vars::urban_tgl
real(rp), dimension(:,:,:), allocatable, public urban_tgl
Definition: mod_urban_vars.f90:68
scale_atmos_hydrometeor
module Hydrometeor
Definition: scale_atmos_sub_hydrometeor.F90:14
scale_grid_index::jag
integer, public jag
of computational grids
Definition: scale_grid_index.F90:82
scale_stdio::io_get_available_fid
integer function, public io_get_available_fid()
search & get available file ID
Definition: scale_stdio.F90:313
mod_mkinit::rect_setup
subroutine rect_setup
Bubble.
Definition: mod_mkinit.f90:626
scale_grid_real
module GRID (real space)
Definition: scale_grid_real.F90:14
scale_landuse
module LANDUSE
Definition: scale_landuse.F90:11
scale_grid_index::ka
integer, public ka
of whole cells: z, local, with HALO
Definition: scale_grid_index.F90:54
mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_sw_dn
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_sflx_sw_dn
Definition: mod_atmos_phy_rd_vars.f90:61
mod_ocean_driver
module OCEAN driver
Definition: mod_ocean_driver.f90:10
mod_mkinit::i_warmbubble
integer, parameter, public i_warmbubble
Definition: mod_mkinit.f90:118
mod_ocean_driver::ocean_surface_set
subroutine, public ocean_surface_set(countup)
Put surface boundary to other model.
Definition: mod_ocean_driver.f90:262
scale_const::const_pre00
real(rp), public const_pre00
pressure reference [Pa]
Definition: scale_const.F90:90
mod_mkinit::read_sounding
subroutine read_sounding(DENS, VELX, VELY, POTT, QV)
Read sounding data from file.
Definition: mod_mkinit.f90:1225
scale_grid::grid_fz
real(rp), dimension(:), allocatable, public grid_fz
face coordinate [m]: z, local=global
Definition: scale_grid_cartesian.F90:66
mod_atmos_driver
module ATMOSPHERE driver
Definition: mod_atmos_driver.f90:18
mod_mkinit::i_mountainwave
integer, parameter, public i_mountainwave
Definition: mod_mkinit.f90:116
scale_landuse::landuse_write
subroutine, public landuse_write
Write landuse data.
Definition: scale_landuse.F90:271
scale_grid_index::jhalo
integer, public jhalo
of halo cells: y
Definition: scale_grid_index.F90:52
mod_land_vars::land_temp
real(rp), dimension(:,:,:), allocatable, public land_temp
temperature of each soil layer [K]
Definition: mod_land_vars.f90:62
scale_comm
module COMMUNICATION
Definition: scale_comm.F90:23
scale_atmos_hydrostatic
module ATMOSPHERE / Hydrostatic barance
Definition: scale_atmos_sub_hydrostatic.F90:11
scale_const::const_grav
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:48
scale_grid_index::js
integer, public js
start point of inner domain: y, local
Definition: scale_grid_index.F90:62
scale_grid_index::iag
integer, public iag
of computational grids
Definition: scale_grid_index.F90:81
mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_lw_dn
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_lw_dn
Definition: mod_atmos_phy_rd_vars.f90:64
mod_urban_vars::urban_roff
real(rp), dimension(:,:), allocatable, public urban_roff
Definition: mod_urban_vars.f90:72
scale_process
module PROCESS
Definition: scale_process.F90:14
mod_atmos_phy_ae_vars::atmos_phy_ae_ccn
real(rp), dimension(:,:,:), allocatable, public atmos_phy_ae_ccn
Definition: mod_atmos_phy_ae_vars.f90:57
mod_realinput::realinput_atmos
subroutine, public realinput_atmos(flg_intrp)
Definition: mod_realinput.f90:162
mod_mkinit::i_real
integer, parameter, public i_real
Definition: mod_mkinit.f90:139
mod_land_vars::land_sfc_albedo
real(rp), dimension(:,:,:), allocatable, public land_sfc_albedo
land surface albedo (0-1)
Definition: mod_land_vars.f90:65
mod_mkinit::i_heatisland
integer, parameter, public i_heatisland
Definition: mod_mkinit.f90:135
mod_mkinit::i_triplecouple
integer, parameter, public i_triplecouple
Definition: mod_mkinit.f90:131
mod_land_vars
module LAND Variables
Definition: mod_land_vars.f90:11
scale_atmos_phy_mp_suzuki10::nccn
integer, public nccn
Definition: scale_atmos_phy_mp_suzuki10.F90:111
mod_admin_restart
module administrator for restart
Definition: mod_admin_restart.f90:10
scale_const
module CONSTANT
Definition: scale_const.F90:14
mod_urban_vars::urban_sfc_albedo
real(rp), dimension(:,:,:), allocatable, public urban_sfc_albedo
Definition: mod_urban_vars.f90:91
scale_atmos_phy_tb
module ATMOSPHERE / Physics Turbulence
Definition: scale_atmos_phy_tb.F90:15
mod_mkinit::i_planestate
integer, parameter, public i_planestate
Definition: mod_mkinit.f90:108
mod_atmos_admin::atmos_phy_mp_type
character(len=h_short), public atmos_phy_mp_type
Definition: mod_atmos_admin.f90:35
mod_mkinit::i_khwave
integer, parameter, public i_khwave
Definition: mod_mkinit.f90:114
scale_grid_index::ks
integer, public ks
start point of inner domain: z, local
Definition: scale_grid_index.F90:58
scale_grid::grid_cx
real(rp), dimension(:), allocatable, public grid_cx
center coordinate [m]: x, local
Definition: scale_grid_cartesian.F90:57
mod_mkinit::i_squallline
integer, parameter, public i_squallline
Definition: mod_mkinit.f90:120
mod_urban_vars::urban_tc
real(rp), dimension(:,:), allocatable, public urban_tc
Definition: mod_urban_vars.f90:63
mod_atmos_vars::momy
real(rp), dimension(:,:,:), allocatable, target, public momy
Definition: mod_atmos_vars.f90:75
mod_admin_time::time_dourban_restart
logical, public time_dourban_restart
execute urban restart output in this step?
Definition: mod_admin_time.f90:60
scale_grid
module GRID (cartesian)
Definition: scale_grid_cartesian.F90:17
mod_urban_vars::urban_rainr
real(rp), dimension(:,:), allocatable, public urban_rainr
Definition: mod_urban_vars.f90:69
scale_random::random_get
subroutine, public random_get(var)
Get random number.
Definition: scale_random.F90:136
scale_const::const_undef8
real(dp), parameter, public const_undef8
undefined value (REAL8)
Definition: scale_const.F90:42
mod_land_vars::land_water
real(rp), dimension(:,:,:), allocatable, public land_water
moisture of each soil layer [m3/m3]
Definition: mod_land_vars.f90:63
mod_atmos_phy_rd_vars::atmos_phy_rd_sflx_downall
real(rp), dimension(:,:,:,:), allocatable, public atmos_phy_rd_sflx_downall
Definition: mod_atmos_phy_rd_vars.f90:68
mod_mkinit::mkinit_type
integer, public mkinit_type
Definition: mod_mkinit.f90:105
scale_rm_process
module RM PROCESS
Definition: scale_rm_process.F90:10
mod_admin_restart::admin_restart_write
subroutine, public admin_restart_write
Write data to restart files.
Definition: mod_admin_restart.f90:293
mod_atmos_phy_mp_vars::atmos_phy_mp_sflx_snow
real(rp), dimension(:,:), allocatable, public atmos_phy_mp_sflx_snow
Definition: mod_atmos_phy_mp_vars.f90:65
mod_mkinit::i_rico
integer, parameter, public i_rico
Definition: mod_mkinit.f90:124
scale_grid::grid_domain_center_y
real(rp), public grid_domain_center_y
center position of global domain [m]: y
Definition: scale_grid_cartesian.F90:54
scale_prof
module profiler
Definition: scale_prof.F90:10
scale_grid_index::ie
integer, public ie
end point of inner domain: x, local
Definition: scale_grid_index.F90:61
mod_mkinit::i_ignore
integer, parameter, public i_ignore
Definition: mod_mkinit.f90:106
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:36
mod_mkinit::i_lambwave
integer, parameter, public i_lambwave
Definition: mod_mkinit.f90:112
mod_mkinit::i_supercell
integer, parameter, public i_supercell
Definition: mod_mkinit.f90:119
scale_atmos_hydrometeor::i_qv
integer, public i_qv
Definition: scale_atmos_sub_hydrometeor.F90:68
mod_ocean_vars::ocean_sfc_z0m
real(rp), dimension(:,:), allocatable, public ocean_sfc_z0m
ocean surface roughness length for momentum [m]
Definition: mod_ocean_vars.f90:62
mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_sw_dn
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_sw_dn
Definition: mod_atmos_phy_rd_vars.f90:66
mod_mkinit::i_bomex
integer, parameter, public i_bomex
Definition: mod_mkinit.f90:147
scale_const::const_i_sw
integer, public const_i_sw
short-wave radiation index
Definition: scale_const.F90:96
scale_precision
module PRECISION
Definition: scale_precision.F90:13
scale_atmos_phy_ae_kajino13
module ATMOSPHERE / Physics Aerosol Microphysics
Definition: scale_atmos_phy_ae_kajino13.F90:15
mod_mkinit::i_dycoms2_rf01
integer, parameter, public i_dycoms2_rf01
Definition: mod_mkinit.f90:122
mod_mkinit::i_interporation
integer, parameter, public i_interporation
Definition: mod_mkinit.f90:126
mod_admin_time::time_doland_restart
logical, public time_doland_restart
execute land restart output in this step?
Definition: mod_admin_time.f90:58
scale_const::const_pi
real(rp), public const_pi
pi
Definition: scale_const.F90:34
mod_urban_vars::urban_tg
real(rp), dimension(:,:), allocatable, public urban_tg
Definition: mod_urban_vars.f90:62
mod_admin_time
module ADMIN TIME
Definition: mod_admin_time.f90:10
mod_ocean_vars::ocean_sfc_temp
real(rp), dimension(:,:), allocatable, public ocean_sfc_temp
ocean surface skin temperature [K]
Definition: mod_ocean_vars.f90:60
mod_urban_vars::urban_sfc_temp
real(rp), dimension(:,:), allocatable, public urban_sfc_temp
Definition: mod_urban_vars.f90:90
mod_urban_vars::urban_trl
real(rp), dimension(:,:,:), allocatable, public urban_trl
Definition: mod_urban_vars.f90:66
scale_stdio::io_fid_conf
integer, public io_fid_conf
Config file ID.
Definition: scale_stdio.F90:55
mod_mkinit::i_warmbubbleaero
integer, parameter, public i_warmbubbleaero
Definition: mod_mkinit.f90:143
mod_mkinit::i_wk1982
integer, parameter, public i_wk1982
Definition: mod_mkinit.f90:121
mod_mkinit
module INITIAL
Definition: mod_mkinit.f90:31
mod_mkinit::i_bubblecouple
integer, parameter, public i_bubblecouple
Definition: mod_mkinit.f90:132
scale_random
module RANDOM
Definition: scale_random.F90:13
mod_mkinit::i_tracerbubble
integer, parameter, public i_tracerbubble
Definition: mod_mkinit.f90:109
scale_stdio::io_fid_log
integer, public io_fid_log
Log file ID.
Definition: scale_stdio.F90:56
mod_mkinit::i_coldbubble
integer, parameter, public i_coldbubble
Definition: mod_mkinit.f90:110
mod_realinput::realinput_surface
subroutine, public realinput_surface
Definition: mod_realinput.f90:354
mod_mkinit::i_cavityflow
integer, parameter, public i_cavityflow
Definition: mod_mkinit.f90:145
scale_atmos_phy_mp_suzuki10
module Spectran Bin Microphysics
Definition: scale_atmos_phy_mp_suzuki10.F90:31
mod_urban_vars::urban_tbl
real(rp), dimension(:,:,:), allocatable, public urban_tbl
Definition: mod_urban_vars.f90:67
mod_atmos_phy_ae_vars
module ATMOSPHERE / Physics Aerosol Microphysics
Definition: mod_atmos_phy_ae_vars.f90:14
mod_mkinit::i_barocwave
integer, parameter, public i_barocwave
Definition: mod_mkinit.f90:146
mod_land_driver
module LAND driver
Definition: mod_land_driver.f90:11
scale_landuse::landuse_frac_land
real(rp), dimension(:,:), allocatable, public landuse_frac_land
land fraction
Definition: scale_landuse.F90:39
mod_mkinit::mkinit_setup
subroutine, public mkinit_setup
Setup.
Definition: mod_mkinit.f90:228
mod_mkinit::i_dycoms2_rf02
integer, parameter, public i_dycoms2_rf02
Definition: mod_mkinit.f90:123
mod_land_driver::land_surface_set
subroutine, public land_surface_set(countup)
Put surface boundary to other model.
Definition: mod_land_driver.f90:262
mod_ocean_vars::ocean_sfc_z0h
real(rp), dimension(:,:), allocatable, public ocean_sfc_z0h
ocean surface roughness length for heat [m]
Definition: mod_ocean_vars.f90:63
mod_ocean_vars
module OCEAN Variables
Definition: mod_ocean_vars.f90:11
scale_stdio::io_fid_nml
integer, public io_fid_nml
Log file ID (only for output namelist)
Definition: scale_stdio.F90:57
scale_const::const_pstd
real(rp), public const_pstd
standard pressure [Pa]
Definition: scale_const.F90:89
scale_grid::grid_fyg
real(rp), dimension(:), allocatable, public grid_fyg
face coordinate [m]: y, global
Definition: scale_grid_cartesian.F90:84
scale_grid::grid_cy
real(rp), dimension(:), allocatable, public grid_cy
center coordinate [m]: y, local
Definition: scale_grid_cartesian.F90:58
mod_urban_vars::urban_rainb
real(rp), dimension(:,:), allocatable, public urban_rainb
Definition: mod_urban_vars.f90:70
mod_atmos_phy_rd_vars::atmos_phy_rd_toaflx_lw_up
real(rp), dimension(:,:), allocatable, public atmos_phy_rd_toaflx_lw_up
Definition: mod_atmos_phy_rd_vars.f90:63
mod_mkinit::urban_setup
subroutine urban_setup
Urban setup.
Definition: mod_mkinit.f90:1074
mod_mkinit::flux_setup
subroutine flux_setup
flux setup
Definition: mod_mkinit.f90:887
mod_atmos_vars::qtrc
real(rp), dimension(:,:,:,:), allocatable, target, public qtrc
Definition: mod_atmos_vars.f90:77
scale_grid_index::ihalo
integer, public ihalo
of halo cells: x
Definition: scale_grid_index.F90:51
gtool_file::filegetshape
subroutine, public filegetshape(dims, basename, varname, myrank, single, error)
Definition: gtool_file.f90:1529
scale_grid_index::ja
integer, public ja
of whole cells: y, local, with HALO
Definition: scale_grid_index.F90:56
scale_grid::grid_fy
real(rp), dimension(:), allocatable, public grid_fy
face coordinate [m]: y, local
Definition: scale_grid_cartesian.F90:68
Generated by   doxygen 1.8.14