SCALE-RM
scale_atmos_phy_rd_mstrnx.F90
Go to the documentation of this file.
1 !-------------------------------------------------------------------------------
13 !-------------------------------------------------------------------------------
14 #include "scalelib.h"
15 module scale_atmos_phy_rd_mstrnx
16  !-----------------------------------------------------------------------------
17  !
18  !++ used modules
19  !
20  use scale_precision
21  use scale_io
22  use scale_prof
23  use scale_cpl_sfc_index
24 
25  use scale_atmos_phy_rd_common, only: &
26  i_sw, &
27  i_lw, &
28  i_dn, &
29  i_up
30 
31  use scale_atmos_hydrometeor, only: &
32  n_hyd
33  use scale_atmos_aerosol, only: &
34  n_ae
35  !-----------------------------------------------------------------------------
36  implicit none
37  private
38  !-----------------------------------------------------------------------------
39  !
40  !++ Public procedure
41  !
42  public :: atmos_phy_rd_mstrnx_setup
43  public :: atmos_phy_rd_mstrnx_finalize
44  public :: atmos_phy_rd_mstrnx_flux
45 
46  !-----------------------------------------------------------------------------
47  !
48  !++ Public parameters & variables
49  !
50  !-----------------------------------------------------------------------------
51  !
52  !++ Private procedure
53  !
54  private :: rd_mstrn_setup
55  private :: rd_mstrn_dtrn3
56  private :: rd_mstrn_two_stream
57 
58  !-----------------------------------------------------------------------------
59  !
60  !++ Private parameters & variables
61  !
62  real(RP), private, parameter :: RD_cosSZA_min = 0.017_rp ! minimum SZA (>89.0)
63  real(RP), private, parameter :: RD_EPS = 1.e-4_rp
64 
65  real(RP), private :: RD_TOA = 100.0_rp
66  integer, private :: RD_KADD = 10
67  integer, private, parameter :: RD_naero = n_hyd + n_ae ! # of cloud/aerosol species
68  integer, private, parameter :: RD_hydro_str = 1 ! start index for cloud
69  integer, private, parameter :: RD_hydro_end = n_hyd ! end index for cloud
70  integer, private, parameter :: RD_aero_str = n_hyd + 1 ! start index for aerosol
71  integer, private, parameter :: RD_aero_end = n_hyd + n_ae ! end index for aerosol
72 
73  integer, private :: RD_KMAX ! # of computational cells: z for radiation scheme
74 
75  real(RP), private, allocatable :: RD_zh (:) ! altitude at the interface [km]
76  real(RP), private, allocatable :: RD_z (:) ! altitude at the center [km]
77  real(RP), private, allocatable :: RD_rhodz (:) ! density * delta z [kg/m2]
78  real(RP), private, allocatable :: RD_pres (:) ! pressure at the center [hPa]
79  real(RP), private, allocatable :: RD_presh (:) ! pressure at the interface [hPa]
80  real(RP), private, allocatable :: RD_temp (:) ! temperature at the center [K]
81  real(RP), private, allocatable :: RD_temph (:) ! temperature at the interface [K]
82  real(RP), private, allocatable :: RD_gas (:,:) ! gas species volume mixing ratio [ppmv]
83  real(RP), private, allocatable :: RD_cfc (:,:) ! CFCs volume mixing ratio [ppmv]
84  real(RP), private, allocatable :: RD_aerosol_conc(:,:) ! cloud/aerosol volume mixing ratio [ppmv]
85  real(RP), private, allocatable :: RD_aerosol_radi(:,:) ! cloud/aerosol effective radius [cm]
86  real(RP), private, allocatable :: RD_cldfrac (:) ! cloud fraction (0-1)
87 
88  integer, private :: I_MPAE2RD(RD_naero) ! look-up table between input aerosol category and MSTRN particle type
89 
90  character(len=H_LONG), private :: MSTRN_GASPARA_INPUTFILE = 'PARAG.29'
91  character(len=H_LONG), private :: MSTRN_AEROPARA_INPUTFILE = 'PARAPC.29'
92  character(len=H_LONG), private :: MSTRN_HYGROPARA_INPUTFILE = 'VARDATA.RM29'
93 
94  integer, private :: MSTRN_nband = 29
95 
96  integer, private, parameter :: MSTRN_nstream = 1
97  integer, private, parameter :: MSTRN_ch_limit = 10
98  integer, private, parameter :: MSTRN_nflag = 7 ! # of optical properties flag
99  integer, private, parameter :: MSTRN_nfitP = 26 ! # of fitting point for log(pressure)
100  integer, private, parameter :: MSTRN_nfitT = 3 ! # of fitting point for temperature
101 
102  integer, private, parameter :: MSTRN_ngas = 7
103  ! 1: H2O
104  ! 2: CO2
105  ! 3: O3
106  ! 4: N2O
107  ! 5: CO
108  ! 6: CH4
109  ! 7: O2
110  integer, private, parameter :: MSTRN_ncfc = 28
111  ! 1: CFC-11
112  ! 2: CFC-12
113  ! 3: CFC-13
114  ! 4: CFC-14
115  ! 5: CFC-113
116  ! 6: CFC-114
117  ! 7: CFC-115
118  ! 8: HCFC-21
119  ! 9: HCFC-22
120  ! 10: HCFC-123
121  ! 11: HCFC-124
122  ! 12: HCFC-141b
123  ! 13: HCFC-142b
124  ! 14: HCFC-225ca
125  ! 15: HCFC-225cb
126  ! 16: HFC-32
127  ! 17: HFC-125
128  ! 18: HFC-134
129  ! 19: HFC-134a
130  ! 20: HFC-143a
131  ! 21: HFC-152a
132  ! 22: SF6
133  ! 23: ClONO2
134  ! 24: CCl4
135  ! 25: N2O5
136  ! 26: C2F6
137  ! 27: HNO4
138  integer, private :: MSTRN_nptype = 9
139  ! 1: water cloud
140  ! 2: ice cloud
141  ! 3: Soil dust
142  ! 4: Carbonacerous (BC/OC=0.3)
143  ! 5: Carbonacerous (BC/OC=0.15)
144  ! 6: Carbonacerous (BC/OC=0.)
145  ! 7: Black carbon
146  ! 8: Sulfate
147  ! 9: Sea salt
148 
149  integer, private, parameter :: MSTRN_nsfc = 7
150  ! 1: ocean
151  ! 2: wet land
152  ! 3: dry land
153  ! 4: low plants
154  ! 5: forest
155  ! 6: snow
156  ! 7: ice
157  integer, private, parameter :: MSTRN_nfitPLK = 5
158  integer, private, parameter :: MSTRN_nplkord = 3
159  integer, private, parameter :: MSTRN_nmoment = 6
160  integer, private :: MSTRN_nradius = 8
161  integer, private, parameter :: MSTRN_ncloud = 2
162 
163  logical, private :: ATMOS_PHY_RD_MSTRN_ONLY_QCI = .false.
164  logical, private :: ATMOS_PHY_RD_MSTRN_ONLY_TROPOCLOUD = .false.
165  logical, private :: ATMOS_PHY_RD_MSTRN_USE_AERO = .false.
166 
167 
168 
169  real(RP), private, allocatable :: waveh (:) ! wavenumbers at band boundary [1/cm]
170 
171  real(RP), private, allocatable :: logfitP (:) ! fitting point for log10(pressure)
172  real(RP), private, allocatable :: fitT (:) ! fitting point for temperature
173  real(RP), private, allocatable :: logfitT (:) ! fitting point for log10(temperature)
174  integer, private, allocatable :: iflgb (:,:) ! optical properties flag in each band
175  integer, private, allocatable :: nch (:) ! number of subintervals in each band
176  real(RP), private, allocatable :: wgtch (:,:) ! weights of subintervals in each band
177  integer, private, allocatable :: ngasabs (:) ! number of absorbers(gas) in each band
178  integer, private, allocatable :: igasabs (:,:) ! index of absorbers(gas) in each band
179 
180  real(RP), private, allocatable :: akd (:,:,:,:,:) ! absorption coefficient table
181  real(RP), private, allocatable :: skd (:,:,:,:) ! absorption coefficient table for H2O self broadening
182  real(RP), private, allocatable :: acfc_pow(:,:) ! 10 to the power of absorption coefficient table for CFC
183 
184  real(RP), private, allocatable :: fitPLK (:,:) ! fitting point for planck function
185  real(RP), private, allocatable :: fsol (:) ! solar insolation in each band
186  real(RP), private :: fsol_tot ! total solar insolation
187  real(RP), private, allocatable :: sfc (:,:) ! surface condition in each band
188  real(RP), private, allocatable :: rayleigh(:) ! rayleigh scattering in each band
189  real(RP), private, allocatable :: qmol (:,:) ! moments for rayleigh scattering phase function
190  real(RP), private, allocatable :: q (:,:,:,:) ! moments for aerosol scattering phase function
191 
192  integer, private, allocatable :: hygro_flag(:) ! flag for hygroscopic enlargement
193  real(RP), private, allocatable :: radmode (:,:) ! radius mode for hygroscopic parameter
194 
195  integer, private, allocatable :: ptype_nradius(:) ! number limit of radius mode for cloud / aerosol particle type
196 
197 
198  ! index for optical flag iflgb
199  integer, private, parameter :: I_SWLW = 4
200  integer, private, parameter :: I_H2O_continuum = 5
201  integer, private, parameter :: I_CFC_continuum = 7
202  ! for cloud type
203  integer, private, parameter :: I_ClearSky = 1
204  integer, private, parameter :: I_Cloud = 2
205 
206  ! pre-calc
207  real(RP), private :: RHO_std ! rho(0C,1atm) [kg/m3]
208 
209  real(RP), private :: M(2) ! discrete quadrature mu for two-stream approximation
210  real(RP), private :: W(2) ! discrete quadrature w for two-stream approximation
211  real(RP), private :: Wmns(2), Wpls(2) ! W-, W+
212  real(RP), private :: Wbar(2), Wscale(2)
213  !$acc declare create(M, W, Wmns, Wpls, Wscale)
214 
215  !-----------------------------------------------------------------------------
216 contains
217  !-----------------------------------------------------------------------------
219  subroutine atmos_phy_rd_mstrnx_setup( &
220  KA, KS, KE, &
221  CZ, FZ )
222  use scale_prc, only: &
223  prc_abort
224  use scale_time, only: &
225  time_nowdate
226  use scale_atmos_grid_cartesc_real, only: &
227  atmos_grid_cartesc_real_basepoint_lat
228  use scale_atmos_phy_rd_profile, only: &
229  rd_profile_setup => atmos_phy_rd_profile_setup, &
230  rd_profile_setup_zgrid => atmos_phy_rd_profile_setup_zgrid, &
231  rd_profile_read => atmos_phy_rd_profile_read
232  use scale_atmos_hydrometeor, only: &
233  n_hyd
234  use scale_atmos_aerosol, only: &
235  n_ae
236  implicit none
237  integer, intent(in) :: ka, ks, ke
238 
239  real(rp), intent(in) :: cz(ka)
240  real(rp), intent(in) :: fz(0:ka)
241 
242  real(rp) :: atmos_phy_rd_mstrn_toa
243  integer :: atmos_phy_rd_mstrn_kadd
244  character(len=H_LONG) :: atmos_phy_rd_mstrn_gaspara_in_filename
245  character(len=H_LONG) :: atmos_phy_rd_mstrn_aeropara_in_filename
246  character(len=H_LONG) :: atmos_phy_rd_mstrn_hygropara_in_filename
247  integer :: atmos_phy_rd_mstrn_nband
248  integer :: atmos_phy_rd_mstrn_nptype
249  integer :: atmos_phy_rd_mstrn_nradius
250  integer :: atmos_phy_rd_mstrn_nradius_cloud
251  integer :: atmos_phy_rd_mstrn_nradius_aero
252 
253  namelist / param_atmos_phy_rd_mstrn / &
254  atmos_phy_rd_mstrn_toa, &
255  atmos_phy_rd_mstrn_kadd, &
256  atmos_phy_rd_mstrn_gaspara_in_filename, &
257  atmos_phy_rd_mstrn_aeropara_in_filename, &
258  atmos_phy_rd_mstrn_hygropara_in_filename, &
259  atmos_phy_rd_mstrn_nband, &
260  atmos_phy_rd_mstrn_nptype, &
261  atmos_phy_rd_mstrn_nradius, &
262  atmos_phy_rd_mstrn_nradius_cloud, &
263  atmos_phy_rd_mstrn_nradius_aero, &
264  atmos_phy_rd_mstrn_only_qci, &
265  atmos_phy_rd_mstrn_only_tropocloud, &
266  atmos_phy_rd_mstrn_use_aero
267 
268  integer :: kmax
269  integer :: ngas, ncfc
270  integer :: ierr
271  !---------------------------------------------------------------------------
272 
273  log_newline
274  log_info("ATMOS_PHY_RD_mstrnx_setup",*) 'Setup'
275  log_info("ATMOS_PHY_RD_mstrnx_setup",*) 'Sekiguchi and Nakajima (2008) mstrnX radiation process'
276 
277  !--- read namelist
278  atmos_phy_rd_mstrn_toa = rd_toa
279  atmos_phy_rd_mstrn_kadd = rd_kadd
280  atmos_phy_rd_mstrn_gaspara_in_filename = mstrn_gaspara_inputfile
281  atmos_phy_rd_mstrn_aeropara_in_filename = mstrn_aeropara_inputfile
282  atmos_phy_rd_mstrn_hygropara_in_filename = mstrn_hygropara_inputfile
283  atmos_phy_rd_mstrn_nband = mstrn_nband
284  atmos_phy_rd_mstrn_nptype = mstrn_nptype
285  atmos_phy_rd_mstrn_nradius = mstrn_nradius
286  atmos_phy_rd_mstrn_nradius_cloud = -1
287  atmos_phy_rd_mstrn_nradius_aero = -1
288 
289  rewind(io_fid_conf)
290  read(io_fid_conf,nml=param_atmos_phy_rd_mstrn,iostat=ierr)
291  if( ierr < 0 ) then !--- missing
292  log_info("ATMOS_PHY_RD_mstrnx_setup",*) 'Not found namelist. Default used.'
293  elseif( ierr > 0 ) then !--- fatal error
294  log_error("ATMOS_PHY_RD_mstrnx_setup",*) 'Not appropriate names in namelist PARAM_ATMOS_PHY_RD_MSTRN. Check!'
295  call prc_abort
296  endif
297  log_nml(param_atmos_phy_rd_mstrn)
298 
299  rd_toa = atmos_phy_rd_mstrn_toa
300  rd_kadd = atmos_phy_rd_mstrn_kadd
301  mstrn_gaspara_inputfile = atmos_phy_rd_mstrn_gaspara_in_filename
302  mstrn_aeropara_inputfile = atmos_phy_rd_mstrn_aeropara_in_filename
303  mstrn_hygropara_inputfile = atmos_phy_rd_mstrn_hygropara_in_filename
304  mstrn_nband = atmos_phy_rd_mstrn_nband
305  mstrn_nptype = atmos_phy_rd_mstrn_nptype
306 
307  if ( atmos_phy_rd_mstrn_nradius_cloud < 0 ) then ! set default
308  atmos_phy_rd_mstrn_nradius_cloud = atmos_phy_rd_mstrn_nradius
309  endif
310  if ( atmos_phy_rd_mstrn_nradius_aero < 0 ) then ! set default
311  atmos_phy_rd_mstrn_nradius_aero = atmos_phy_rd_mstrn_nradius
312  endif
313 
314  mstrn_nradius = max( atmos_phy_rd_mstrn_nradius, &
315  atmos_phy_rd_mstrn_nradius_cloud, &
316  atmos_phy_rd_mstrn_nradius_aero )
317 
318  !--- setup MP/AE - RD(particle) coupling parameter
319  allocate( ptype_nradius(mstrn_nptype) )
320 
321  if ( mstrn_nptype == 9 ) then
322 
323  i_mpae2rd( 1) = 1 ! cloud water Qc -> Water
324  i_mpae2rd( 2) = 1 ! rain Qr -> Water
325  i_mpae2rd( 3) = 2 ! cloud ice Qi -> Ice
326  i_mpae2rd( 4) = 2 ! snow Qs -> Ice
327  i_mpae2rd( 5) = 2 ! graupel Qg -> Ice
328  i_mpae2rd( 6) = 2 ! hail Qh -> Ice
329  i_mpae2rd( 7) = 3 ! aerosol type1 -> Soil dust
330  i_mpae2rd( 8) = 4 ! aerosol type2 -> Carbonacerous (BC/OC=0.3)
331  i_mpae2rd( 9) = 5 ! aerosol type3 -> Carbonacerous (BC/OC=0.15)
332  i_mpae2rd(10) = 6 ! aerosol type4 -> Carbonacerous (BC/OC=0.)
333  i_mpae2rd(11) = 7 ! aerosol type5 -> Black carbon
334  i_mpae2rd(12) = 8 ! aerosol type6 -> Sulfate
335  i_mpae2rd(13) = 9 ! aerosol type7 -> Sea salt
336 
337  ptype_nradius( 1) = atmos_phy_rd_mstrn_nradius_cloud
338  ptype_nradius( 2) = atmos_phy_rd_mstrn_nradius_cloud
339  ptype_nradius( 3) = atmos_phy_rd_mstrn_nradius_aero
340  ptype_nradius( 4) = atmos_phy_rd_mstrn_nradius_aero
341  ptype_nradius( 5) = atmos_phy_rd_mstrn_nradius_aero
342  ptype_nradius( 6) = atmos_phy_rd_mstrn_nradius_aero
343  ptype_nradius( 7) = atmos_phy_rd_mstrn_nradius_aero
344  ptype_nradius( 8) = atmos_phy_rd_mstrn_nradius_aero
345  ptype_nradius( 9) = atmos_phy_rd_mstrn_nradius_aero
346 
347  elseif( mstrn_nptype == 12 ) then
348 
349  i_mpae2rd( 1) = 1 ! cloud water Qc -> CLOUD
350  i_mpae2rd( 2) = 2 ! rain Qr -> RAIN
351  i_mpae2rd( 3) = 3 ! cloud ice Qi -> ICE
352  i_mpae2rd( 4) = 4 ! snow Qs -> SNOW
353  i_mpae2rd( 5) = 5 ! graupel Qg -> GRAUPEL
354  i_mpae2rd( 6) = 5 ! hail Qh -> GRAUPEL
355  i_mpae2rd( 7) = 6 ! aerosol type1 -> Soil dust
356  i_mpae2rd( 8) = 7 ! aerosol type2 -> Carbonacerous (BC/OC=0.3)
357  i_mpae2rd( 9) = 8 ! aerosol type3 -> Carbonacerous (BC/OC=0.15)
358  i_mpae2rd(10) = 9 ! aerosol type4 -> Carbonacerous (BC/OC=0.)
359  i_mpae2rd(11) = 10 ! aerosol type5 -> Black carbon
360  i_mpae2rd(12) = 11 ! aerosol type6 -> Sulfate
361  i_mpae2rd(13) = 12 ! aerosol type7 -> Sea salt
362 
363  ptype_nradius( 1) = atmos_phy_rd_mstrn_nradius_cloud
364  ptype_nradius( 2) = atmos_phy_rd_mstrn_nradius_cloud
365  ptype_nradius( 3) = atmos_phy_rd_mstrn_nradius_cloud
366  ptype_nradius( 4) = atmos_phy_rd_mstrn_nradius_cloud
367  ptype_nradius( 5) = atmos_phy_rd_mstrn_nradius_cloud
368  ptype_nradius( 6) = atmos_phy_rd_mstrn_nradius_aero
369  ptype_nradius( 7) = atmos_phy_rd_mstrn_nradius_aero
370  ptype_nradius( 8) = atmos_phy_rd_mstrn_nradius_aero
371  ptype_nradius( 9) = atmos_phy_rd_mstrn_nradius_aero
372  ptype_nradius(10) = atmos_phy_rd_mstrn_nradius_aero
373  ptype_nradius(11) = atmos_phy_rd_mstrn_nradius_aero
374  ptype_nradius(12) = atmos_phy_rd_mstrn_nradius_aero
375 
376  endif
377  !$acc enter data copyin(I_MPAE2RD, ptype_nradius)
378 
379  !--- setup MSTRN parameter
380  call rd_mstrn_setup( ngas, & ! [OUT]
381  ncfc ) ! [OUT]
382 
383  !--- setup climatological profile
384  call rd_profile_setup
385 
386  kmax = ke - ks + 1
387  rd_kmax = kmax + rd_kadd
388 
389  !--- allocate arrays
390  ! input
391  allocate( rd_zh(rd_kmax+1) )
392  allocate( rd_z(rd_kmax ) )
393 
394  allocate( rd_rhodz(rd_kmax ) )
395  allocate( rd_pres(rd_kmax ) )
396  allocate( rd_presh(rd_kmax+1) )
397  allocate( rd_temp(rd_kmax ) )
398  allocate( rd_temph(rd_kmax+1) )
399 
400  allocate( rd_gas(rd_kmax,ngas ) )
401  allocate( rd_cfc(rd_kmax,ncfc ) )
402  allocate( rd_aerosol_conc(rd_kmax,rd_naero) )
403  allocate( rd_aerosol_radi(rd_kmax,rd_naero) )
404  allocate( rd_cldfrac(rd_kmax ) )
405 
406  !--- setup vartical grid for radiation (larger TOA than Model domain)
407  call rd_profile_setup_zgrid( &
408  ka, ks, ke, &
409  rd_kmax, rd_kadd, & ! [IN]
410  rd_toa, cz(:), fz(:), & ! [IN]
411  rd_zh(:), rd_z(:) ) ! [OUT]
412 
413  !--- read climatological profile
414  call rd_profile_read( rd_kmax, & ! [IN]
415  ngas, & ! [IN]
416  ncfc, & ! [IN]
417  rd_naero, & ! [IN]
418  atmos_grid_cartesc_real_basepoint_lat, & ! [IN]
419  time_nowdate(:), & ! [IN]
420  rd_zh(:), & ! [IN]
421  rd_z(:), & ! [IN]
422  rd_rhodz(:), & ! [OUT]
423  rd_pres(:), & ! [OUT]
424  rd_presh(:), & ! [OUT]
425  rd_temp(:), & ! [OUT]
426  rd_temph(:), & ! [OUT]
427  rd_gas(:,:), & ! [OUT]
428  rd_cfc(:,:), & ! [OUT]
429  rd_aerosol_conc(:,:), & ! [OUT]
430  rd_aerosol_radi(:,:), & ! [OUT]
431  rd_cldfrac(:) ) ! [OUT]
432 
433  !$acc enter data copyin(RD_rhodz,RD_pres,RD_presh,RD_temp,RD_temph,RD_gas,RD_cfc,RD_aerosol_conc,RD_aerosol_radi,RD_cldfrac)
434 
435  return
436  end subroutine atmos_phy_rd_mstrnx_setup
437 
438  !-----------------------------------------------------------------------------
440  subroutine atmos_phy_rd_mstrnx_finalize
441 
442  !$acc exit data delete(I_MPAE2RD, ptype_nradius)
443 
444  !$acc exit data &
445  !$acc delete(waveh, &
446  !$acc wgtch, fitPLK, logfitP, logfitT, fitT, &
447  !$acc radmode, ngasabs, igasabs, &
448  !$acc fsol, q, qmol, rayleigh, acfc_pow, nch, AKD, SKD )
449 
450  !$acc exit data delete(RD_rhodz,RD_pres,RD_presh,RD_temp,RD_temph,RD_gas,RD_cfc,RD_aerosol_conc,RD_aerosol_radi,RD_cldfrac)
451 
452  deallocate( ptype_nradius )
453  deallocate( rd_zh )
454  deallocate( rd_z )
455 
456  deallocate( rd_rhodz )
457  deallocate( rd_pres )
458  deallocate( rd_presh )
459  deallocate( rd_temp )
460  deallocate( rd_temph )
461 
462  deallocate( rd_gas )
463  deallocate( rd_cfc )
464  deallocate( rd_aerosol_conc )
465  deallocate( rd_aerosol_radi )
466  deallocate( rd_cldfrac )
467 
468  deallocate( waveh )
469  deallocate( logfitp )
470  deallocate( fitt )
471  deallocate( logfitt )
472 
473  deallocate( iflgb )
474  deallocate( nch )
475  deallocate( wgtch )
476  deallocate( ngasabs )
477  deallocate( igasabs )
478 
479  deallocate( akd )
480  deallocate( skd )
481  deallocate( acfc_pow )
482 
483  deallocate( fitplk )
484  deallocate( fsol )
485  deallocate( sfc )
486  deallocate( rayleigh )
487 
488  deallocate( qmol )
489  deallocate( q )
490 
491  deallocate( hygro_flag )
492  deallocate( radmode )
493 
494  return
495  end subroutine atmos_phy_rd_mstrnx_finalize
496 
497  !-----------------------------------------------------------------------------
499  subroutine atmos_phy_rd_mstrnx_flux( &
500  KA, KS, KE, IA, IS, IE, JA, JS, JE, &
501  DENS, TEMP, PRES, QV, &
502  CZ, FZ, &
503  fact_ocean, &
504  fact_land, &
505  fact_urban, &
506  temp_sfc, albedo_sfc, &
507  solins, cosSZA, &
508  CLDFRAC, MP_Re, MP_Qe, &
509  AE_Re, AE_Qe, &
510  flux_rad, &
511  flux_rad_top, &
512  flux_rad_sfc_dn, &
513  dtau_s, dem_s )
514  !Jval )
515 #ifdef _OPENACC
516  use openacc
517 #endif
518  use scale_const, only: &
519  eps => const_eps, &
520  mdry => const_mdry, &
521  mvap => const_mvap, &
522  ppm => const_ppm
523  use scale_time, only: &
524  time_nowdate
525  use scale_atmos_grid_cartesc_real, only: &
526  atmos_grid_cartesc_real_basepoint_lat
527  use scale_atmos_hydrometeor, only: &
528  n_hyd, &
529  i_hc, &
530  i_hi, &
531  hyd_dens
532  use scale_atmos_aerosol, only: &
533  n_ae, &
534  ae_dens
535  use scale_atmos_phy_rd_profile, only: &
536  rd_profile_read => atmos_phy_rd_profile_read, &
537  rd_profile_use_climatology => atmos_phy_rd_profile_use_climatology
538  implicit none
539  integer, intent(in) :: ka, ks, ke
540  integer, intent(in) :: ia, is, ie
541  integer, intent(in) :: ja, js, je
542 
543  real(rp), intent(in) :: dens (ka,ia,ja)
544  real(rp), intent(in) :: temp (ka,ia,ja)
545  real(rp), intent(in) :: pres (ka,ia,ja)
546  real(rp), intent(in) :: qv (ka,ia,ja)
547  real(rp), intent(in) :: cz ( ka,ia,ja)
548  real(rp), intent(in) :: fz (0:ka,ia,ja)
549  real(rp), intent(in) :: fact_ocean (ia,ja)
550  real(rp), intent(in) :: fact_land (ia,ja)
551  real(rp), intent(in) :: fact_urban (ia,ja)
552  real(rp), intent(in) :: temp_sfc (ia,ja)
553  real(rp), intent(in) :: albedo_sfc (ia,ja,n_rad_dir,n_rad_rgn)
554  real(rp), intent(in) :: solins (ia,ja)
555  real(rp), intent(in) :: cossza (ia,ja)
556  real(rp), intent(in) :: cldfrac (ka,ia,ja)
557  real(rp), intent(in) :: mp_re (ka,ia,ja,n_hyd)
558  real(rp), intent(in) :: mp_qe (ka,ia,ja,n_hyd)
559  real(rp), intent(in) :: ae_re (ka,ia,ja,n_ae)
560  real(rp), intent(in) :: ae_qe (ka,ia,ja,n_ae)
561  real(rp), intent(out) :: flux_rad (ka,ia,ja,2,2,2)
562  real(rp), intent(out) :: flux_rad_top (ia,ja,2,2,2)
563  real(rp), intent(out) :: flux_rad_sfc_dn(ia,ja,n_rad_dir,n_rad_rgn)
564 
565  real(rp), intent(out), optional :: dtau_s(ka,ia,ja) ! 0.67 micron cloud optical depth
566  real(rp), intent(out), optional :: dem_s (ka,ia,ja) ! 10.5 micron cloud emissivity
567 ! real(RP), intent(out), optional :: Jval (KA,IA,JA,CH_QA_photo)
568 
569  integer :: tropopause(ia,ja)
570  real(rp) :: gamma
571 
572  real(rp), parameter :: min_cldfrac = 1.e-8_rp
573 
574  real(rp) :: rhodz_merge (rd_kmax,ia,ja)
575  real(rp) :: pres_merge (rd_kmax,ia,ja)
576  real(rp) :: temp_merge (rd_kmax,ia,ja)
577  real(rp) :: temph_merge (rd_kmax+1,ia,ja)
578 
579  real(rp) :: gas_merge (rd_kmax,ia,ja,mstrn_ngas)
580  real(rp) :: cfc_merge (rd_kmax,ia,ja,mstrn_ncfc)
581  real(rp) :: aerosol_conc_merge(rd_kmax,ia,ja,rd_naero )
582  real(rp) :: aerosol_radi_merge(rd_kmax,ia,ja,rd_naero )
583  real(rp) :: cldfrac_merge (rd_kmax,ia,ja)
584 
585  ! output
586  real(rp) :: flux_rad_merge(rd_kmax+1,ia,ja,2,2,mstrn_ncloud)
587  real(rp) :: taucld_067u (rd_kmax,ia,ja) ! 0.67 micron cloud optical depth
588  real(rp) :: emiscld_105u (rd_kmax,ia,ja) ! 10.5 micron cloud emissivity
589 
590  real(rp) :: zerosw
591 
592  integer :: ihydro, iaero
593  integer :: rd_k, k, i, j, v, ic
594  !---------------------------------------------------------------------------
595 
596  log_progress(*) 'atmosphere / physics / radiation / mstrnX'
597 
598  call prof_rapstart('RD_Profile', 3)
599 
600  !$acc data copyin(DENS,TEMP,PRES,QV,CZ,FZ,fact_ocean,fact_land,fact_urban,temp_sfc,albedo_sfc,solins,cosSZA,cldfrac,MP_Re,AE_Re,AE_Qe) &
601  !$acc copyout(flux_rad,flux_rad_top,flux_rad_sfc_dn) &
602  !$acc create(tropopause,rhodz_merge,pres_merge,temp_merge,temph_merge,gas_merge,cfc_merge,aerosol_conc_merge,aerosol_radi_merge,cldfrac_merge,flux_rad_merge,tauCLD_067u,emisCLD_105u)
603 
604  !$acc data copyout(dtau_s) if(acc_is_present(dtau_s))
605  !$acc data copyout(dem_s) if(acc_is_present(dem_s))
606 
607  if ( atmos_phy_rd_mstrn_only_tropocloud ) then
608  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
609  !$omp private(k,i,j, &
610  !$omp gamma) &
611  !$omp shared (tropopause,pres,temp,CZ, &
612  !$omp KS,KE,IS,IE,JS,JE)
613  !$acc kernels
614  do j = js, je
615  do i = is, ie
616  tropopause(i,j) = ke+1
617  do k = ke, ks, -1
618  if ( pres(k,i,j) >= 300.e+2_rp ) then
619  exit
620  elseif( pres(k,i,j) < 50.e+2_rp ) then
621  tropopause(i,j) = k
622  else
623  gamma = ( temp(k+1,i,j) - temp(k,i,j) ) / ( cz(k+1,i,j) - cz(k,i,j) )
624  if( gamma > 1.e-4_rp ) tropopause(i,j) = k
625  endif
626  enddo
627  enddo
628  enddo
629  !$acc end kernels
630  else
631  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
632  !$omp private(i,j) &
633  !$omp shared (tropopause, &
634  !$omp KE,IS,IE,JS,JE)
635 !OCL XFILL
636  !$acc kernels
637  do j = js, je
638  do i = is, ie
639  tropopause(i,j) = ke+1
640  end do
641  end do
642  !$acc end kernels
643  endif
644 
645  ! marge basic profile and value in model domain
646 
647  if ( rd_profile_use_climatology ) then
648  call rd_profile_read( rd_kmax, & ! [IN]
649  mstrn_ngas, & ! [IN]
650  mstrn_ncfc, & ! [IN]
651  rd_naero, & ! [IN]
652  atmos_grid_cartesc_real_basepoint_lat, & ! [IN]
653  time_nowdate(:), & ! [IN]
654  rd_zh(:), & ! [IN]
655  rd_z(:), & ! [IN]
656  rd_rhodz(:), & ! [OUT]
657  rd_pres(:), & ! [OUT]
658  rd_presh(:), & ! [OUT]
659  rd_temp(:), & ! [OUT]
660  rd_temph(:), & ! [OUT]
661  rd_gas(:,:), & ! [OUT]
662  rd_cfc(:,:), & ! [OUT]
663  rd_aerosol_conc(:,:), & ! [OUT]
664  rd_aerosol_radi(:,:), & ! [OUT]
665  rd_cldfrac(:) ) ! [OUT]
666  endif
667 
668 !OCL XFILL
669  !$omp parallel do default(none) &
670  !$omp shared(JS,JE,IS,IE,RD_KADD,temph_merge,RD_temph,KE,RD_KMAX,KS,temp,CZ,FZ) &
671  !$omp private(i,j,k,RD_k) OMP_SCHEDULE_ collapse(2)
672  !$acc kernels
673  do j = js, je
674  do i = is, ie
675  do rd_k = 1, rd_kadd
676  temph_merge(rd_k,i,j) = rd_temph(rd_k)
677  enddo
678 
679  temph_merge(rd_kadd+1,i,j) = temp(ke,i,j)
680  do rd_k = rd_kadd+2, rd_kmax
681  k = ks + rd_kmax - rd_k ! reverse axis
682 
683  temph_merge(rd_k,i,j) = ( temp(k ,i,j) * (cz(k+1,i,j)-fz(k,i,j)) / (cz(k+1,i,j)-cz(k,i,j)) &
684  + temp(k+1,i,j) * (fz(k ,i,j)-cz(k,i,j)) / (cz(k+1,i,j)-cz(k,i,j)) )
685  enddo
686  temph_merge(rd_kmax+1,i,j) = temp(ks,i,j)
687  enddo
688  enddo
689  !$acc end kernels
690 
691 !OCL XFILL
692  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
693  !$omp private(k,i,j,RD_k) &
694  !$omp shared (RD_temp,dens,FZ,pres,temp, &
695  !$omp rhodz_merge,RD_rhodz,pres_merge,RD_pres,temp_merge, &
696  !$omp KS,JS,JE,IS,IE,RD_KMAX,RD_KADD)
697  !$acc kernels
698  do j = js, je
699  do i = is, ie
700  do rd_k = 1, rd_kadd
701  rhodz_merge(rd_k,i,j) = rd_rhodz(rd_k)
702  pres_merge(rd_k,i,j) = rd_pres(rd_k)
703  temp_merge(rd_k,i,j) = rd_temp(rd_k)
704  enddo
705 
706  do rd_k = rd_kadd+1, rd_kmax
707  k = ks + rd_kmax - rd_k ! reverse axis
708 
709  rhodz_merge(rd_k,i,j) = dens(k,i,j) * ( fz(k,i,j)-fz(k-1,i,j) ) ! [kg/m2]
710  pres_merge(rd_k,i,j) = pres(k,i,j) * 1.e-2_rp ! [hPa]
711  temp_merge(rd_k,i,j) = temp(k,i,j)
712  enddo
713  enddo
714  enddo
715  !$acc end kernels
716 
717 !OCL XFILL
718  !$omp parallel default(none) &
719  !$omp private(v,i,j,RD_k) &
720  !$omp shared (gas_merge,RD_gas, &
721  !$omp IS,IE,JS,JE,RD_KMAX)
722  !$acc kernels
723  do v = 1, mstrn_ngas
724  !$omp do OMP_SCHEDULE_ collapse(2)
725  do j = js, je
726  do i = is, ie
727  do rd_k = 1, rd_kmax
728  gas_merge(rd_k,i,j,v) = rd_gas(rd_k,v)
729  enddo
730  enddo
731  enddo
732  !$omp end do nowait
733  enddo
734  !$acc end kernels
735  !$omp end parallel
736 
737  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
738  !$omp private(k,i,j,RD_k, &
739  !$omp zerosw ) &
740  !$omp shared (gas_merge,QV,EPS,Mvap,Mdry, &
741  !$omp KS,IS,IE,JS,JE,RD_KADD,RD_KMAX)
742  !$acc kernels
743  do j = js, je
744  do i = is, ie
745  do rd_k = rd_kadd+1, rd_kmax
746  k = ks + rd_kmax - rd_k ! reverse axis
747  zerosw = sign(0.5_rp, qv(k,i,j)-eps) + 0.5_rp
748  gas_merge(rd_k,i,j,1) = qv(k,i,j) / mvap * mdry / ppm * zerosw ! [PPM]
749  enddo
750  enddo
751  enddo
752  !$acc end kernels
753 
754 !OCL XFILL
755  !$omp parallel default(none) &
756  !$omp private(v,i,j,RD_k) &
757  !$omp shared (cfc_merge,RD_cfc, &
758  !$omp IS,IE,JS,JE,RD_KMAX)
759  !$acc kernels
760  do v = 1, mstrn_ncfc
761  !$omp do OMP_SCHEDULE_ collapse(2)
762  do j = js, je
763  do i = is, ie
764  do rd_k = 1, rd_kmax
765  cfc_merge(rd_k,i,j,v) = rd_cfc(rd_k,v)
766  enddo
767  enddo
768  enddo
769  !$omp end do nowait
770  enddo
771  !$acc end kernels
772  !$omp end parallel
773 
774  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
775  !$omp private(k,i,j,RD_k) &
776  !$omp shared (cldfrac_merge,RD_cldfrac,cldfrac, &
777  !$omp KS,IS,IE,JS,JE,RD_KMAX,RD_KADD)
778  !$acc kernels
779  do j = js, je
780  do i = is, ie
781  do rd_k = 1, rd_kadd
782  cldfrac_merge(rd_k,i,j) = rd_cldfrac(rd_k)
783  enddo
784 
785  do rd_k = rd_kadd+1, rd_kmax
786  k = ks + rd_kmax - rd_k ! reverse axis
787 
788  cldfrac_merge(rd_k,i,j) = min( max( cldfrac(k,i,j), 0.0_rp ), 1.0_rp )
789 ! cldfrac_merge(RD_k,i,j) = 0.5_RP + sign( 0.5_RP, cldfrac(k,i,j)-min_cldfrac )
790  enddo
791  enddo
792  enddo
793  !$acc end kernels
794 
795 !OCL XFILL
796  !$omp parallel default(none) &
797  !$omp private(v,i,j,RD_k) &
798  !$omp shared (aerosol_conc_merge,RD_aerosol_conc,aerosol_radi_merge,RD_aerosol_radi, &
799  !$omp IS,IE,JS,JE,RD_KADD)
800  !$acc kernels
801  do v = 1, rd_naero
802  !$omp do OMP_SCHEDULE_ collapse(2)
803  do j = js, je
804  do i = is, ie
805  do rd_k = 1, rd_kadd
806  aerosol_conc_merge(rd_k,i,j,v) = rd_aerosol_conc(rd_k,v)
807  aerosol_radi_merge(rd_k,i,j,v) = rd_aerosol_radi(rd_k,v)
808  enddo
809  enddo
810  enddo
811  !$omp end do nowait
812  enddo
813  !$acc end kernels
814  !$omp end parallel
815 
816  !$omp parallel default(none) &
817  !$omp private(ihydro,k,i,j,RD_k) &
818  !$omp shared (aerosol_conc_merge,tropopause,MP_Qe,HYD_DENS,RHO_std, &
819  !$omp ATMOS_PHY_RD_MSTRN_ONLY_QCI, &
820  !$omp KS,KE,IS,IE,JS,JE,RD_KMAX,RD_KADD)
821  do ihydro = 1, n_hyd
822  if ( atmos_phy_rd_mstrn_only_qci .and. &
823  ( ihydro /= i_hc .and. ihydro /= i_hi ) ) then
824 !OCL XFILL
825  !$omp do OMP_SCHEDULE_ collapse(2)
826  !$acc kernels
827  do j = js, je
828  do i = is, ie
829  do rd_k = rd_kadd+1, rd_kmax
830  aerosol_conc_merge(rd_k,i,j,ihydro) = 0.0_rp
831  end do
832  end do
833  end do
834  !$acc end kernels
835  !$omp end do nowait
836  else
837  !$omp do OMP_SCHEDULE_ collapse(2)
838  !$acc kernels
839  do j = js, je
840  do i = is, ie
841  do rd_k = rd_kadd+1, rd_kadd+1 + ke - tropopause(i,j)
842  aerosol_conc_merge(rd_k,i,j,ihydro) = 0.0_rp
843  end do
844  !$acc loop independent
845  do rd_k = rd_kadd+1 + ke - tropopause(i,j) + 1, rd_kmax
846  k = ks + rd_kmax - rd_k ! reverse axis
847  aerosol_conc_merge(rd_k,i,j,ihydro) = max( mp_qe(k,i,j,ihydro), 0.0_rp ) &
848  / hyd_dens(ihydro) * rho_std / ppm ! [PPM to standard air]
849  enddo
850  enddo
851  enddo
852  !$acc end kernels
853  !$omp end do nowait
854  end if
855  enddo
856  !$omp end parallel
857 
858  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(3) &
859  !$omp private(ihydro,k,i,j,RD_k) &
860  !$omp shared (aerosol_radi_merge,MP_Re, &
861  !$omp KS,IS,IE,JS,JE,RD_KMAX,RD_KADD)
862  !$acc kernels
863  do ihydro = 1, n_hyd
864  do j = js, je
865  do i = is, ie
866  do rd_k = rd_kadd+1, rd_kmax
867  k = ks + rd_kmax - rd_k ! reverse axis
868  aerosol_radi_merge(rd_k,i,j,ihydro) = mp_re(k,i,j,ihydro)
869  enddo
870  enddo
871  enddo
872  enddo
873  !$acc end kernels
874 
875  !$omp parallel default(none) &
876  !$omp private(iaero,k,i,j,RD_k) &
877  !$omp shared (aerosol_conc_merge,RD_aerosol_conc,aerosol_radi_merge,RD_aerosol_radi, &
878  !$omp ATMOS_PHY_RD_MSTRN_USE_AERO, &
879  !$omp AE_Qe,RHO_std,AE_Re, &
880  !$omp KS,IS,IE,JS,JE,RD_KMAX,RD_KADD)
881  do iaero = 1, n_ae
882 
883  if ( atmos_phy_rd_mstrn_use_aero ) then
884  !$omp do OMP_SCHEDULE_ collapse(2)
885  !$acc kernels
886  do j = js, je
887  do i = is, ie
888  do rd_k = rd_kadd+1, rd_kmax
889  k = ks + rd_kmax - rd_k ! reverse axis
890  aerosol_conc_merge(rd_k,i,j,n_hyd+iaero) = max( ae_qe(k,i,j,iaero), 0.0_rp ) &
891  / ae_dens(iaero) * rho_std / ppm ! [PPM to standard air]
892  aerosol_radi_merge(rd_k,i,j,n_hyd+iaero) = ae_re(k,i,j,iaero)
893  enddo
894  enddo
895  enddo
896  !$acc end kernels
897  !$omp end do nowait
898  else
899  !$omp do OMP_SCHEDULE_ collapse(2)
900  !$acc kernels
901  do j = js, je
902  do i = is, ie
903  do rd_k = rd_kadd+1, rd_kmax
904  aerosol_conc_merge(rd_k,i,j,n_hyd+iaero) = rd_aerosol_conc(rd_k,n_hyd+iaero)
905  aerosol_radi_merge(rd_k,i,j,n_hyd+iaero) = rd_aerosol_radi(rd_k,n_hyd+iaero)
906  enddo
907  enddo
908  enddo
909  !$acc end kernels
910  !$omp end do nowait
911  endif
912 
913  enddo
914  !$omp end parallel
915 
916  call prof_rapend ('RD_Profile', 3)
917  call prof_rapstart('RD_MSTRN_DTRN3', 3)
918 
919  ! calc radiative transfer
920 #ifndef _OPENACC
921  !$omp parallel do
922  do j = js, je
923  do i = is, ie
924 #endif
925  call rd_mstrn_dtrn3( rd_kmax, &
926 #ifdef _OPENACC
927  ia, is, ie, ja, js, je, &
928 #else
929  ia, i, i, ja, j, j, &
930  i, j, &
931 #endif
932  mstrn_ngas, mstrn_ncfc, rd_naero, & ! [IN]
933  rd_hydro_str, rd_hydro_end, rd_aero_str, rd_aero_end, & ! [IN]
934  solins(:,:), cossza(:,:), & ! [IN]
935  rhodz_merge(:,:,:), pres_merge(:,:,:), & ! [IN]
936  temp_merge(:,:,:), temph_merge(:,:,:), temp_sfc(:,:), & ! [IN]
937  gas_merge(:,:,:,:), cfc_merge(:,:,:,:), & ! [IN]
938  aerosol_conc_merge(:,:,:,:), aerosol_radi_merge(:,:,:,:), & ! [IN]
939  i_mpae2rd(:), & ! [IN]
940  cldfrac_merge(:,:,:), & ! [IN]
941  albedo_sfc(:,:,:,:), & ! [IN]
942  fact_ocean(:,:), fact_land(:,:), fact_urban(:,:), & ! [IN]
943  flux_rad_merge(:,:,:,:,:,:), flux_rad_sfc_dn(:,:,:,:), & ! [OUT]
944  taucld_067u(:,:,:), emiscld_105u(:,:,:) ) ! [OUT]
945 #ifndef _OPENACC
946  end do
947  end do
948 #endif
949 
950 
951  call prof_rapend ('RD_MSTRN_DTRN3', 3)
952 
953  ! return to grid coordinate of model domain
954  !$omp parallel default(none) &
955  !$omp private(ic,k,i,j,RD_k) &
956  !$omp shared (flux_rad,flux_rad_merge, &
957  !$omp KS,IS,IE,JS,JE,RD_KMAX,RD_KADD)
958  !$acc kernels
959  !$acc loop independent collapse(4)
960  do ic = 1, 2
961  !$omp do OMP_SCHEDULE_ collapse(2)
962  do j = js, je
963  do i = is, ie
964  do rd_k = rd_kadd+1, rd_kmax+1
965  k = ks + rd_kmax - rd_k ! reverse axis
966 
967  flux_rad(k,i,j,i_lw,i_up,ic) = flux_rad_merge(rd_k,i,j,i_lw,i_up,ic)
968  flux_rad(k,i,j,i_lw,i_dn,ic) = flux_rad_merge(rd_k,i,j,i_lw,i_dn,ic)
969  flux_rad(k,i,j,i_sw,i_up,ic) = flux_rad_merge(rd_k,i,j,i_sw,i_up,ic)
970  flux_rad(k,i,j,i_sw,i_dn,ic) = flux_rad_merge(rd_k,i,j,i_sw,i_dn,ic)
971  enddo
972  enddo
973  enddo
974  !$omp end do nowait
975  enddo
976  !$acc end kernels
977  !$omp end parallel
978 
979 !OCL XFILL
980  !$omp parallel default(none) &
981  !$omp private(ic,i,j) &
982  !$omp shared (flux_rad_top,flux_rad_merge, &
983  !$omp IS,IE,JS,JE)
984  !$acc kernels
985  do ic = 1, 2
986  !$omp do OMP_SCHEDULE_
987  do j = js, je
988  do i = is, ie
989  flux_rad_top(i,j,i_lw,i_up,ic) = flux_rad_merge(1,i,j,i_lw,i_up,ic)
990  flux_rad_top(i,j,i_lw,i_dn,ic) = flux_rad_merge(1,i,j,i_lw,i_dn,ic)
991  flux_rad_top(i,j,i_sw,i_up,ic) = flux_rad_merge(1,i,j,i_sw,i_up,ic)
992  flux_rad_top(i,j,i_sw,i_dn,ic) = flux_rad_merge(1,i,j,i_sw,i_dn,ic)
993  enddo
994  enddo
995  !$omp end do nowait
996  enddo
997  !$acc end kernels
998  !$omp end parallel
999 
1000  if ( present( dtau_s ) ) then
1001  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
1002  !$omp private(k,i,j,RD_k) &
1003  !$omp shared (dtau_s,tauCLD_067u, &
1004  !$omp KS,IS,IE,JS,JE,RD_KMAX,RD_KADD)
1005  !$acc kernels
1006  !$acc loop independent collapse(3)
1007  do j = js, je
1008  do i = is, ie
1009  do rd_k = rd_kadd+1, rd_kmax
1010  k = ks + rd_kmax - rd_k ! reverse axis
1011  dtau_s(k,i,j) = taucld_067u(rd_k,i,j)
1012  enddo
1013  enddo
1014  enddo
1015  !$acc end kernels
1016  end if
1017 
1018  if ( present( dem_s ) ) then
1019  !$omp parallel do default(none) OMP_SCHEDULE_ collapse(2) &
1020  !$omp private(k,i,j,RD_k) &
1021  !$omp shared (dem_s,emisCLD_105u, &
1022  !$omp KS,IS,IE,JS,JE,RD_KMAX,RD_KADD)
1023  !$acc kernels
1024  !$acc loop independent collapse(3)
1025  do j = js, je
1026  do i = is, ie
1027  do rd_k = rd_kadd+1, rd_kmax
1028  k = ks + rd_kmax - rd_k ! reverse axis
1029  dem_s(k,i,j) = emiscld_105u(rd_k,i,j)
1030  enddo
1031  enddo
1032  enddo
1033  !$acc end kernels
1034  end if
1035 
1036  !$acc end data
1037  !$acc end data
1038  !$acc end data
1039 
1040  return
1041  end subroutine atmos_phy_rd_mstrnx_flux
1042 
1043  !-----------------------------------------------------------------------------
1045  subroutine rd_mstrn_setup( &
1046  ngas, &
1047  ncfc )
1048  use scale_prc, only: &
1049  prc_abort
1050  use scale_const, only: &
1051  grav => const_grav, &
1052  rdry => const_rdry, &
1053  pstd => const_pstd, &
1054  tem00 => const_tem00
1055  implicit none
1056 
1057  integer, intent(out) :: ngas
1058  integer, intent(out) :: ncfc
1059 
1060  real(rp) :: acfc(mstrn_ncfc)
1061 
1062  integer :: nband, nstream, nfitp, nfitt, nflag
1063  integer :: nsfc, nptype, nplkord, nfitplk
1064  integer :: nradius
1065 
1066  character(len=H_LONG) :: dummy, fname
1067 
1068  integer :: fid, ierr
1069  integer :: iw, ich, ip, it, igas, icfc, iptype, im
1070  !---------------------------------------------------------------------------
1071 
1072  !---< gas absorption parameter input >---
1073  ngas = mstrn_ngas
1074 
1075  ! allocate arrays
1076  allocate( waveh(mstrn_nband+1) )
1077  allocate( logfitp(mstrn_nfitp) )
1078  allocate( fitt(mstrn_nfitt) )
1079  allocate( logfitt(mstrn_nfitt) )
1080 
1081  allocate( iflgb(mstrn_nflag, mstrn_nband) )
1082  allocate( nch( mstrn_nband) )
1083  allocate( wgtch(mstrn_ch_limit,mstrn_nband) )
1084  allocate( ngasabs( mstrn_nband) )
1085  allocate( igasabs(mstrn_ngas, mstrn_nband) )
1086 
1087  allocate( akd(mstrn_nfitp,mstrn_ch_limit,mstrn_nfitt,mstrn_ngas,mstrn_nband) )
1088  allocate( skd(mstrn_nfitp,mstrn_ch_limit,mstrn_nfitt, mstrn_nband) )
1089  allocate( acfc_pow(mstrn_ncfc,mstrn_nband) )
1090 
1091  fid = io_get_available_fid()
1092  call io_get_fname(fname, mstrn_gaspara_inputfile)
1093  open( fid, &
1094  file = fname, &
1095  form = 'formatted', &
1096  status = 'old', &
1097  iostat = ierr )
1098 
1099  if ( ierr /= 0 ) then
1100  log_error("RD_MSTRN_setup",*) 'Input data file does not found! ', trim(fname)
1101  call prc_abort
1102  endif
1103 
1104  ! read gas parameters for check
1105  read(fid,*) nband, nstream, nfitp, nfitt, nflag, ncfc
1106 
1107  if ( nband /= mstrn_nband ) then
1108  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nband(given,file)=', mstrn_nband, nband
1109  call prc_abort
1110  endif
1111  if ( nstream /= mstrn_nstream ) then
1112  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nstream(given,file)=', mstrn_nstream, nstream
1113  call prc_abort
1114  endif
1115  if ( nfitp /= mstrn_nfitp ) then
1116  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nfitP(given,file)=', mstrn_nfitp, nfitp
1117  call prc_abort
1118  endif
1119  if ( nfitt /= mstrn_nfitt ) then
1120  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nfitT(given,file)=', mstrn_nfitt, nfitt
1121  call prc_abort
1122  endif
1123  if ( nflag /= mstrn_nflag ) then
1124  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nflag(given,file)=', mstrn_nflag, nflag
1125  call prc_abort
1126  endif
1127  if ( ncfc /= mstrn_ncfc ) then
1128  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! ncfc(given,file)=', mstrn_ncfc, ncfc
1129  call prc_abort
1130  endif
1131 
1132  ! wave band boundaries
1133  read(fid,*) dummy
1134  read(fid,*) waveh(:)
1135  ! fitting point for log(pressure)
1136  read(fid,*) dummy
1137  read(fid,*) logfitp(:)
1138  ! fitting point for temperature
1139  read(fid,*) dummy
1140  read(fid,*) fitt(:)
1141 
1142  logfitt(:) = log10( fitt(:) )
1143 
1144  ! for each band
1145  do iw = 1, mstrn_nband
1146 
1147  ! optical properties flag
1148  read(fid,*) dummy
1149  read(fid,*) iflgb(:,iw)
1150  ! number of subintervals
1151  read(fid,*) dummy
1152  read(fid,*) nch(iw)
1153  ! weights for subintervals
1154  read(fid,*) dummy
1155  read(fid,*) (wgtch(ich,iw),ich=1,nch(iw))
1156  ! number of considering gases
1157  read(fid,*) dummy
1158  read(fid,*) ngasabs(iw)
1159 
1160  ! major gas absorption
1161  if ( ngasabs(iw) > 0 ) then
1162  do igas = 1, ngasabs(iw)
1163  read(fid,*) igasabs(igas,iw)
1164  do it = 1, mstrn_nfitt
1165  do ip = 1, mstrn_nfitp
1166  read(fid,*) (akd(ip,ich,it,igasabs(igas,iw),iw),ich=1,nch(iw))
1167  enddo
1168  enddo
1169  enddo
1170  endif
1171 
1172  ! H2O continuum
1173  if ( iflgb(i_h2o_continuum,iw) > 0 ) then
1174  read(fid,*) dummy
1175  do it = 1, mstrn_nfitt
1176  do ip = 1, mstrn_nfitp
1177  read(fid,*) (skd(ip,ich,it,iw),ich=1,nch(iw))
1178  enddo
1179  enddo
1180  endif
1181 
1182  ! CFC absorption
1183  if ( iflgb(i_cfc_continuum,iw) > 0 ) then
1184  read(fid,*) dummy
1185  read(fid,*) (acfc(icfc),icfc=1,mstrn_ncfc)
1186 
1187  do icfc = 1, mstrn_ncfc
1188  acfc_pow(icfc,iw) = 10.0_rp**acfc(icfc)
1189  end do
1190  endif
1191 
1192  enddo ! band loop
1193 
1194  close(fid)
1195 
1196  !---< aerosol(particle) parameter input >---
1197 
1198  ! allocate arrays
1199  allocate( fitplk(mstrn_nfitplk,mstrn_nband) )
1200  allocate( fsol( mstrn_nband) )
1201  allocate( sfc(mstrn_nsfc, mstrn_nband) )
1202  allocate( rayleigh( mstrn_nband) )
1203 
1204  allocate( qmol( mstrn_nmoment,mstrn_nband) )
1205  allocate( q(-1:mstrn_nradius+1,mstrn_nptype,mstrn_nmoment,mstrn_nband) )
1206  do iptype = 1, mstrn_nptype
1207  nradius = ptype_nradius(iptype)
1208 
1209  q( -1:0 ,iptype,:,:) = 0.d0 ! dummy for NaN
1210  q(nradius+1:mstrn_nradius+1,iptype,:,:) = 0.d0 ! dummy for extrapolation
1211  enddo
1212 
1213 
1214  call io_get_fname(fname, mstrn_aeropara_inputfile)
1215  open( fid, &
1216  file = fname, &
1217  form = 'formatted', &
1218  status = 'old', &
1219  iostat = ierr )
1220 
1221  if ( ierr /= 0 ) then
1222  log_error("RD_MSTRN_setup",*) 'Input data file does not found! ', trim(fname)
1223  call prc_abort
1224  endif
1225 
1226  ! read aerosol/surface parameters for check
1227  read(fid,*) nband, nsfc, nptype, nstream, nplkord, nfitplk
1228 
1229  if ( nband /= mstrn_nband ) then
1230  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nband(given,file)=', mstrn_nband, nband
1231  call prc_abort
1232  endif
1233  if ( nsfc /= mstrn_nsfc ) then
1234  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nsfc(given,file)=', mstrn_nsfc, nsfc
1235  call prc_abort
1236  endif
1237  if ( nptype /= mstrn_nptype ) then
1238  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nptype(given,file)=', mstrn_nptype, nptype
1239  call prc_abort
1240  endif
1241  if ( nstream /= mstrn_nstream ) then
1242  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nstream(given,file)=', mstrn_nstream, nstream
1243  call prc_abort
1244  endif
1245  if ( nplkord /= mstrn_nplkord ) then
1246  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nplkord(given,file)=', mstrn_nplkord, nplkord
1247  call prc_abort
1248  endif
1249  if ( nfitplk /= mstrn_nfitplk ) then
1250  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nfitPLK(given,file)=', mstrn_nfitplk, nfitplk
1251  call prc_abort
1252  endif
1253 
1254  ! wave band boundaries
1255  read(fid,*) dummy
1256  read(fid,*) waveh(:)
1257 
1258  ! for each band
1259  do iw = 1, mstrn_nband
1260 
1261  ! fitting point for planck functions
1262  read(fid,*) dummy
1263  read(fid,*) fitplk(:,iw)
1264  ! solar insolation
1265  read(fid,*) dummy
1266  read(fid,*) fsol(iw)
1267  ! surface properties
1268  read(fid,*) dummy
1269  read(fid,*) sfc(:,iw)
1270  ! rayleigh scattering
1271  read(fid,*) dummy
1272  read(fid,*) rayleigh(iw)
1273 
1274  ! moments
1275  read(fid,*) dummy
1276  do im = 1, mstrn_nmoment
1277  ! for rayleigh scattering phase function
1278  read(fid,*) qmol(im,iw)
1279  ! for aerosol scattering phase function
1280  do iptype = 1, mstrn_nptype
1281  nradius = ptype_nradius(iptype)
1282  read(fid,*) q(1:nradius,iptype,im,iw)
1283  enddo
1284  enddo
1285 
1286  enddo
1287 
1288  close(fid)
1289 
1290  fsol_tot = 0.0_rp
1291  do iw = 1, mstrn_nband
1292  fsol_tot = fsol_tot + fsol(iw)
1293  enddo
1294 
1295  !---< radius mode & hygroscopic parameter input >---
1296 
1297  ! allocate arrays
1298  allocate( hygro_flag(mstrn_nptype) )
1299  allocate( radmode(mstrn_nptype,mstrn_nradius) )
1300 
1301  call io_get_fname(fname, mstrn_hygropara_inputfile)
1302  open( fid, &
1303  file = fname, &
1304  form = 'formatted', &
1305  status = 'old', &
1306  iostat = ierr )
1307 
1308  if ( ierr /= 0 ) then
1309  log_error("RD_MSTRN_setup",*) 'Input data file does not found! ', trim(fname)
1310  call prc_abort
1311  endif
1312 
1313  read(fid,*) nptype
1314 
1315  if ( nptype /= mstrn_nptype ) then
1316  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nptype(given,file)=', mstrn_nptype, nptype
1317  call prc_abort
1318  endif
1319 
1320  do iptype = 1, mstrn_nptype
1321  read(fid,*) dummy
1322  read(fid,*) hygro_flag(iptype), nradius
1323 
1324  if ( nradius /= ptype_nradius(iptype) ) then
1325  log_error("RD_MSTRN_setup",*) 'Inconsistent parameter value! nradius(given,file)=', ptype_nradius(iptype), nradius
1326  call prc_abort
1327  endif
1328 
1329  read(fid,*) radmode(iptype,1:nradius)
1330  enddo
1331 
1332  close(fid)
1333 
1334  !----- report data -----
1335  log_newline
1336  log_info("RD_MSTRN_setup",'(1x,A,F12.7)') 'Baseline of total solar insolation : ', fsol_tot
1337 
1338  !---< constant parameter for main scheme >---
1339  rho_std = pstd / ( rdry * tem00 ) ! [kg/m3]
1340 
1341  m(i_sw) = 1.0_rp / sqrt(3.0_rp)
1342  w(i_sw) = 1.0_rp
1343  wbar(i_sw) = 1.0_rp
1344 
1345  m(i_lw) = 1.0_rp / 1.66_rp
1346  w(i_lw) = 1.0_rp
1347  wbar(i_lw) = 2.0_rp * m(i_lw)
1348 
1349  do im = 1, 2
1350  wmns(im) = sqrt( w(im) / m(im) )
1351  wpls(im) = sqrt( w(im) * m(im) )
1352  wscale(im) = wpls(im) / wbar(im)
1353  end do
1354  !$acc update device(M, W, Wmns, Wpls, Wscale)
1355 
1356  !$acc enter data &
1357  !$acc copyin(waveh, &
1358  !$acc wgtch, fitPLK, logfitP, logfitT, fitT, &
1359  !$acc radmode, ngasabs, igasabs, ptype_nradius, &
1360  !$acc fsol, q, qmol, rayleigh, acfc_pow, nch, AKD, SKD )
1361 
1362  return
1363  end subroutine rd_mstrn_setup
1364 
1365  !-----------------------------------------------------------------------------
1367 !OCL SERIAL
1368  subroutine rd_mstrn_dtrn3( &
1369  rd_kmax, IA, IS, IE, JA, JS, JE, &
1370 #ifndef _OPENACC
1371  i, j, &
1372 #endif
1373  ngas, &
1374  ncfc, &
1375  naero, &
1376  hydro_str, &
1377  hydro_end, &
1378  aero_str, &
1379  aero_end, &
1380  solins, &
1381  cosSZA0, &
1382  rhodz, &
1383  pres, &
1384  temp, &
1385  temph, &
1386  temp_sfc, &
1387  gas, &
1388  cfc, &
1389  aerosol_conc, &
1390  aerosol_radi, &
1391  aero2ptype, &
1392  cldfrac, &
1393  albedo_sfc, &
1394  fact_ocean, &
1395  fact_land, &
1396  fact_urban, &
1397  rflux, &
1398  rflux_sfc_dn, &
1399  tauCLD_067u, &
1400  emisCLD_105u )
1401  use scale_prc, only: &
1402  prc_abort
1403  use scale_const, only: &
1404  grav => const_grav, &
1405  pstd => const_pstd, &
1406  ppm => const_ppm
1407  implicit none
1408  integer, intent(in) :: rd_kmax
1409  integer, intent(in) :: ia, is, ie
1410  integer, intent(in) :: ja, js, je
1411 #ifndef _OPENACC
1412  integer, intent(in) :: i, j
1413 #endif
1414 
1415  integer, intent(in) :: ngas
1416  integer, intent(in) :: ncfc
1417  integer, intent(in) :: naero
1418  integer, intent(in) :: hydro_str
1419  integer, intent(in) :: hydro_end
1420  integer, intent(in) :: aero_str
1421  integer, intent(in) :: aero_end
1422  real(rp), intent(in) :: solins (ia,ja)
1423  real(rp), intent(in) :: cossza0 (ia,ja)
1424  real(rp), intent(in) :: rhodz (rd_kmax ,ia,ja)
1425  real(rp), intent(in) :: pres (rd_kmax ,ia,ja)
1426  real(rp), intent(in) :: temp (rd_kmax ,ia,ja)
1427  real(rp), intent(in) :: temph (rd_kmax+1,ia,ja)
1428  real(rp), intent(in) :: temp_sfc (ia,ja)
1429  real(rp), intent(in) :: gas (rd_kmax,ia,ja,ngas )
1430  real(rp), intent(in) :: cfc (rd_kmax,ia,ja,ncfc )
1431  real(rp), intent(in) :: aerosol_conc(rd_kmax,ia,ja,naero)
1432  real(rp), intent(in) :: aerosol_radi(rd_kmax,ia,ja,naero)
1433  integer, intent(in) :: aero2ptype (naero)
1434  real(rp), intent(in) :: cldfrac (rd_kmax,ia,ja)
1435  real(rp), intent(in) :: albedo_sfc (ia,ja,n_rad_dir,n_rad_rgn)
1436  real(rp), intent(in) :: fact_ocean (ia,ja)
1437  real(rp), intent(in) :: fact_land (ia,ja)
1438  real(rp), intent(in) :: fact_urban (ia,ja)
1439  real(rp), intent(inout) :: rflux (rd_kmax+1,ia,ja,2,2,mstrn_ncloud)
1440  real(rp), intent(inout) :: rflux_sfc_dn(ia,ja,n_rad_dir,n_rad_rgn) ! surface downward radiation flux (direct/diffuse,IR/NIR/VIS)
1441  real(rp), intent(inout) :: taucld_067u (rd_kmax,ia,ja) ! 0.67 micron cloud optical depth
1442  real(rp), intent(inout) :: emiscld_105u(rd_kmax,ia,ja) ! 10.5 micron cloud emissivity
1443 
1444 #ifdef _OPENACC
1445 #define LOOP_INNER do j=JS,JE; do i=IS,IE
1446 #define LOOP_END_INNER end do; end do
1447 #else
1448 #define LOOP_INNER
1449 #define LOOP_END_INNER
1450 #endif
1451 
1452  ! for P-T fitting
1453  real(rp) :: dz_std (rd_kmax,is:ie,js:je) ! layer thickness at 0C, 1atm [cm]
1454  real(rp) :: logp ! log10(pres)
1455  real(rp) :: logt ! log10(temp)
1456  integer :: indexp (rd_kmax,is:ie,js:je) ! index for interpolation in P-fitting
1457  real(rp) :: factp (rd_kmax,is:ie,js:je) ! interpolation factor in P-fitting
1458  real(rp) :: factt32(rd_kmax,is:ie,js:je) ! interpolation factor in T-fitting
1459  real(rp) :: factt21(rd_kmax,is:ie,js:je) ! interpolation factor in T-fitting
1460  integer :: indexr (rd_kmax,naero,is:ie,js:je) ! index for interpolation in R-fitting
1461  real(rp) :: factr (rd_kmax,naero,is:ie,js:je) ! interpolation factor in R-fitting
1462 
1463  ! for optical thickness by gas
1464  real(rp) :: taugas(rd_kmax,mstrn_ch_limit,is:ie,js:je) ! optical thickness by gas absorption (total)
1465  real(rp) :: a1, a2, a3, factpt
1466  real(rp) :: qv, length
1467  integer :: gasno
1468 
1469  ! for optical thickness by particles
1470  real(rp) :: taupr (rd_kmax,mstrn_ncloud) ! optical thickness by Rayleigh/cloud/aerosol
1471  real(rp) :: omgpr (rd_kmax,mstrn_ncloud) ! single scattering albedo by Rayleigh/cloud/aerosol
1472  real(dp) :: optparam(rd_kmax,mstrn_nmoment,mstrn_ncloud,is:ie,js:je) ! optical parameters
1473  real(rp) :: q_fit, dp_p
1474 
1475  ! for planck functions
1476  real(rp) :: bbar (rd_kmax ) ! planck functions for thermal source at the interface
1477  real(rp) :: bbarh(rd_kmax+1) ! planck functions for thermal source at the center
1478  real(rp) :: b_sfc ! planck functions for thermal source at the surface
1479  real(rp) :: wl, beta
1480 
1481  ! for two-stream
1482  real(rp) :: tau(rd_kmax, mstrn_ncloud,mstrn_ch_limit) ! total optical thickness
1483  real(rp) :: omg(rd_kmax, mstrn_ncloud,mstrn_ch_limit) ! single scattering albedo
1484  real(rp) :: g (rd_kmax,0:2,mstrn_ncloud ) ! two-stream approximation factors
1485  ! 0: always 1
1486  ! 1: asymmetry factor
1487  ! 2: truncation factor
1488  real(rp) :: b (rd_kmax,0:2,mstrn_ncloud,mstrn_ch_limit) ! planck expansion coefficients (zero if SW)
1489  real(rp) :: fsol_rgn ! solar insolation (zero if LW)
1490 
1491  real(rp) :: flux (rd_kmax+1,2,mstrn_ncloud,mstrn_ch_limit) ! upward/downward flux
1492  real(rp) :: flux_direct(rd_kmax+1, mstrn_ncloud,mstrn_ch_limit) ! downward flux (direct solar)
1493 
1494  ! for satellite simulator
1495  real(rp) :: emiscld(rd_kmax,mstrn_ch_limit) ! cloud emissivity
1496 
1497  ! work
1498  real(rp) :: cossza(ia,ja)
1499 
1500  ! works for two_stream
1501  ! main factors
1502  real(rp) :: tdir0(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! transmission factor for solar direct (clear-sky/cloud)
1503  real(rp) :: r0 (rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! reflection factor (clear-sky/cloud)
1504  real(rp) :: t0 (rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! transmission factor (clear-sky/cloud)
1505  real(rp) :: em_lw(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! thermal source (sfc->TOA) (clear-sky/cloud)
1506  real(rp) :: em_sw(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! solar source (sfc->TOA) (clear-sky/cloud)
1507  real(rp) :: ep_lw(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! thermal source (TOA->sfc) (clear-sky/cloud)
1508  real(rp) :: ep_sw(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! solar source (TOA->sfc) (clear-sky/cloud)
1509  ! Averaged factors, considering cloud overwrap
1510  real(rp) :: cf (rd_kmax ,mstrn_ncloud*mstrn_ch_limit) ! cloud fraction
1511  real(rp) :: tau_bar_sol(rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! solar insolation through accumulated optical thickness at each layer
1512  real(rp) :: r (rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! reflection factor
1513  real(rp) :: t (rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! transmission factor
1514  real(rp) :: em (rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! source (sfc->TOA)
1515  real(rp) :: ep (rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! source (TOA->sfc)
1516  real(rp) :: em0 (mstrn_ncloud)
1517  ! Doubling-Adding
1518  real(rp) :: r12mns(mstrn_ncloud*mstrn_ch_limit,rd_kmax+1) ! reflection factor in doubling method
1519  real(rp) :: r12pls(mstrn_ncloud*mstrn_ch_limit,rd_kmax+1) ! reflection factor in doubling method
1520  real(rp) :: e12mns(mstrn_ncloud*mstrn_ch_limit,rd_kmax+1) ! source function in doubling method
1521  real(rp) :: e12pls(mstrn_ncloud*mstrn_ch_limit,rd_kmax+1) ! source function in doubling method
1522 
1523 
1524  real(rp) :: zerosw
1525  real(rp) :: valsum(rd_kmax)
1526  integer :: chmax
1527  integer :: ip, ir, irgn, irgn_alb
1528  integer :: igas, icfc, iaero, iptype, nradius
1529  integer :: iw, ich, iplk, icloud, im, idir
1530  integer :: k
1531 #ifdef _OPENACC
1532  integer :: i, j
1533 #endif
1534  !---------------------------------------------------------------------------
1535 
1536  !$acc data &
1537  !$acc copyin(solins, cosSZA0, rhodz, pres, temp, temph, temp_sfc, gas, cfc, &
1538  !$acc aerosol_conc, aerosol_radi, aero2ptype, cldfrac, albedo_sfc, &
1539  !$acc fact_ocean, fact_land, fact_urban) &
1540  !$acc copyout(rflux, rflux_sfc_dn, tauCLD_067u, emisCLD_105u) &
1541  !$acc create(dz_std, indexP, factP, factT32, factT21, indexR, factR, cosSZA, optparam, tauGAS)
1542 
1543  !$acc kernels
1544  loop_inner
1545  cossza(i,j) = max( cossza0(i,j), rd_cossza_min )
1546  loop_end_inner
1547  !$acc end kernels
1548 
1549  !$acc parallel loop collapse(2)
1550  loop_inner
1551  do k = 1, rd_kmax
1552  dz_std(k,i,j) = rhodz(k,i,j) / rho_std * 100.0_rp ! [cm]
1553 
1554  logp = log10( pres(k,i,j) )
1555  logt = log10( temp(k,i,j) )
1556 
1557  indexp(k,i,j) = mstrn_nfitp
1558  !$acc loop seq
1559  do ip = mstrn_nfitp, 2, -1
1560  if( logp >= logfitp(ip) ) indexp(k,i,j) = ip
1561  enddo
1562 
1563  ip = indexp(k,i,j)
1564 
1565  factp(k,i,j) = ( logp - logfitp(ip-1) ) / ( logfitp(ip) - logfitp(ip-1) )
1566 
1567  factt32(k,i,j) = ( logt - logfitt(2) ) / ( logfitt(3) - logfitt(2) ) &
1568  * ( temp(k,i,j) - fitt(1) ) / ( fitt(3) - fitt(1) )
1569  factt21(k,i,j) = ( logt - logfitt(2) ) / ( logfitt(2) - logfitt(1) ) &
1570  * ( fitt(3) - temp(k,i,j) ) / ( fitt(3) - fitt(1) )
1571  enddo
1572  loop_end_inner
1573  !$acc end parallel
1574 
1575  !---< interpolation of mode radius & hygroscopic parameter (R-fitting) >---
1576  !$acc parallel loop collapse(3)
1577  loop_inner
1578  do iaero = 1, naero
1579  do k = 1, rd_kmax
1580  iptype = aero2ptype(iaero)
1581  nradius = ptype_nradius(iptype)
1582 
1583  ! [Note] Extrapolation sometimes makes unexpected value
1584  if ( aerosol_radi(k,i,j,iaero) <= radmode(iptype,1) ) then
1585  indexr(k,iaero,i,j) = 1
1586  factr(k,iaero,i,j) = 0.0_rp
1587  elseif( aerosol_radi(k,i,j,iaero) > radmode(iptype,nradius) ) then
1588  indexr(k,iaero,i,j) = nradius
1589  factr(k,iaero,i,j) = 1.0_rp
1590  else
1591  indexr(k,iaero,i,j) = -1
1592  !$acc loop seq
1593  do ir = 1, nradius-1
1594  if ( aerosol_radi(k,i,j,iaero) <= radmode(iptype,ir+1) &
1595  .AND. aerosol_radi(k,i,j,iaero) > radmode(iptype,ir ) ) then ! interpolation
1596  indexr(k,iaero,i,j) = ir
1597  factr(k,iaero,i,j) = ( aerosol_radi(k,i,j,iaero) - radmode(iptype,ir) ) &
1598  / ( radmode(iptype,ir+1) - radmode(iptype,ir) )
1599  exit
1600  endif
1601  enddo
1602  endif
1603  enddo
1604  enddo
1605  loop_end_inner
1606  !$acc end parallel
1607 
1608  ! initialize
1609  !$acc parallel loop collapse(5)
1610  do icloud = 1, mstrn_ncloud
1611  do idir = 1, 2
1612  do irgn = 1, 2
1613  loop_inner
1614  do k = 1, rd_kmax+1
1615  rflux(k,i,j,irgn,idir,icloud) = 0.0_rp
1616  end do
1617  loop_end_inner
1618  end do
1619  end do
1620  end do
1621  !$acc end parallel
1622  !$acc kernels
1623  do irgn_alb = 1, n_rad_rgn
1624  do idir = 1, n_rad_dir
1625  loop_inner
1626  rflux_sfc_dn(i,j,idir,irgn_alb) = 0.0_rp
1627  loop_end_inner
1628  end do
1629  end do
1630  !$acc end kernels
1631 
1632  do iw = 1, mstrn_nband
1633  if ( iflgb(i_swlw,iw) == 0 ) then ! Long wave region
1634  irgn = i_lw
1635  irgn_alb = i_r_ir ! IR
1636  else ! Short wave region
1637  irgn = i_sw
1638  if ( waveh(iw) >= 10000.0_rp / 0.7_rp ) then
1639  irgn_alb = i_r_vis ! Visible
1640  else
1641  irgn_alb = i_r_nir ! Near-IR
1642  endif
1643  endif
1644 
1645  chmax = nch(iw)
1646 
1647  !---< interpolation of gas parameters (P-T fitting) >---
1648  !$acc parallel loop collapse(2)
1649  loop_inner
1650  !$acc loop collapse(2)
1651  do ich = 1, chmax
1652  do k = 1, rd_kmax
1653  taugas(k,ich,i,j) = 0.0_rp
1654  enddo
1655  enddo
1656  loop_end_inner
1657  !$acc end parallel
1658 
1659  !--- Gas line absorption
1660  !$acc parallel loop collapse(2)
1661  loop_inner
1662  !$acc loop seq
1663  do igas = 1, ngasabs(iw)
1664  gasno = igasabs(igas,iw)
1665  !$acc loop collapse(2)
1666  do ich = 1, chmax
1667  do k = 1, rd_kmax
1668  ip = indexp(k,i,j)
1669  a1 = akd(ip-1,ich,1,gasno,iw) * ( 1.0_rp - factp(k,i,j) )&
1670  + akd(ip ,ich,1,gasno,iw) * ( factp(k,i,j) )
1671  a2 = akd(ip-1,ich,2,gasno,iw) * ( 1.0_rp - factp(k,i,j) )&
1672  + akd(ip ,ich,2,gasno,iw) * ( factp(k,i,j) )
1673  a3 = akd(ip-1,ich,3,gasno,iw) * ( 1.0_rp - factp(k,i,j) )&
1674  + akd(ip ,ich,3,gasno,iw) * ( factp(k,i,j) )
1675  factpt = factt32(k,i,j)*(a3-a2) + a2 + factt21(k,i,j)*(a2-a1)
1676 
1677  length = gas(k,i,j,igasabs(igas,iw)) * ppm * dz_std(k,i,j)
1678  taugas(k,ich,i,j) = taugas(k,ich,i,j) + 10.0_rp**factpt * length
1679  end do
1680  end do
1681  end do
1682  loop_end_inner
1683  !$acc end parallel
1684 
1685  !--- Gas broad absorption
1686  if ( iflgb(i_h2o_continuum,iw) == 1 ) then
1687  !$acc parallel loop collapse(2)
1688  loop_inner
1689  !$acc loop collapse(2)
1690  do ich = 1, chmax
1691  do k = 1, rd_kmax
1692  ip = indexp(k,i,j)
1693  a1 = skd(ip-1,ich,1,iw) * ( 1.0_rp-factp(k,i,j) )&
1694  + skd(ip ,ich,1,iw) * ( factp(k,i,j) )
1695  a2 = skd(ip-1,ich,2,iw) * ( 1.0_rp-factp(k,i,j) )&
1696  + skd(ip ,ich,2,iw) * ( factp(k,i,j) )
1697  a3 = skd(ip-1,ich,3,iw) * ( 1.0_rp-factp(k,i,j) )&
1698  + skd(ip ,ich,3,iw) * ( factp(k,i,j) )
1699  factpt = factt32(k,i,j)*(a3-a2) + a2 + factt21(k,i,j)*(a2-a1)
1700 
1701  qv = gas(k,i,j,1) * ppm * dz_std(k,i,j)
1702  length = qv*qv / ( qv + dz_std(k,i,j) )
1703  taugas(k,ich,i,j) = taugas(k,ich,i,j) + 10.0_rp**factpt * length
1704  enddo
1705  enddo
1706  loop_end_inner
1707  !$acc end parallel
1708  endif
1709 
1710  if ( iflgb(i_cfc_continuum,iw) == 1 ) then
1711  !$acc parallel loop collapse(2) private(valsum)
1712  loop_inner
1713  do k = 1, rd_kmax
1714  valsum(k) = 0.0_rp
1715  end do
1716  !$acc loop seq
1717  do icfc = 1, ncfc
1718  do k = 1, rd_kmax
1719  valsum(k) = valsum(k) + acfc_pow(icfc,iw) * cfc(k,i,j,icfc)
1720  end do
1721  end do
1722  do k = 1, rd_kmax
1723  valsum(k) = valsum(k) * ppm * dz_std(k,i,j)
1724  end do
1725 
1726  do ich = 1, chmax
1727  do k = 1, rd_kmax
1728  taugas(k,ich,i,j) = taugas(k,ich,i,j) + valsum(k)
1729  enddo
1730  enddo
1731  loop_end_inner
1732  !$acc end parallel
1733  endif
1734 
1735  !---< particle (Rayleigh/Cloud/Aerosol) >---
1736 
1737  ! optical thickness, phase function
1738  ! im=1: extinction coefficient
1739  ! im=2,3,4: moments of the volume scattering phase function
1740 
1741  !--- Rayleigh scattering
1742  !$acc parallel loop collapse(2)
1743  loop_inner
1744  !$acc loop collapse(2)
1745  do im = 1, mstrn_nstream*2+2
1746  do k = 1, rd_kmax
1747  dp_p = rhodz(k,i,j) * grav / pstd
1748  length = rayleigh(iw) * dp_p
1749 
1750  optparam(k,im,i_cloud ,i,j) = qmol(im,iw) * length
1751  optparam(k,im,i_clearsky,i,j) = qmol(im,iw) * length
1752  enddo
1753  enddo
1754  loop_end_inner
1755  !$acc end parallel
1756 
1757  !--- Cloud scattering
1758  !$acc parallel loop collapse(2)
1759  loop_inner
1760  do k = 1, rd_kmax
1761  taucld_067u(k,i,j) = 0.0_rp
1762  enddo
1763  loop_end_inner
1764  !$acc end parallel
1765  !$acc kernels
1766  loop_inner
1767  !$acc loop seq
1768  do iaero = hydro_str, hydro_end
1769  iptype = aero2ptype(iaero)
1770 
1771  do im = 1, mstrn_nstream*2+2
1772  do k = 1, rd_kmax
1773  ir = indexr(k,iaero,i,j)
1774  length = aerosol_conc(k,i,j,iaero) * ppm * dz_std(k,i,j)
1775 
1776  q_fit = q(ir ,iptype,im,iw) * ( 1.0_rp-factr(k,iaero,i,j) ) &
1777  + q(ir+1,iptype,im,iw) * ( factr(k,iaero,i,j) )
1778  q_fit = max( q_fit, 0.0_rp )
1779 
1780  optparam(k,im,i_cloud,i,j) = optparam(k,im,i_cloud,i,j) + q_fit * length
1781  enddo
1782  enddo
1783 
1784  if ( waveh(iw) <= 1.493e+4_rp .AND. 1.493e+4_rp < waveh(iw+1) ) then ! 0.67 micron
1785  do k = 1, rd_kmax
1786  ir = indexr(k,iaero,i,j)
1787  length = aerosol_conc(k,i,j,iaero) * ppm * dz_std(k,i,j)
1788 
1789  q_fit = q(ir ,iptype,1,iw) * ( 1.0_rp-factr(k,iaero,i,j) ) &
1790  + q(ir+1,iptype,1,iw) * ( factr(k,iaero,i,j) )
1791  q_fit = max( q_fit, 0.0_rp )
1792 
1793  taucld_067u(k,i,j) = taucld_067u(k,i,j) + q_fit * length
1794  end do
1795  end if
1796 
1797  enddo
1798  loop_end_inner
1799  !$acc end kernels
1800 
1801  !--- Aerosol scattering
1802  !$acc kernels
1803  loop_inner
1804  !$acc loop seq
1805  do iaero = aero_str, aero_end
1806  iptype = aero2ptype(iaero)
1807 
1808  do im = 1, mstrn_nstream*2+2
1809  do k = 1, rd_kmax
1810  ir = indexr(k,iaero,i,j)
1811  length = aerosol_conc(k,i,j,iaero) * ppm * dz_std(k,i,j)
1812 
1813  q_fit = q(ir ,iptype,im,iw) * ( 1.0_rp-factr(k,iaero,i,j) ) &
1814  + q(ir+1,iptype,im,iw) * ( factr(k,iaero,i,j) )
1815  q_fit = max( q_fit, 0.0_rp )
1816 
1817  optparam(k,im,i_cloud ,i,j) = optparam(k,im,i_cloud ,i,j) + q_fit * length
1818  optparam(k,im,i_clearsky,i,j) = optparam(k,im,i_clearsky,i,j) + q_fit * length
1819  enddo
1820  enddo
1821  enddo
1822  loop_end_inner
1823  !$acc end kernels
1824 
1825  !$acc kernels
1826  !$acc loop independent collapse(2) &
1827  !$acc private(tauPR,omgPR,g,b,b_sfc,fsol_rgn,bbar,bbarh,tau,omg, &
1828  !$acc flux,flux_direct,emisCLD, &
1829  !$acc Tdir0,R0,T0,Em_LW,Em_SW,Ep_LW,Ep_SW,cf,tau_bar_sol,R,T,Em,Ep,Em0, &
1830  !$acc R12mns,R12pls,E12mns,E12pls )
1831  loop_inner
1832 
1833  !$acc loop collapse(2)
1834  do icloud = 1, mstrn_ncloud
1835  do k = 1, rd_kmax
1836  taupr(k,icloud) = optparam(k,1,icloud,i,j)
1837  omgpr(k,icloud) = optparam(k,1,icloud,i,j) - optparam(k,2,icloud,i,j)
1838 
1839  !--- g
1840  zerosw = 0.5_rp - sign(0.5_rp,omgpr(k,icloud)-rd_eps)
1841 
1842  g(k,0,icloud) = 1.0_rp
1843  g(k,1,icloud) = optparam(k,3,icloud,i,j) * ( 1.0_rp-zerosw ) / ( omgpr(k,icloud)+zerosw )
1844  g(k,2,icloud) = optparam(k,4,icloud,i,j) * ( 1.0_rp-zerosw ) / ( omgpr(k,icloud)+zerosw )
1845  !g(k,1,icloud) = max( optparam(k,3,icloud,i,j) * ( 1.0_RP-zerosw ) / ( omgPR(k,icloud)+zerosw ), 0.0_RP )
1846  !g(k,2,icloud) = max( optparam(k,4,icloud,i,j) * ( 1.0_RP-zerosw ) / ( omgPR(k,icloud)+zerosw ), 0.0_RP )
1847  enddo
1848  enddo
1849 
1850  if ( irgn == i_sw ) then ! solar
1851  b(:,:,:,:) = 0.0_rp
1852 
1853  b_sfc = 0.0_rp
1854  fsol_rgn = fsol(iw) / fsol_tot * solins(i,j)
1855 
1856  elseif( irgn == i_lw ) then ! IR
1857  !--- set planck functions
1858  wl = 10000.0_rp / sqrt( waveh(iw) * waveh(iw+1) )
1859 
1860  ! from temp at cell center
1861  do k = 1, rd_kmax
1862  beta = 0.0_rp
1863  !$acc loop seq
1864  do iplk = mstrn_nfitplk, 1, -1
1865  beta = beta / ( wl*temp(k,i,j) ) + fitplk(iplk,iw)
1866  enddo
1867  bbar(k) = exp(-beta) * temp(k,i,j) / (wl*wl)
1868  enddo
1869 
1870  ! from temp at cell wall
1871  do k = 1, rd_kmax+1
1872  beta = 0.0_rp
1873  !$acc loop seq
1874  do iplk = mstrn_nfitplk, 1, -1
1875  beta = beta / ( wl*temph(k,i,j) ) + fitplk(iplk,iw)
1876  enddo
1877  bbarh(k) = exp(-beta) * temph(k,i,j) / (wl*wl)
1878  enddo
1879 
1880  ! from temp_sfc
1881  beta = 0.0_rp
1882  !$acc loop seq
1883  do iplk = mstrn_nfitplk, 1, -1
1884  beta = beta / ( wl*temp_sfc(i,j) ) + fitplk(iplk,iw)
1885  enddo
1886 
1887  b_sfc = exp(-beta) * temp_sfc(i,j) / (wl*wl)
1888  fsol_rgn = 0.0_rp
1889 
1890  endif
1891 
1892  !--- total tau & omega
1893  !$acc loop collapse(3)
1894  do ich = 1, chmax
1895  do icloud = 1, mstrn_ncloud
1896  do k = 1, rd_kmax
1897  tau(k,icloud,ich) = taugas(k,ich,i,j) + taupr(k,icloud)
1898  !tau(k,icloud,ich) = max( tauGAS(k,ich,i,j) + tauPR(k,icloud), 0.0_RP )
1899 
1900  zerosw = 0.5_rp - sign( 0.5_rp, tau(k,icloud,ich)-rd_eps ) ! if tau < EPS, zerosw = 1
1901 
1902  omg(k,icloud,ich) = ( 1.0_rp-zerosw ) * omgpr(k,icloud) / ( tau(k,icloud,ich)-zerosw ) &
1903  + ( zerosw ) * 1.0_rp
1904 
1905  !omg(k,icloud,ich) = min( max( omg(k,icloud,ich), 0.0_RP ), 1.0_RP )
1906  enddo
1907  enddo
1908  enddo
1909 
1910  if( irgn == i_lw ) then ! IR
1911  !$acc loop collapse(3)
1912  do ich = 1, chmax
1913  do icloud = 1, mstrn_ncloud
1914  do k = 1, rd_kmax
1915  zerosw = 0.5_rp - sign( 0.5_rp, tau(k,icloud,ich)-rd_eps ) ! if tau < EPS, zerosw = 1
1916 
1917  b(k,0,icloud,ich) = bbarh(k)
1918  b(k,1,icloud,ich) = ( 1.0_rp-zerosw ) &
1919  * ( - bbarh(k+1) &
1920  + 4.0_rp * bbar(k ) &
1921  - 3.0_rp * bbarh(k ) &
1922  ) / ( tau(k,icloud,ich)-zerosw )
1923  b(k,2,icloud,ich) = ( 1.0_rp-zerosw ) &
1924  * ( + bbarh(k+1) &
1925  - 2.0_rp * bbar(k ) &
1926  + bbarh(k ) &
1927  ) / ( tau(k,icloud,ich)**2-zerosw ) * 2.0_rp
1928  enddo
1929  enddo
1930  enddo
1931 
1932  endif ! solar/IR switch
1933 
1934 
1935  ! two-stream transfer
1936  !call PROF_rapstart('RD_MSTRN_twst', 3)
1937  call rd_mstrn_two_stream( ia, is, ie, ja, js, je, &
1938  rd_kmax, chmax, & ! [IN]
1939  cossza(i,j), & ! [IN]
1940  fsol_rgn, irgn, & ! [IN]
1941  tau(:,:,:), omg(:,:,:), & ! [IN]
1942  g(:,:,:), b(:,:,:,:), b_sfc, & ! [IN]
1943  albedo_sfc(i,j,i_r_direct, irgn_alb), & ! [IN]
1944  albedo_sfc(i,j,i_r_diffuse,irgn_alb), & ! [IN]
1945  cldfrac(:,i,j), & ! [IN]
1946  flux(:,:,:,:), flux_direct(:,:,:), & ! [OUT]
1947  emiscld(:,:), & ! [OUT]
1948  tdir0(:,:,:), r0(:,:,:), t0(:,:,:), & ! [WORK]
1949  em_lw(:,:,:), em_sw(:,:,:), & ! [WORK]
1950  ep_lw(:,:,:), ep_sw(:,:,:), & ! [WORK]
1951  cf(:,:), tau_bar_sol(:,:), & ! [WORK]
1952  r(:,:), t(:,:), & ! [WORK]
1953  em(:,:), ep(:,:), em0(:), & ! [WORK]
1954  r12mns(:,:), r12pls(:,:), & ! [WORK]
1955  e12mns(:,:), e12pls(:,:) ) ! [WORK]
1956  !call PROF_rapend ('RD_MSTRN_twst', 3)
1957 
1958  !$acc loop seq
1959  do ich = 1, chmax
1960  !$acc loop collapse(2)
1961  do icloud = 1, mstrn_ncloud
1962  do k = 1, rd_kmax+1
1963  rflux(k,i,j,irgn,i_up,icloud) = rflux(k,i,j,irgn,i_up,icloud) + flux(k,i_up,icloud,ich) * wgtch(ich,iw)
1964  rflux(k,i,j,irgn,i_dn,icloud) = rflux(k,i,j,irgn,i_dn,icloud) + flux(k,i_dn,icloud,ich) * wgtch(ich,iw)
1965  enddo
1966  enddo
1967  enddo
1968 
1969  !$acc loop seq
1970  do ich = 1, chmax
1971  rflux_sfc_dn(i,j,i_r_direct ,irgn_alb) = rflux_sfc_dn(i,j,i_r_direct ,irgn_alb) &
1972  + ( flux_direct(rd_kmax+1,i_cloud,ich) ) * wgtch(ich,iw)
1973  rflux_sfc_dn(i,j,i_r_diffuse,irgn_alb) = rflux_sfc_dn(i,j,i_r_diffuse,irgn_alb) &
1974  + ( flux(rd_kmax+1,i_dn,i_cloud,ich) - flux_direct(rd_kmax+1,i_cloud,ich) ) * wgtch(ich,iw)
1975  end do
1976 
1977  if ( waveh(iw) <= 952.0_rp .AND. 952.0_rp < waveh(iw+1) ) then ! 10.5 micron
1978  ! 10.5 micron emissivity for resolved clouds
1979 
1980  do k = 1, rd_kmax
1981  emiscld_105u(k,i,j) = 0.0_rp
1982  end do
1983  !$acc loop seq
1984  do ich = 1, chmax
1985  do k = 1, rd_kmax
1986  emiscld_105u(k,i,j) = emiscld_105u(k,i,j) + emiscld(k,ich) * wgtch(ich,iw)
1987  enddo
1988  end do
1989  endif
1990 
1991  loop_end_inner
1992  !$acc end kernels
1993 
1994  enddo ! IW loop
1995 
1996  !$acc end data
1997 
1998  return
1999  end subroutine rd_mstrn_dtrn3
2000 
2001  !-----------------------------------------------------------------------------
2003 !OCL SERIAL
2004  subroutine rd_mstrn_two_stream( &
2005  IA, IS, IE, JA, JS, JE, &
2006  RD_KMAX, chmax, &
2007  cosSZA, &
2008  fsol, irgn, &
2009  tau, omg, &
2010  g, b, b_sfc, &
2011  albedo_sfc_direct, &
2012  albedo_sfc_diffuse, &
2013  cldfrac, &
2014  flux, flux_direct, &
2015  emisCLD, &
2016  Tdir0, R0, T0, &
2017  Em_LW, Em_SW, &
2018  Ep_LW, Ep_SW, &
2019  cf, tau_bar_sol, &
2020  R, T, &
2021  Em, Ep, Em0, &
2022  R12mns, R12pls, &
2023  E12mns, E12pls )
2024  !$acc routine vector
2025  use scale_const, only: &
2026  pi => const_pi, &
2027  eps => const_eps, &
2028  eps1 => const_eps1
2029  implicit none
2030 
2031  integer, intent(in) :: ia, is, ie, ja, js, je
2032  integer, intent(in) :: rd_kmax
2033  integer, intent(in), value :: chmax
2034  real(rp), intent(in) :: cossza ! cos(SZA) = mu0
2035  real(rp), intent(in), value :: fsol ! solar radiation intensity
2036  integer, intent(in), value :: irgn ! 1:LW 2:SW
2037  real(rp), intent(in) :: tau (rd_kmax, mstrn_ncloud,mstrn_ch_limit) ! total optical thickness (clear-sky/cloud)
2038  real(rp), intent(in) :: omg (rd_kmax, mstrn_ncloud,mstrn_ch_limit) ! single scattering albedo (clear-sky/cloud)
2039  real(rp), intent(in) :: g (rd_kmax,0:2,mstrn_ncloud ) ! two-stream approximation factors (clear-sky/cloud)
2040  real(rp), intent(in) :: b (rd_kmax,0:2,mstrn_ncloud,mstrn_ch_limit) ! planck expansion coefficients (clear-sky/cloud)
2041  real(rp), intent(in), value :: b_sfc ! planck function at surface
2042  real(rp), intent(in) :: albedo_sfc_direct ! surface albedo (DIRECT)
2043  real(rp), intent(in) :: albedo_sfc_diffuse ! surface albedo (DIFFUSE)
2044  real(rp), intent(in) :: cldfrac (rd_kmax) ! cloud fraction
2045  real(rp), intent(out) :: flux (rd_kmax+1,2,mstrn_ncloud,mstrn_ch_limit) ! upward(sfc->TOA)/downward(TOA->sfc) flux (clear-sky/cloud)
2046  real(rp), intent(out) :: flux_direct(rd_kmax+1, mstrn_ncloud,mstrn_ch_limit) ! downward(TOA->sfc) flux, solar direct (clear-sky/cloud)
2047  real(rp), intent(out) :: emiscld (rd_kmax, mstrn_ch_limit) ! cloud emissivity factor (cloud)
2048 
2049  ! work
2050  ! main factors
2051  real(rp), intent(out) :: tdir0(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! transmission factor for solar direct (clear-sky/cloud)
2052  real(rp), intent(out) :: r0 (rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! reflection factor (clear-sky/cloud)
2053  real(rp), intent(out) :: t0 (rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! transmission factor (clear-sky/cloud)
2054  real(rp), intent(out) :: em_lw(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! thermal source (sfc->TOA) (clear-sky/cloud)
2055  real(rp), intent(out) :: em_sw(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! solar source (sfc->TOA) (clear-sky/cloud)
2056  real(rp), intent(out) :: ep_lw(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! thermal source (TOA->sfc) (clear-sky/cloud)
2057  real(rp), intent(out) :: ep_sw(rd_kmax,mstrn_ncloud,mstrn_ch_limit) ! solar source (TOA->sfc) (clear-sky/cloud)
2058  ! Averaged factors, considering cloud overwrap
2059  real(rp), intent(out) :: cf (rd_kmax ,mstrn_ncloud*mstrn_ch_limit) ! cloud fraction
2060  real(rp), intent(out) :: tau_bar_sol(rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! solar insolation through accumulated optical thickness at each layer
2061  real(rp), intent(out) :: r (rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! reflection factor
2062  real(rp), intent(out) :: t (rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! transmission factor
2063  real(rp), intent(out) :: em (rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! source (sfc->TOA)
2064  real(rp), intent(out) :: ep (rd_kmax+1,mstrn_ncloud*mstrn_ch_limit) ! source (TOA->sfc)
2065  real(rp), intent(out) :: em0 (mstrn_ncloud)
2066  ! Doubling-Adding
2067  real(rp), intent(out) :: r12mns(mstrn_ncloud*mstrn_ch_limit,rd_kmax+1) ! reflection factor in doubling method
2068  real(rp), intent(out) :: r12pls(mstrn_ncloud*mstrn_ch_limit,rd_kmax+1) ! reflection factor in doubling method
2069  real(rp), intent(out) :: e12mns(mstrn_ncloud*mstrn_ch_limit,rd_kmax+1) ! source function in doubling method
2070  real(rp), intent(out) :: e12pls(mstrn_ncloud*mstrn_ch_limit,rd_kmax+1) ! source function in doubling method
2071 
2072 
2073  ! parameters with two-stream truncation
2074  real(rp) :: tau_new ! optical thickness : two-stream truncation
2075  real(rp) :: omg_new ! single scattering albedo : two-stream truncation
2076  real(rp) :: g_new ! asymmetric factor : two-stream truncation
2077  real(rp) :: b_new0 ! planck function : two-stream truncation
2078  real(rp) :: b_new1
2079  real(rp) :: b_new2
2080  real(rp) :: c0
2081  real(rp) :: c1
2082  real(rp) :: c2
2083  real(rp) :: pmns, ppls ! Phase function : two-stream truncation
2084  real(rp) :: smns, spls ! Source function : two-stream truncation
2085 
2086  ! working
2087  real(rp) :: x, y ! X-, X+
2088  real(rp) :: lamda ! eigenvalue of X-, X+
2089  real(rp) :: e
2090  real(rp) :: apls_mns, bpls_mns ! A+/A-, B+/B-
2091  real(dp) :: v0mns, v0pls ! V0-, V0+
2092  real(dp) :: v1mns, v1pls ! V1-, V1+
2093  real(dp) :: dmns0, dmns1, dmns2 ! D0-, D1-, D2-
2094  real(dp) :: dpls0, dpls1, dpls2 ! D0+, D1+, D2+
2095  real(rp) :: sigmns, sigpls ! sigma-, sigma+
2096  real(rp) :: qgamma ! Q * gamma
2097  real(rp) :: zerosw, tmp
2098 
2099  real(rp) :: tdir ! transmission factor for solar direct
2100 
2101  real(rp) :: umns, upls ! flux intensity
2102 
2103  real(rp) :: factor
2104  real(rp) :: wmns_irgn, m_irgn, w_irgn, wpls_irgn, wscale_irgn
2105 
2106  integer, parameter :: i_sfc2toa = 1
2107  integer, parameter :: i_toa2sfc = 2
2108 
2109  real(rp) :: sw
2110  integer :: k, kk, icloud, ich, ic
2111  !---------------------------------------------------------------------------
2112 
2113  m_irgn = m(irgn)
2114  w_irgn = w(irgn)
2115  wmns_irgn = wmns(irgn)
2116  wpls_irgn = wpls(irgn)
2117  wscale_irgn = wscale(irgn)
2118 
2119  !$acc loop collapse(3)
2120  do ich = 1, chmax
2121  do icloud = 1, mstrn_ncloud
2122 !OCL LOOP_FISSION_TARGET(LS)
2123  do k = 1, rd_kmax
2124 
2125  !---< two-stream truncation >---
2126  tau_new = ( 1.0_rp - omg(k,icloud,ich)*g(k,2,icloud) ) * tau(k,icloud,ich)
2127 
2128  omg_new = ( 1.0_rp - g(k,2,icloud) ) / ( 1.0_rp - omg(k,icloud,ich)*g(k,2,icloud) ) * omg(k,icloud,ich)
2129  omg_new = min( omg_new, eps1 )
2130 
2131  g_new = ( g(k,1,icloud) - g(k,2,icloud) ) / ( 1.0_rp - g(k,2,icloud) )
2132 
2133 #if defined(NVIDIA) || defined(SX)
2134  tdir0(k,icloud,ich) = exp( -min( tau_new/cossza, 1.e+3_rp ) ) ! apply exp limiter
2135 #else
2136  tdir0(k,icloud,ich) = exp(-tau_new/cossza)
2137 #endif
2138 
2139  factor = ( 1.0_rp - omg(k,icloud,ich)*g(k,2,icloud) )
2140  b_new0 = b(k,0,icloud,ich)
2141  b_new1 = b(k,1,icloud,ich) / factor
2142  b_new2 = b(k,2,icloud,ich) / (factor*factor)
2143  c0 = wmns_irgn * 2.0_rp * pi * ( 1.0_rp - omg_new ) * b_new0
2144  c1 = wmns_irgn * 2.0_rp * pi * ( 1.0_rp - omg_new ) * b_new1
2145  c2 = wmns_irgn * 2.0_rp * pi * ( 1.0_rp - omg_new ) * b_new2
2146 
2147  !--- P+, P-
2148  pmns = omg_new * 0.5_rp * ( 1.0_rp - 3.0_rp * g_new * m_irgn*m_irgn )
2149  ppls = omg_new * 0.5_rp * ( 1.0_rp + 3.0_rp * g_new * m_irgn*m_irgn )
2150 
2151  !--- S+, S-
2152  smns = omg_new * 0.5_rp * ( 1.0_rp - 3.0_rp * g_new * m_irgn*cossza )
2153  spls = omg_new * 0.5_rp * ( 1.0_rp + 3.0_rp * g_new * m_irgn*cossza )
2154 
2155  !---< calculate R, T, e+, e- >---
2156  sw = 0.5_rp + sign(0.5_rp,tau_new-rd_eps)
2157 
2158  !--- X, Y
2159  x = ( 1.0_rp - w_irgn * ( ppls - pmns ) ) / m_irgn
2160  y = ( 1.0_rp - w_irgn * ( ppls + pmns ) ) / m_irgn
2161  !X = max( ( 1.0_RP - W_irgn * ( Ppls - Pmns ) ) / M_irgn, 1.E-30 )
2162  !Y = max( ( 1.0_RP - W_irgn * ( Ppls + Pmns ) ) / M_irgn, 1.E-30 )
2163  lamda = sqrt(x*y)
2164 #if defined(NVIDIA) || defined(SX)
2165  e = exp( -min( lamda*tau_new, 1.e+3_rp ) ) ! apply exp limiter
2166 #else
2167  e = exp(-lamda*tau_new)
2168 #endif
2169 
2170  !--- A+/A-, B+/B-
2171  apls_mns = ( x * ( 1.0_dp+e ) - lamda * ( 1.0_dp-e ) ) &
2172  / ( x * ( 1.0_dp+e ) + lamda * ( 1.0_dp-e ) )
2173  bpls_mns = ( x * ( 1.0_dp-e ) - lamda * ( 1.0_dp+e ) ) &
2174  / ( x * ( 1.0_dp-e ) + lamda * ( 1.0_dp+e ) )
2175 
2176  !--- R, T
2177  r0(k,icloud,ich) = ( sw ) * 0.5_rp * ( apls_mns + bpls_mns ) &
2178  + ( 1.0_rp-sw ) * ( tau_new * ( pmns ) / m_irgn )
2179  t0(k,icloud,ich) = ( sw ) * 0.5_rp * ( apls_mns - bpls_mns ) &
2180  + ( 1.0_rp-sw ) * ( 1.0_rp - tau_new * ( 1.0_rp - ppls ) / m_irgn )
2181 
2182  !--- thermal source
2183  dmns0 = c0 / y + 2.0_rp * c2 / (x*y*y) + c1 / (x*y)
2184  dpls0 = c0 / y + 2.0_rp * c2 / (x*y*y) - c1 / (x*y)
2185  dmns1 = c1 / y + 2.0_rp * c2 / (x*y)
2186  dpls1 = c1 / y - 2.0_rp * c2 / (x*y)
2187  dmns2 = c2 / y
2188  dpls2 = c2 / y
2189 
2190  v0mns = dmns0
2191  v0pls = dpls0
2192  v1mns = dmns0 + dmns1*tau_new + dmns2*tau_new*tau_new
2193  v1pls = dpls0 + dpls1*tau_new + dpls2*tau_new*tau_new
2194 
2195  em_lw(k,icloud,ich) = ( sw ) * ( v0mns - r0(k,icloud,ich) * v0pls - t0(k,icloud,ich) * v1mns ) &
2196  + ( 1.0_rp-sw ) * 0.5_rp * tau_new * ( 2.0_rp*c0 + c1*tau_new + c2*tau_new*tau_new )
2197  ep_lw(k,icloud,ich) = ( sw ) * ( v1pls - t0(k,icloud,ich) * v0pls - r0(k,icloud,ich) * v1mns ) &
2198  + ( 1.0_rp-sw ) * 0.5_rp * tau_new * ( 2.0_rp*c0 + c1*tau_new + c2*tau_new*tau_new )
2199 
2200  !--- solar source
2201  sigmns = wmns_irgn * ( spls - smns )
2202  sigpls = wmns_irgn * ( spls + smns )
2203 
2204  tmp = x*y*cossza - 1.0/cossza
2205  zerosw = 1.0_rp - sign(1.0_rp,abs(tmp)-eps) ! if abs(tmp)<EPS then 2, otherwise 0
2206  qgamma = ( sigpls*x*cossza + sigmns ) / ( tmp + zerosw*eps )
2207 
2208  v0pls = 0.5_rp * ( ( 1.0_rp + 1.0_rp/(x*cossza) ) * qgamma + sigmns / x )
2209  v0mns = 0.5_rp * ( ( 1.0_rp - 1.0_rp/(x*cossza) ) * qgamma - sigmns / x )
2210 
2211  v1pls = v0pls * tdir0(k,icloud,ich)
2212  v1mns = v0mns * tdir0(k,icloud,ich)
2213 
2214  em_sw(k,icloud,ich) = ( sw ) * ( v0mns - r0(k,icloud,ich) * v0pls - t0(k,icloud,ich) * v1mns ) &
2215  + ( 1.0_rp-sw ) * wmns_irgn * smns * tau_new * sqrt( tdir0(k,icloud,ich) )
2216  ep_sw(k,icloud,ich) = ( sw ) * ( v1pls - t0(k,icloud,ich) * v0pls - r0(k,icloud,ich) * v1mns ) &
2217  + ( 1.0_rp-sw ) * wmns_irgn * spls * tau_new * sqrt( tdir0(k,icloud,ich) )
2218  enddo
2219  enddo
2220  enddo
2221 
2222 
2223  !---< consider partial cloud layer: semi-random over-wrapping >---
2224 
2225  do icloud = 1, mstrn_ncloud
2226  if ( icloud == i_clearsky ) then
2227  !$acc loop collapse(2)
2228  do ich = 1, chmax
2229  do k = 1, rd_kmax
2230  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2231  cf(k,ic) = 0.0_rp
2232  enddo
2233  enddo
2234  else
2235  !$acc loop independent collapse(2)
2236  do ich = 1, chmax
2237  do k = 1, rd_kmax
2238  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2239  cf(k,ic) = cldfrac(k)
2240  enddo
2241  enddo
2242  endif
2243  enddo
2244 
2245 
2246  !$acc loop collapse(2)
2247  do ich = 1, chmax
2248  do icloud = 1, mstrn_ncloud
2249  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2250  tau_bar_sol(1,ic) = fsol
2251  !$acc loop seq
2252  do k = 2, rd_kmax+1 ! k-recurrence
2253  tdir = ( cf(k-1,ic) ) * tdir0(k-1,i_cloud ,ich) &
2254  + ( 1.0_rp-cf(k-1,ic) ) * tdir0(k-1,i_clearsky,ich)
2255  tau_bar_sol(k,ic) = tau_bar_sol(k-1,ic) * tdir
2256  enddo
2257  enddo
2258  enddo
2259 
2260  !$acc loop collapse(3)
2261  do ich = 1, chmax
2262  do icloud = 1, mstrn_ncloud
2263  do k = 1, rd_kmax+1
2264  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2265  flux_direct(k,icloud,ich) = cossza * tau_bar_sol(k,ic)
2266  enddo
2267  enddo
2268  enddo
2269 
2270  !$acc loop collapse(3)
2271  do ich = 1, chmax
2272  do icloud = 1, mstrn_ncloud
2273  do k = 1, rd_kmax
2274  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2275  em(k,ic) = ( cf(k,ic) ) * ( em_lw(k,i_cloud, ich) &
2276  + em_sw(k,i_cloud, ich) * tau_bar_sol(k,ic) ) &
2277  + ( 1.0_rp-cf(k,ic) ) * ( em_lw(k,i_clearsky,ich) &
2278  + em_sw(k,i_clearsky,ich) * tau_bar_sol(k,ic) )
2279 
2280  ep(k,ic) = ( cf(k,ic) ) * ( ep_lw(k,i_cloud, ich) &
2281  + ep_sw(k,i_cloud, ich) * tau_bar_sol(k,ic) ) &
2282  + ( 1.0_rp-cf(k,ic) ) * ( ep_lw(k,i_clearsky,ich) &
2283  + ep_sw(k,i_clearsky,ich) * tau_bar_sol(k,ic) )
2284 
2285  r(k,ic) = ( cf(k,ic) ) * r0(k,i_cloud ,ich) &
2286  + ( 1.0_rp-cf(k,ic) ) * r0(k,i_clearsky,ich)
2287  t(k,ic) = ( cf(k,ic) ) * t0(k,i_cloud ,ich) &
2288  + ( 1.0_rp-cf(k,ic) ) * t0(k,i_clearsky,ich)
2289  enddo
2290  enddo
2291  enddo
2292 
2293  !---< Adding-Doubling method >---
2294  ! [note] TOA->Surface is positive direction. "pls" means upper to lower altitude.
2295  !$acc loop collapse(2)
2296  do ich = 1, chmax
2297  do icloud = 1, mstrn_ncloud
2298  ! at lambert surface
2299  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2300  r(rd_kmax+1,ic) = ( cf(rd_kmax,ic) ) * albedo_sfc_diffuse &
2301  + ( 1.0_rp-cf(rd_kmax,ic) ) * albedo_sfc_diffuse
2302  t(rd_kmax+1,ic) = 0.0_rp
2303 
2304  ! currently, Em0(I_Cloud) and Em0(I_ClearSky) is the same
2305  em0(i_cloud ) = wpls_irgn * ( flux_direct(rd_kmax+1,icloud,ich) * albedo_sfc_direct / (w_irgn*m_irgn) &
2306  + 2.0_rp * pi * ( 1.0_rp-albedo_sfc_diffuse ) * b_sfc )
2307  em0(i_clearsky) = wpls_irgn * ( flux_direct(rd_kmax+1,icloud,ich) * albedo_sfc_direct / (w_irgn*m_irgn) &
2308  + 2.0_rp * pi * ( 1.0_rp-albedo_sfc_diffuse ) * b_sfc )
2309 
2310  em(rd_kmax+1,ic) = ( cf(rd_kmax,ic) ) * em0(i_cloud ) &
2311  + ( 1.0_rp-cf(rd_kmax,ic) ) * em0(i_clearsky)
2312  ep(rd_kmax+1,ic) = 0.0_rp
2313 
2314 
2315  r12pls(ic,rd_kmax+1) = r(rd_kmax+1,ic)
2316  e12mns(ic,rd_kmax+1) = em(rd_kmax+1,ic)
2317  r12mns(ic,1) = r(1,ic)
2318  e12pls(ic,1) = ep(1,ic)
2319  enddo
2320  enddo
2321 
2322  ! cloud emissivity
2323  ich = chmax / 2 + 1
2324  icloud = 2
2325  !$acc loop collapse(2)
2326  do ich = 1, chmax
2327  do k = 1, rd_kmax
2328  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2329  emiscld(k,ich) = 1.0_rp - r(k,ic) - t(k,ic)
2330  enddo
2331  enddo
2332 
2333  !$acc loop seq
2334  do kk = 1, rd_kmax
2335 !OCL ITERATIONS(max=20)
2336  do ic = 1, chmax * mstrn_ncloud
2337  ! adding: surface to TOA
2338  k = rd_kmax - kk + 1
2339  r12pls(ic,k) = r(k,ic) + t(k,ic) / ( 1.0_rp - r12pls(ic,k+1) * r(k,ic) ) &
2340  * ( r12pls(ic,k+1) * t(k,ic) )
2341  e12mns(ic,k) = em(k,ic) + t(k,ic) / ( 1.0_rp - r12pls(ic,k+1) * r(k,ic) ) &
2342  * ( r12pls(ic,k+1) * ep(k,ic) + e12mns(ic,k+1) )
2343  ! adding: TOA to surface
2344  k = kk + 1
2345  r12mns(ic,k) = r(k,ic) + t(k,ic) / ( 1.0_rp - r12mns(ic,k-1) * r(k,ic) ) &
2346  * ( r12mns(ic,k-1) * t(k,ic) )
2347  e12pls(ic,k) = ep(k,ic) + t(k,ic) / ( 1.0_rp - r12mns(ic,k-1) * r(k,ic) ) &
2348  * ( r12mns(ic,k-1) * em(k,ic) + e12pls(ic,k-1) )
2349  enddo
2350  enddo
2351 
2352  !--- radiative flux at cell wall
2353  !$acc loop collapse(2)
2354  do ich = 1, chmax
2355  do icloud = 1, mstrn_ncloud
2356  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2357  ! TOA boundary
2358  upls = 0.0_rp
2359  umns = e12mns(ic,1) + r12pls(ic,1) * upls
2360  flux(1,i_up,icloud,ich) = wscale_irgn * umns
2361  flux(1,i_dn,icloud,ich) = wscale_irgn * upls + flux_direct(1,icloud,ich)
2362  enddo
2363  enddo
2364 
2365  !$acc loop collapse(3)
2366  do ich = 1, chmax
2367  do icloud = 1, mstrn_ncloud
2368  do k = 2, rd_kmax+1
2369  ic = icloud + ( ich - 1 ) * mstrn_ncloud
2370  upls = ( e12pls(ic,k-1) + r12mns(ic,k-1)*e12mns(ic,k) ) / ( 1.0_rp - r12mns(ic,k-1)*r12pls(ic,k) )
2371  umns = e12mns(ic,k) + r12pls(ic,k) * upls
2372  flux(k,i_up,icloud,ich) = wscale_irgn * umns
2373  flux(k,i_dn,icloud,ich) = wscale_irgn * upls + flux_direct(k,icloud,ich)
2374  enddo
2375  enddo
2376  enddo
2377 
2378  return
2379  end subroutine rd_mstrn_two_stream
2380 
2381 end module scale_atmos_phy_rd_mstrnx
scale_const::const_grav
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:49
scale_cpl_sfc_index::n_rad_dir
integer, parameter, public n_rad_dir
Definition: scale_cpl_sfc_index.F90:36
scale_prc::prc_abort
subroutine, public prc_abort
Abort Process.
Definition: scale_prc.F90:350
scale_cpl_sfc_index::i_r_direct
integer, parameter, public i_r_direct
Definition: scale_cpl_sfc_index.F90:37
scale_atmos_phy_rd_profile::atmos_phy_rd_profile_use_climatology
logical, public atmos_phy_rd_profile_use_climatology
use climatology?
Definition: scale_atmos_phy_rd_profile.F90:39
scale_cpl_sfc_index::i_r_diffuse
integer, parameter, public i_r_diffuse
Definition: scale_cpl_sfc_index.F90:38
scale_precision
module PRECISION
Definition: scale_precision.F90:14
scale_atmos_aerosol::n_ae
integer, parameter, public n_ae
Definition: scale_atmos_aerosol.F90:33
scale_atmos_phy_rd_common::i_dn
integer, parameter, public i_dn
Definition: scale_atmos_phy_rd_common.F90:35
scale_atmos_phy_rd_common::i_lw
integer, parameter, public i_lw
Definition: scale_atmos_phy_rd_common.F90:37
scale_cpl_sfc_index::i_r_ir
integer, parameter, public i_r_ir
Definition: scale_cpl_sfc_index.F90:29
scale_const::const_eps
real(rp), public const_eps
small number
Definition: scale_const.F90:35
scale_prof::prof_rapstart
subroutine, public prof_rapstart(rapname_base, level, disable_barrier)
Start raptime.
Definition: scale_prof.F90:174
scale_atmos_phy_rd_common::i_up
integer, parameter, public i_up
Definition: scale_atmos_phy_rd_common.F90:34
scale_atmos_hydrometeor
module atmosphere / hydrometeor
Definition: scale_atmos_hydrometeor.F90:12
scale_atmos_phy_rd_mstrnx::atmos_phy_rd_mstrnx_finalize
subroutine, public atmos_phy_rd_mstrnx_finalize
finalize
Definition: scale_atmos_phy_rd_mstrnx.F90:441
scale_io::io_get_available_fid
integer function, public io_get_available_fid()
search & get available file ID
Definition: scale_io.F90:373
scale_atmos_grid_cartesc_real
module Atmosphere GRID CartesC Real(real space)
Definition: scale_atmos_grid_cartesC_real.F90:11
scale_atmos_grid_cartesc_real::atmos_grid_cartesc_real_basepoint_lat
real(rp), public atmos_grid_cartesc_real_basepoint_lat
position of base point in real world [rad,-pi,pi]
Definition: scale_atmos_grid_cartesC_real.F90:37
scale_prc
module PROCESS
Definition: scale_prc.F90:11
scale_atmos_hydrometeor::i_hi
integer, parameter, public i_hi
ice water cloud
Definition: scale_atmos_hydrometeor.F90:99
scale_io::io_get_fname
subroutine, public io_get_fname(outstr, instr, rank, ext, len)
generate process specific filename
Definition: scale_io.F90:421
scale_io
module STDIO
Definition: scale_io.F90:10
scale_const::const_pstd
real(rp), public const_pstd
standard pressure [Pa]
Definition: scale_const.F90:96
scale_cpl_sfc_index::i_r_nir
integer, parameter, public i_r_nir
Definition: scale_cpl_sfc_index.F90:30
scale_const
module CONSTANT
Definition: scale_const.F90:11
scale_const::const_mdry
real(rp), public const_mdry
mass weight (dry air) [g/mol]
Definition: scale_const.F90:58
scale_const::const_eps1
real(rp), public const_eps1
small number
Definition: scale_const.F90:36
scale_const::const_ppm
real(rp), parameter, public const_ppm
parts par million
Definition: scale_const.F90:100
scale_atmos_phy_rd_mstrnx::atmos_phy_rd_mstrnx_setup
subroutine, public atmos_phy_rd_mstrnx_setup(KA, KS, KE, CZ, FZ)
Setup.
Definition: scale_atmos_phy_rd_mstrnx.F90:222
scale_atmos_phy_rd_profile::atmos_phy_rd_profile_setup_zgrid
subroutine, public atmos_phy_rd_profile_setup_zgrid(KA, KS, KE, KMAX, KADD, toa, CZ, FZ, zh, z)
Setup vertical grid for radiation.
Definition: scale_atmos_phy_rd_profile.F90:1099
scale_prof
module profiler
Definition: scale_prof.F90:11
scale_precision::dp
integer, parameter, public dp
Definition: scale_precision.F90:32
scale_time
module TIME
Definition: scale_time.F90:11
scale_const::const_pi
real(rp), parameter, public const_pi
pi
Definition: scale_const.F90:32
scale_atmos_hydrometeor::i_hc
integer, parameter, public i_hc
liquid water cloud
Definition: scale_atmos_hydrometeor.F90:97
scale_atmos_aerosol::ae_dens
real(rp), dimension(n_ae), parameter, public ae_dens
Definition: scale_atmos_aerosol.F90:53
scale_atmos_aerosol
module atmosphere / aerosol
Definition: scale_atmos_aerosol.F90:12
scale_atmos_phy_rd_common::i_sw
integer, parameter, public i_sw
Definition: scale_atmos_phy_rd_common.F90:38
scale_const::const_tem00
real(rp), parameter, public const_tem00
temperature reference (0C) [K]
Definition: scale_const.F90:99
scale_cpl_sfc_index
module coupler / surface-atmospehre
Definition: scale_cpl_sfc_index.F90:11
scale_time::time_nowdate
integer, dimension(6), public time_nowdate
current time [YYYY MM DD HH MM SS]
Definition: scale_time.F90:68
scale_cpl_sfc_index::i_r_vis
integer, parameter, public i_r_vis
Definition: scale_cpl_sfc_index.F90:31
scale_atmos_phy_rd_profile::atmos_phy_rd_profile_read
subroutine, public atmos_phy_rd_profile_read(kmax, ngas, ncfc, naero, real_lat, now_date, zh, z, rhodz, pres, presh, temp, temph, gas, cfc, aerosol_conc, aerosol_radi, cldfrac)
Read profile for radiation.
Definition: scale_atmos_phy_rd_profile.F90:546
scale_const::const_rdry
real(rp), public const_rdry
specific gas constant (dry air) [J/kg/K]
Definition: scale_const.F90:59
scale_atmos_phy_rd_common
module atmosphere / physics / radiation / common
Definition: scale_atmos_phy_rd_common.F90:12
scale_prof::prof_rapend
subroutine, public prof_rapend(rapname_base, level, disable_barrier)
Save raptime.
Definition: scale_prof.F90:246
scale_const::const_mvap
real(rp), public const_mvap
mass weight (water vapor) [g/mol]
Definition: scale_const.F90:67
scale_atmos_phy_rd_profile::atmos_phy_rd_profile_setup
subroutine, public atmos_phy_rd_profile_setup
Setup.
Definition: scale_atmos_phy_rd_profile.F90:135
scale_cpl_sfc_index::n_rad_rgn
integer, parameter, public n_rad_rgn
Definition: scale_cpl_sfc_index.F90:28
scale_atmos_phy_rd_mstrnx::atmos_phy_rd_mstrnx_flux
subroutine, public atmos_phy_rd_mstrnx_flux(KA, KS, KE, IA, IS, IE, JA, JS, JE, DENS, TEMP, PRES, QV, CZ, FZ, fact_ocean, fact_land, fact_urban, temp_sfc, albedo_sfc, solins, cosSZA, CLDFRAC, MP_Re, MP_Qe, AE_Re, AE_Qe, flux_rad, flux_rad_top, flux_rad_sfc_dn, dtau_s, dem_s)
Radiation main.
Definition: scale_atmos_phy_rd_mstrnx.F90:514
scale_atmos_hydrometeor::hyd_dens
real(rp), dimension(n_hyd), public hyd_dens
Definition: scale_atmos_hydrometeor.F90:111
scale_io::io_fid_conf
integer, public io_fid_conf
Config file ID.
Definition: scale_io.F90:57
scale_atmos_hydrometeor::n_hyd
integer, parameter, public n_hyd
Definition: scale_atmos_hydrometeor.F90:95
scale_atmos_phy_rd_mstrnx
module atmosphere / physics / radiation / mstrnX
Definition: scale_atmos_phy_rd_mstrnx.F90:15
scale_atmos_phy_rd_profile
module atmosphere / physics/ radiation / profile
Definition: scale_atmos_phy_rd_profile.F90:15
Generated by   doxygen 1.8.17