Line data Source code
1 : ! path: $Source: /storm/rc1/cvsroot/rc/rrtmg_sw/src/mcica_subcol_gen_sw.f90,v $
2 : ! author: $Author: mike $
3 : ! revision: $Revision: 1.4 $
4 : ! created: $Date: 2008/01/03 21:35:35 $
5 : !
6 :
7 : module mcica_subcol_gen_sw
8 :
9 : ! --------------------------------------------------------------------------
10 : ! | |
11 : ! | Copyright 2006-2007, Atmospheric & Environmental Research, Inc. (AER). |
12 : ! | This software may be used, copied, or redistributed as long as it is |
13 : ! | not sold and this copyright notice is reproduced on each copy made. |
14 : ! | This model is provided as is without any express or implied warranties. |
15 : ! | (http://www.rtweb.aer.com/) |
16 : ! | |
17 : ! --------------------------------------------------------------------------
18 :
19 : ! Purpose: Create McICA stochastic arrays for cloud physical or optical properties.
20 : ! Two options are possible:
21 : ! 1) Input cloud physical properties: cloud fraction, ice and liquid water
22 : ! paths, ice fraction, and particle sizes. Output will be stochastic
23 : ! arrays of these variables. (inflag = 1)
24 : ! 2) Input cloud optical properties directly: cloud optical depth, single
25 : ! scattering albedo and asymmetry parameter. Output will be stochastic
26 : ! arrays of these variables. (inflag = 0)
27 :
28 : ! --------- Modules ----------
29 :
30 : use shr_kind_mod, only: r8 => shr_kind_r8
31 : use cam_abortutils, only: endrun
32 :
33 : use parrrsw, only: ngptsw
34 : use rrsw_wvn, only: ngb
35 :
36 : implicit none
37 : private
38 :
39 : public :: mcica_subcol_sw
40 :
41 : !=========================================================================================
42 : contains
43 : !=========================================================================================
44 :
45 115200 : subroutine mcica_subcol_sw(lchnk, ncol, nlay, icld, permuteseed, play, &
46 38400 : cldfrac, ciwp, clwp, rei, rel, tauc, ssac, asmc, fsfc, &
47 76800 : cldfmcl, ciwpmcl, clwpmcl, reicmcl, relqmcl, &
48 38400 : taucmcl, ssacmcl, asmcmcl, fsfcmcl)
49 :
50 : ! ----- Input -----
51 : ! Control
52 : integer, intent(in) :: lchnk ! chunk identifier
53 : integer, intent(in) :: ncol ! number of columns
54 : integer, intent(in) :: nlay ! number of model layers
55 : integer, intent(in) :: icld ! clear/cloud, cloud overlap flag
56 : integer, intent(in) :: permuteseed ! if the cloud generator is called multiple times,
57 : ! permute the seed between each call;
58 : ! between calls for LW and SW, recommended
59 : ! permuteseed differs by 'ngpt'
60 :
61 : ! Atmosphere
62 : real(kind=r8), intent(in) :: play(:,:) ! layer pressures (mb)
63 : ! Dimensions: (ncol,nlay)
64 :
65 : ! Atmosphere/clouds - cldprop
66 : real(kind=r8), intent(in) :: cldfrac(:,:) ! layer cloud fraction
67 : ! Dimensions: (ncol,nlay)
68 : real(kind=r8), intent(in) :: tauc(:,:,:) ! cloud optical depth
69 : ! Dimensions: (nbndsw,ncol,nlay)
70 : real(kind=r8), intent(in) :: ssac(:,:,:) ! cloud single scattering albedo (non-delta scaled)
71 : ! Dimensions: (nbndsw,ncol,nlay)
72 : real(kind=r8), intent(in) :: asmc(:,:,:) ! cloud asymmetry parameter (non-delta scaled)
73 : ! Dimensions: (nbndsw,ncol,nlay)
74 : real(kind=r8), intent(in) :: fsfc(:,:,:) ! cloud forward scattering fraction (non-delta scaled)
75 : ! Dimensions: (nbndsw,ncol,nlay)
76 : real(kind=r8), intent(in) :: ciwp(:,:) ! cloud ice water path
77 : ! Dimensions: (ncol,nlay)
78 : real(kind=r8), intent(in) :: clwp(:,:) ! cloud liquid water path
79 : ! Dimensions: (ncol,nlay)
80 : real(kind=r8), intent(in) :: rei(:,:) ! cloud ice particle size
81 : ! Dimensions: (ncol,nlay)
82 : real(kind=r8), intent(in) :: rel(:,:) ! cloud liquid particle size
83 : ! Dimensions: (ncol,nlay)
84 :
85 : ! ----- Output -----
86 : ! Atmosphere/clouds - cldprmc [mcica]
87 : real(kind=r8), intent(out) :: cldfmcl(:,:,:) ! cloud fraction [mcica]
88 : ! Dimensions: (ngptsw,ncol,nlay)
89 : real(kind=r8), intent(out) :: ciwpmcl(:,:,:) ! cloud ice water path [mcica]
90 : ! Dimensions: (ngptsw,ncol,nlay)
91 : real(kind=r8), intent(out) :: clwpmcl(:,:,:) ! cloud liquid water path [mcica]
92 : ! Dimensions: (ngptsw,ncol,nlay)
93 : real(kind=r8), intent(out) :: relqmcl(:,:) ! liquid particle size (microns)
94 : ! Dimensions: (ncol,nlay)
95 : real(kind=r8), intent(out) :: reicmcl(:,:) ! ice partcle size (microns)
96 : ! Dimensions: (ncol,nlay)
97 : real(kind=r8), intent(out) :: taucmcl(:,:,:) ! cloud optical depth [mcica]
98 : ! Dimensions: (ngptsw,ncol,nlay)
99 : real(kind=r8), intent(out) :: ssacmcl(:,:,:) ! cloud single scattering albedo [mcica]
100 : ! Dimensions: (ngptsw,ncol,nlay)
101 : real(kind=r8), intent(out) :: asmcmcl(:,:,:) ! cloud asymmetry parameter [mcica]
102 : ! Dimensions: (ngptsw,ncol,nlay)
103 : real(kind=r8), intent(out) :: fsfcmcl(:,:,:) ! cloud forward scattering fraction [mcica]
104 : ! Dimensions: (ngptsw,ncol,nlay)
105 :
106 : ! ----- Local -----
107 :
108 : ! Stochastic cloud generator variables [mcica]
109 : integer, parameter :: nsubcsw = ngptsw ! number of sub-columns (g-point intervals)
110 : integer :: km, im, nm ! loop indices
111 :
112 0 : real(kind=r8) :: pmid(ncol,nlay) ! layer pressures (Pa)
113 : !----------------------------------------------------------------------------
114 :
115 : ! Return if clear sky; or stop if icld out of range
116 38400 : if (icld.eq.0) return
117 38400 : if (icld.lt.0.or.icld.gt.3) then
118 0 : call endrun('MCICA_SUBCOL: INVALID ICLD')
119 : endif
120 :
121 : ! NOTE: For GCM mode, permuteseed must be offset between LW and SW by at
122 : ! least number of subcolumns
123 :
124 : ! Pass particle sizes to new arrays, no subcolumns for these properties yet
125 : ! Convert pressures from mb to Pa
126 :
127 17041920 : reicmcl(:ncol,:nlay) = rei(:ncol,:nlay)
128 17041920 : relqmcl(:ncol,:nlay) = rel(:ncol,:nlay)
129 17041920 : pmid(:ncol,:nlay) = play(:ncol,:nlay)*1.e2_r8
130 :
131 : ! Generate the stochastic subcolumns of cloud optical properties for the shortwave;
132 : call generate_stochastic_clouds_sw( &
133 : ncol, nlay, nsubcsw, icld, pmid, &
134 : cldfrac, clwp, ciwp, tauc, ssac, &
135 : asmc, fsfc, cldfmcl, clwpmcl, ciwpmcl, &
136 38400 : taucmcl, ssacmcl, asmcmcl, fsfcmcl, permuteseed)
137 :
138 38400 : end subroutine mcica_subcol_sw
139 :
140 : !=========================================================================================
141 :
142 38400 : subroutine generate_stochastic_clouds_sw( &
143 76800 : ncol, nlay, nsubcol, icld, pmid, &
144 38400 : cld, clwp, ciwp, tauc, ssac, &
145 38400 : asmc, fsfc, cld_stoch, clwp_stoch, ciwp_stoch, &
146 38400 : tauc_stoch, ssac_stoch, asmc_stoch, fsfc_stoch, changeSeed)
147 :
148 : !----------------------------------------------------------------------------------------------------------------
149 : ! ---------------------
150 : ! Contact: Cecile Hannay (hannay@ucar.edu)
151 : !
152 : ! Original code: Based on Raisanen et al., QJRMS, 2004.
153 : !
154 : ! Modifications: Generalized for use with RRTMG and added Mersenne Twister as the default
155 : ! random number generator, which can be changed to the optional kissvec random number generator
156 : ! with flag 'irnd' below . Some extra functionality has been commented or removed.
157 : ! Michael J. Iacono, AER, Inc., February 2007
158 : !
159 : ! Given a profile of cloud fraction, cloud water and cloud ice, we produce a set of subcolumns.
160 : ! Each layer within each subcolumn is homogeneous, with cloud fraction equal to zero or one
161 : ! and uniform cloud liquid and cloud ice concentration.
162 : ! The ensemble as a whole reproduces the probability function of cloud liquid and ice within each layer
163 : ! and obeys an overlap assumption in the vertical.
164 : !
165 : ! Overlap assumption:
166 : ! The cloud are consistent with 4 overlap assumptions: random, maximum, maximum-random and exponential.
167 : ! The default option is maximum-random (option 3)
168 : ! The options are: 1=random overlap, 2=max/random, 3=maximum overlap, 4=exponential overlap
169 : ! This is set with the variable "overlap"
170 : !mji - Exponential overlap option (overlap=4) has been deactivated in this version
171 : ! The exponential overlap uses also a length scale, Zo. (real, parameter :: Zo = 2500. )
172 : !
173 : ! Seed:
174 : ! If the stochastic cloud generator is called several times during the same timestep,
175 : ! one should change the seed between the call to insure that the subcolumns are different.
176 : ! This is done by changing the argument 'changeSeed'
177 : ! For example, if one wants to create a set of columns for the shortwave and another set for the longwave ,
178 : ! use 'changeSeed = 1' for the first call and'changeSeed = 2' for the second call
179 : !
180 : ! PDF assumption:
181 : ! We can use arbitrary complicated PDFS.
182 : ! In the present version, we produce homogeneuous clouds (the simplest case).
183 : ! Future developments include using the PDF scheme of Ben Johnson.
184 : !
185 : ! History file:
186 : ! Option to add diagnostics variables in the history file. (using FINCL in the namelist)
187 : ! nsubcol = number of subcolumns
188 : ! overlap = overlap type (1-3)
189 : ! Zo = length scale
190 : ! CLOUD_S = mean of the subcolumn cloud fraction ('_S" means Stochastic)
191 : ! CLDLIQ_S = mean of the subcolumn cloud water
192 : ! CLDICE_S = mean of the subcolumn cloud ice
193 : !
194 : ! Note:
195 : ! Here: we force that the cloud condensate to be consistent with the cloud fraction
196 : ! i.e we only have cloud condensate when the cell is cloudy.
197 : ! In CAM: The cloud condensate and the cloud fraction are obtained from 2 different equations
198 : ! and the 2 quantities can be inconsistent (i.e. CAM can produce cloud fraction
199 : ! without cloud condensate or the opposite).
200 : !---------------------------------------------------------------------------------------------------------------
201 :
202 : use shr_RandNum_mod, only: ShrIntrinsicRandGen, ShrKissRandGen, &
203 : ShrF95MtRandGen, ShrDsfmtRandGen
204 :
205 38400 : type(ShrDsfmtRandGen) :: dsfmt_gen
206 38400 : type(ShrKissRandGen) :: kiss_gen
207 :
208 : ! -- Arguments
209 :
210 : integer, intent(in) :: ncol ! number of layers
211 : integer, intent(in) :: nlay ! number of layers
212 : integer, intent(in) :: icld ! clear/cloud, cloud overlap flag
213 : integer, intent(in) :: nsubcol ! number of sub-columns (g-point intervals)
214 : integer, optional, intent(in) :: changeSeed ! allows permuting seed
215 :
216 : real(kind=r8), intent(in) :: pmid(:,:) ! layer pressure (Pa)
217 : ! Dimensions: (ncol,nlay)
218 : real(kind=r8), intent(in) :: cld(:,:) ! cloud fraction
219 : ! Dimensions: (ncol,nlay)
220 : real(kind=r8), intent(in) :: clwp(:,:) ! cloud liquid water path (g/m2)
221 : ! Dimensions: (ncol,nlay)
222 : real(kind=r8), intent(in) :: ciwp(:,:) ! cloud ice water path (g/m2)
223 : ! Dimensions: (ncol,nlay)
224 : real(kind=r8), intent(in) :: tauc(:,:,:) ! cloud optical depth (non-delta scaled)
225 : ! Dimensions: (nbndsw,ncol,nlay)
226 : real(kind=r8), intent(in) :: ssac(:,:,:) ! cloud single scattering albedo (non-delta scaled)
227 : ! Dimensions: (nbndsw,ncol,nlay)
228 : real(kind=r8), intent(in) :: asmc(:,:,:) ! cloud asymmetry parameter (non-delta scaled)
229 : ! Dimensions: (nbndsw,ncol,nlay)
230 : real(kind=r8), intent(in) :: fsfc(:,:,:) ! cloud forward scattering fraction (non-delta scaled)
231 : ! Dimensions: (nbndsw,ncol,nlay)
232 :
233 : real(kind=r8), intent(out) :: cld_stoch(:,:,:) ! subcolumn cloud fraction
234 : ! Dimensions: (ngptsw,ncol,nlay)
235 : real(kind=r8), intent(out) :: clwp_stoch(:,:,:) ! subcolumn cloud liquid water path
236 : ! Dimensions: (ngptsw,ncol,nlay)
237 : real(kind=r8), intent(out) :: ciwp_stoch(:,:,:) ! subcolumn cloud ice water path
238 : ! Dimensions: (ngptsw,ncol,nlay)
239 : real(kind=r8), intent(out) :: tauc_stoch(:,:,:) ! subcolumn cloud optical depth
240 : ! Dimensions: (ngptsw,ncol,nlay)
241 : real(kind=r8), intent(out) :: ssac_stoch(:,:,:) ! subcolumn cloud single scattering albedo
242 : ! Dimensions: (ngptsw,ncol,nlay)
243 : real(kind=r8), intent(out) :: asmc_stoch(:,:,:) ! subcolumn cloud asymmetry parameter
244 : ! Dimensions: (ngptsw,ncol,nlay)
245 : real(kind=r8), intent(out) :: fsfc_stoch(:,:,:) ! subcolumn cloud forward scattering fraction
246 : ! Dimensions: (ngptsw,ncol,nlay)
247 :
248 : ! -- Local variables
249 76800 : real(kind=r8) :: cldf(ncol,nlay) ! cloud fraction
250 : ! Dimensions: (ncol,nlay)
251 :
252 : ! Constants (min value for cloud fraction and cloud water and ice)
253 : real(kind=r8), parameter :: cldmin = 1.0e-80_r8 ! min cloud fraction
254 :
255 : ! Variables related to random number and seed
256 : integer :: irnd ! flag for random number generator
257 : ! 0 = kissvec
258 : ! 1 = Mersenne Twister
259 :
260 76800 : real(kind=r8), dimension(nsubcol, ncol, nlay) :: CDF, CDF2 ! random numbers
261 76800 : integer, dimension(ncol,4) :: kiss_seed
262 76800 : real(kind=r8), dimension(ncol,1) :: rand_num_1d ! random number
263 76800 : real(kind=r8), dimension(ncol,nlay) :: rand_num ! random number
264 38400 : integer, dimension(ncol) :: iseed ! seed to create random number (Mersenne Teister)
265 :
266 : ! Flag to identify cloud fraction in subcolumns
267 76800 : logical, dimension(nsubcol, ncol, nlay) :: isCloudy ! flag that says whether a gridbox is cloudy
268 :
269 : ! Indices
270 : integer :: ilev, isubcol, i, n, ngbm ! indices
271 : !----------------------------------------------------------------------------
272 :
273 : ! Set randum number generator to use (0 = kissvec; 1 = mersennetwister)
274 38400 : irnd = 0
275 :
276 : ! ensure that cloud fractions are in bounds
277 17041920 : cldf(:,:) = cld(:ncol,:nlay)
278 17041920 : where (cldf(:,:) < cldmin)
279 : cldf(:,:) = 0._r8
280 : end where
281 :
282 : ! ----- Create seed --------
283 :
284 : ! Advance randum number generator by changeseed values
285 : if (irnd == 0) then
286 :
287 : ! For kissvec, create a seed that depends on the state of the columns. Maybe not the best way, but it works.
288 : ! Must use pmid from bottom four layers.
289 314880 : do i = 1, ncol
290 276480 : if (pmid(i,nlay) < pmid(i,nlay-1)) then
291 0 : call endrun('MCICA_SUBCOL: KISSVEC SEED GENERATOR REQUIRES PMID FROM BOTTOM FOUR LAYERS.')
292 : end if
293 276480 : kiss_seed(i,1) = (pmid(i,nlay) - int(pmid(i,nlay))) * 1000000000
294 276480 : kiss_seed(i,2) = (pmid(i,nlay-1) - int(pmid(i,nlay-1))) * 1000000000
295 276480 : kiss_seed(i,3) = (pmid(i,nlay-2) - int(pmid(i,nlay-2))) * 1000000000
296 314880 : kiss_seed(i,4) = (pmid(i,nlay-3) - int(pmid(i,nlay-3))) * 1000000000
297 : end do
298 :
299 38400 : kiss_gen = ShrKissRandGen(kiss_seed)
300 :
301 76800 : do i = 1, changeSeed
302 76800 : call kiss_gen%random(rand_num_1d)
303 : end do
304 :
305 : elseif (irnd == 1) then
306 :
307 : do i = 1, ncol
308 : if (pmid(i,nlay) < pmid(i,nlay-1)) then
309 : call endrun('MCICA_SUBCOL: MT SEED GENERATOR REQUIRES PMID FROM BOTTOM FOUR LAYERS.')
310 : end if
311 : kiss_seed(i,1) = (pmid(i,nlay) - int(pmid(i,nlay))) * 1000000000
312 : kiss_seed(i,2) = (pmid(i,nlay-1) - int(pmid(i,nlay-1))) * 1000000000
313 : kiss_seed(i,3) = (pmid(i,nlay-2) - int(pmid(i,nlay-2))) * 1000000000
314 : kiss_seed(i,4) = (pmid(i,nlay-3) - int(pmid(i,nlay-3))) * 1000000000
315 : end do
316 :
317 : iseed = kiss_seed(:,1) + kiss_seed(:,2) + kiss_seed(:,3) + kiss_seed(:,4)
318 : dsfmt_gen =ShrDsfmtRandGen(iseed,1)
319 :
320 : end if
321 :
322 : ! ------ Apply overlap assumption --------
323 :
324 : ! generate the random numbers
325 :
326 38400 : select case (icld)
327 :
328 : case(1)
329 : ! Random overlap
330 : ! i) pick a random value at every level
331 : if (irnd == 0) then
332 0 : do isubcol = 1,nsubcol
333 0 : call kiss_gen%random(rand_num)
334 0 : CDF(isubcol,:,:) = rand_num(:,:)
335 : end do
336 : else if (irnd == 1) then
337 : do isubcol = 1, nsubcol
338 : call dsfmt_gen%random(rand_num)
339 : CDF(isubcol,:,:) = rand_num(:,:)
340 : end do
341 : end if
342 :
343 : case(2)
344 : ! Maximum-Random overlap
345 : ! i) pick a random number for top layer.
346 : ! ii) walk down the column:
347 : ! - if the layer above is cloudy, we use the same random number than in the layer above
348 : ! - if the layer above is clear, we use a new random number
349 :
350 : if (irnd == 0) then
351 4339200 : do isubcol = 1, nsubcol
352 4300800 : call kiss_gen%random(rand_num)
353 1908733440 : CDF(isubcol,:,:) = rand_num(:,:)
354 : end do
355 : else if (irnd == 1) then
356 : do isubcol = 1, nsubcol
357 : call dsfmt_gen%random(rand_num)
358 : CDF(isubcol,:,:) = rand_num(:,:)
359 : end do
360 : end if
361 :
362 2073600 : do ilev = 2, nlay
363 16727040 : do i = 1, ncol
364 1657873920 : do isubcol = 1, nsubcol
365 1655838720 : if (CDF(isubcol, i, ilev-1) > 1._r8 - cldf(i,ilev-1) ) then
366 140094220 : CDF(isubcol,i,ilev) = CDF(isubcol,i,ilev-1)
367 : else
368 1501091060 : CDF(isubcol,i,ilev) = CDF(isubcol,i,ilev) * (1._r8 - cldf(i,ilev-1))
369 : end if
370 : end do
371 : end do
372 : end do
373 :
374 : case(3)
375 : ! Maximum overlap
376 : ! i) pick same random numebr at every level
377 :
378 38400 : if (irnd == 0) then
379 0 : do isubcol = 1, nsubcol
380 0 : call kiss_gen%random(rand_num_1d)
381 0 : do ilev = 1,nlay
382 0 : CDF(isubcol,:,ilev) = rand_num_1d(:,1)
383 : end do
384 : end do
385 : else if (irnd == 1) then
386 : do isubcol = 1, nsubcol
387 : call dsfmt_gen%random(rand_num_1d)
388 : do ilev = 1, nlay
389 : CDF(isubcol,:,ilev) = rand_num_1d(:,1)
390 : end do
391 : end do
392 : end if
393 :
394 : end select
395 :
396 : ! -- generate subcolumns for homogeneous clouds -----
397 2112000 : do ilev = 1, nlay
398 1689192960 : isCloudy(:,:,ilev) = (CDF(:,:,ilev) >= 1._r8 - spread(cldf(:,ilev), dim=1, nCopies=nsubcol) )
399 : end do
400 :
401 : ! where the subcolumn is cloudy, the subcolumn cloud fraction is 1;
402 : ! where the subcolumn is not cloudy, the subcolumn cloud fraction is 0
403 :
404 2112000 : do ilev = 1, nlay
405 17041920 : do i = 1, ncol
406 1689154560 : do isubcol = 1, nsubcol
407 1687080960 : if (iscloudy(isubcol,i,ilev) ) then
408 145606279 : cld_stoch(isubcol,i,ilev) = 1._r8
409 : else
410 1526544761 : cld_stoch(isubcol,i,ilev) = 0._r8
411 : end if
412 : end do
413 : end do
414 : end do
415 :
416 : ! where there is a cloud, set the subcolumn cloud properties;
417 : ! Incoming clwp, ciwp and tauc should be in-cloud quantites and not grid-averaged quantities
418 :
419 2112000 : do ilev = 1, nlay
420 17041920 : do i = 1, ncol
421 1689154560 : do isubcol = 1, nsubcol
422 1687080960 : if ( iscloudy(isubcol,i,ilev) .and. (cldf(i,ilev) > 0._r8) ) then
423 145606279 : clwp_stoch(isubcol,i,ilev) = clwp(i,ilev)
424 145606279 : ciwp_stoch(isubcol,i,ilev) = ciwp(i,ilev)
425 : else
426 1526544761 : clwp_stoch(isubcol,i,ilev) = 0._r8
427 1526544761 : ciwp_stoch(isubcol,i,ilev) = 0._r8
428 : end if
429 : end do
430 : end do
431 : end do
432 :
433 38400 : ngbm = ngb(1) - 1
434 2112000 : do ilev = 1,nlay
435 17041920 : do i = 1, ncol
436 1689154560 : do isubcol = 1, nsubcol
437 1687080960 : if ( iscloudy(isubcol,i,ilev) .and. (cldf(i,ilev) > 0._r8) ) then
438 145606279 : n = ngb(isubcol) - ngbm
439 145606279 : tauc_stoch(isubcol,i,ilev) = tauc(n,i,ilev)
440 145606279 : ssac_stoch(isubcol,i,ilev) = ssac(n,i,ilev)
441 145606279 : asmc_stoch(isubcol,i,ilev) = asmc(n,i,ilev)
442 145606279 : fsfc_stoch(isubcol,i,ilev) = fsfc(n,i,ilev)
443 : else
444 1526544761 : tauc_stoch(isubcol,i,ilev) = 0._r8
445 1526544761 : ssac_stoch(isubcol,i,ilev) = 1._r8
446 1526544761 : asmc_stoch(isubcol,i,ilev) = 0._r8
447 1526544761 : fsfc_stoch(isubcol,i,ilev) = 0._r8
448 : end if
449 : end do
450 : end do
451 : end do
452 :
453 : if (irnd == 0) then
454 38400 : call kiss_gen%finalize()
455 : else if (irnd == 1) then
456 : call dsfmt_gen%finalize()
457 : end if
458 :
459 38400 : end subroutine generate_stochastic_clouds_sw
460 :
461 : end module mcica_subcol_gen_sw
|