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