Line data Source code
1 : module mcica_subcol_gen
2 :
3 : !----------------------------------------------------------------------------------------
4 : !
5 : ! Purpose: Create McICA stochastic arrays for cloud optical properties.
6 : ! Input cloud optical properties directly: cloud optical depth, single
7 : ! scattering albedo and asymmetry parameter. Output will be stochastic
8 : ! arrays of these variables. (longwave scattering is not yet available)
9 : !
10 : ! Original code: From RRTMG, with the following copyright notice,
11 : ! based on Raisanen et al., QJRMS, 2004:
12 : ! --------------------------------------------------------------------------
13 : ! | |
14 : ! | Copyright 2006-2007, Atmospheric & Environmental Research, Inc. (AER). |
15 : ! | This software may be used, copied, or redistributed as long as it is |
16 : ! | not sold and this copyright notice is reproduced on each copy made. |
17 : ! | This model is provided as is without any express or implied warranties. |
18 : ! | (http://www.rtweb.aer.com/) |
19 : ! | |
20 : ! --------------------------------------------------------------------------
21 : ! This code is a refactored version of code originally in the files
22 : ! mcica_subcol_gen_lw.F90 and mcica_subcol_gen_sw.F90
23 : !
24 : ! Uses the KISS random number generator.
25 : !
26 : ! Overlap assumption: maximum-random.
27 : !
28 : !----------------------------------------------------------------------------------------
29 :
30 : use shr_kind_mod, only: r8 => shr_kind_r8
31 : use ppgrid, only: pcols, pver
32 : use shr_RandNum_mod, only: ShrKissRandGen
33 : use mo_gas_optics_rrtmgp, only: ty_gas_optics_rrtmgp
34 :
35 : implicit none
36 : private
37 : save
38 :
39 : public :: mcica_subcol_lw, mcica_subcol_sw
40 :
41 : !========================================================================================
42 : contains
43 : !========================================================================================
44 :
45 746136 : subroutine mcica_subcol_lw( &
46 : kdist, nbnd, ngpt, ncol, nver, &
47 746136 : changeseed, pmid, cldfrac, tauc, taucmcl )
48 :
49 : ! Arrays use CAM vertical index convention: index increases from top to bottom.
50 : ! This index ordering is assumed in the maximum-random overlap algorithm which starts
51 : ! at the top of a column and marches down, with each layer depending on the state
52 : ! of the layer above it.
53 : !
54 : ! For GCM mode, changeseed must be offset between LW and SW by at least the
55 : ! number of subcolumns
56 :
57 : ! arguments
58 : class(ty_gas_optics_rrtmgp), intent(in) :: kdist ! spectral information
59 : integer, intent(in) :: nbnd ! number of spectral bands
60 : integer, intent(in) :: ngpt ! number of subcolumns (g-point intervals)
61 : integer, intent(in) :: ncol ! number of columns
62 : integer, intent(in) :: nver ! number of layers
63 : integer, intent(in) :: changeseed ! if the subcolumn generator is called multiple times,
64 : ! permute the seed between each call.
65 : real(r8), intent(in) :: pmid(pcols,pver) ! layer pressures (Pa)
66 : real(r8), intent(in) :: cldfrac(ncol,nver) ! layer cloud fraction
67 : real(r8), intent(in) :: tauc(nbnd,ncol,nver) ! cloud optical depth
68 : real(r8), intent(out) :: taucmcl(ngpt,ncol,nver) ! subcolumn cloud optical depth [mcica]
69 :
70 : ! Local variables
71 :
72 : integer :: i, isubcol, k, n
73 :
74 : real(r8), parameter :: cldmin = 1.0e-80_r8 ! min cloud fraction
75 1492272 : real(r8) :: cldf(ncol,pver) ! cloud fraction clipped to cldmin
76 :
77 746136 : type(ShrKissRandGen) :: kiss_gen ! KISS RNG object
78 1492272 : integer :: kiss_seed(ncol,4)
79 1492272 : real(r8) :: rand_num_1d(ncol,1) ! random number (kissvec)
80 1492272 : real(r8) :: rand_num(ncol,nver) ! random number (kissvec)
81 :
82 1492272 : real(r8) :: cdf(ngpt,ncol,nver) ! random numbers
83 1492272 : logical :: iscloudy(ngpt,ncol,nver) ! flag that says whether a gridbox is cloudy
84 : !------------------------------------------------------------------------------------------
85 :
86 : ! clip cloud fraction
87 1159408584 : cldf(:,:) = cldfrac(:ncol,:)
88 1159408584 : where (cldf(:,:) < cldmin)
89 : cldf(:,:) = 0._r8
90 : end where
91 :
92 : ! Create a seed that depends on the state of the columns.
93 : ! Use pmid from bottom four layers.
94 12458736 : do i = 1, ncol
95 11712600 : kiss_seed(i,1) = (pmid(i,pver) - int(pmid(i,pver))) * 1000000000
96 11712600 : kiss_seed(i,2) = (pmid(i,pver-1) - int(pmid(i,pver-1))) * 1000000000
97 11712600 : kiss_seed(i,3) = (pmid(i,pver-2) - int(pmid(i,pver-2))) * 1000000000
98 12458736 : kiss_seed(i,4) = (pmid(i,pver-3) - int(pmid(i,pver-3))) * 1000000000
99 : end do
100 :
101 : ! create the RNG object
102 746136 : kiss_gen = ShrKissRandGen(kiss_seed)
103 :
104 : ! Advance randum number generator by changeseed values
105 96251544 : do i = 1, changeSeed
106 96251544 : call kiss_gen%random(rand_num_1d)
107 : end do
108 :
109 : ! Generate random numbers in each subcolumn at every level
110 96251544 : do isubcol = 1,ngpt
111 95505408 : call kiss_gen%random(rand_num)
112 >14840*10^7 : cdf(isubcol,:,:) = rand_num(:,:)
113 : enddo
114 :
115 : ! Maximum-Random overlap
116 : ! i) pick a random number for top layer.
117 : ! ii) walk down the column:
118 : ! - if the layer above is cloudy, use the same random number as in the layer above
119 : ! - if the layer above is clear, use a new random number
120 :
121 69390648 : do k = 2, nver
122 1146949848 : do i = 1, ncol
123 >13907*10^7 : do isubcol = 1, ngpt
124 >13900*10^7 : if (cdf(isubcol,i,k-1) > 1._r8 - cldf(i,k-1) ) then
125 11091485592 : cdf(isubcol,i,k) = cdf(isubcol,i,k-1)
126 : else
127 >12683*10^7 : cdf(isubcol,i,k) = cdf(isubcol,i,k) * (1._r8 - cldf(i,k-1))
128 : end if
129 : end do
130 : end do
131 : end do
132 :
133 70136784 : do k = 1, nver
134 >14058*10^7 : iscloudy(:,:,k) = (cdf(:,:,k) >= 1._r8 - spread(cldf(:,k), dim=1, nCopies=ngpt) )
135 : end do
136 :
137 : ! -- generate subcolumns for homogeneous clouds -----
138 : ! where there is a cloud, set the subcolumn cloud properties;
139 : ! incoming tauc should be in-cloud quantites and not grid-averaged quantities
140 70136784 : do k = 1,nver
141 1159408584 : do i = 1,ncol
142 >14058*10^7 : do isubcol = 1,ngpt
143 >14051*10^7 : if (iscloudy(isubcol,i,k) .and. (cldf(i,k) > 0._r8) ) then
144 11261941296 : n = kdist%convert_gpt2band(isubcol)
145 11261941296 : taucmcl(isubcol,i,k) = tauc(n,i,k)
146 : else
147 >12816*10^7 : taucmcl(isubcol,i,k) = 0._r8
148 : end if
149 : end do
150 : end do
151 : end do
152 :
153 746136 : call kiss_gen%finalize()
154 :
155 746136 : end subroutine mcica_subcol_lw
156 :
157 : !========================================================================================
158 :
159 389263 : subroutine mcica_subcol_sw( &
160 : kdist, nbnd, ngpt, ncol, nlay, nver, changeseed, &
161 389263 : pmid, cldfrac, tauc, ssac, asmc, &
162 389263 : taucmcl, ssacmcl, asmcmcl)
163 :
164 : ! Arrays use CAM vertical index convention: index increases from top to bottom.
165 : ! This index ordering is assumed in the maximum-random overlap algorithm which starts
166 : ! at the top of a column and marches down, with each layer depending on the state
167 : ! of the layer above it.
168 : !
169 : ! For GCM mode, changeseed must be offset between LW and SW by at least the
170 : ! number of subcolumns
171 :
172 : ! arguments
173 : class(ty_gas_optics_rrtmgp), intent(in) :: kdist ! spectral information
174 : integer, intent(in) :: nbnd ! number of spectral bands
175 : integer, intent(in) :: ngpt ! number of subcolumns (g-point intervals)
176 : integer, intent(in) :: ncol ! number of columns
177 : integer, intent(in) :: nlay ! number of vertical layers in radiation calc;
178 : ! may include an "extra layer"
179 : integer, intent(in) :: nver ! number of CAM's vertical layers in rad calc
180 : integer, intent(in) :: changeseed ! if the subcolumn generator is called multiple times,
181 : ! permute the seed between each call.
182 : real(r8), intent(in) :: pmid(ncol,nlay) ! layer midpoint pressures (Pa)
183 : real(r8), intent(in) :: cldfrac(ncol,nver) ! layer cloud fraction
184 : real(r8), intent(in) :: tauc(nbnd,ncol,nver) ! cloud optical depth
185 : real(r8), intent(in) :: ssac(nbnd,ncol,nver) ! cloud single scattering albedo (non-delta scaled)
186 : real(r8), intent(in) :: asmc(nbnd,ncol,nver) ! cloud asymmetry parameter (non-delta scaled)
187 :
188 :
189 : real(r8), intent(out) :: taucmcl(ngpt,ncol,nver) ! subcolumn cloud optical depth [mcica]
190 : real(r8), intent(out) :: ssacmcl(ngpt,ncol,nver) ! subcolumn cloud single scattering albedo [mcica]
191 : real(r8), intent(out) :: asmcmcl(ngpt,ncol,nver) ! subcolumn cloud asymmetry parameter [mcica]
192 :
193 : ! Local vars
194 :
195 : integer :: i, isubcol, k, n
196 :
197 : real(r8), parameter :: cldmin = 1.0e-80_r8 ! min cloud fraction
198 778526 : real(r8) :: cldf(ncol,nver) ! cloud fraction clipped to cldmin
199 :
200 389263 : type(ShrKissRandGen) :: kiss_gen ! KISS RNG object
201 778526 : integer :: kiss_seed(ncol,4)
202 778526 : real(r8) :: rand_num_1d(ncol,1) ! random number (kissvec)
203 778526 : real(r8) :: rand_num(ncol,nver) ! random number (kissvec)
204 :
205 778526 : real(r8) :: cdf(ngpt,ncol,nver) ! random numbers
206 389263 : logical :: iscloudy(ngpt,ncol,nver) ! flag that says whether a gridbox is cloudy
207 : !------------------------------------------------------------------------------------------
208 :
209 : ! clip cloud fraction
210 581226622 : cldf(:,:) = cldfrac(:ncol,:)
211 581226622 : where (cldf(:,:) < cldmin)
212 : cldf(:,:) = 0._r8
213 : end where
214 :
215 : ! Create a seed that depends on the state of the columns.
216 : ! Use pmid from bottom four layers.
217 6245563 : do i = 1, ncol
218 5856300 : kiss_seed(i,1) = (pmid(i,nlay) - int(pmid(i,nlay))) * 1000000000
219 5856300 : kiss_seed(i,2) = (pmid(i,nlay-1) - int(pmid(i,nlay-1))) * 1000000000
220 5856300 : kiss_seed(i,3) = (pmid(i,nlay-2) - int(pmid(i,nlay-2))) * 1000000000
221 6245563 : kiss_seed(i,4) = (pmid(i,nlay-3) - int(pmid(i,nlay-3))) * 1000000000
222 : end do
223 :
224 : ! create the RNG object
225 389263 : kiss_gen = ShrKissRandGen(kiss_seed)
226 :
227 : ! Advance randum number generator by changeseed values
228 778526 : do i = 1, changeSeed
229 778526 : call kiss_gen%random(rand_num_1d)
230 : end do
231 :
232 : ! Generate random numbers in each subcolumn at every level
233 43986719 : do isubcol = 1,ngpt
234 43597456 : call kiss_gen%random(rand_num)
235 65097770927 : cdf(isubcol,:,:) = rand_num(:,:)
236 : enddo
237 :
238 : ! Maximum-Random overlap
239 : ! i) pick a random number for top layer.
240 : ! ii) walk down the column:
241 : ! - if the layer above is cloudy, use the same random number as in the layer above
242 : ! - if the layer above is clear, use a new random number
243 :
244 36201459 : do k = 2, nver
245 574981059 : do i = 1, ncol
246 60917906996 : do isubcol = 1, ngpt
247 60882094800 : if (cdf(isubcol,i,k-1) > 1._r8 - cldf(i,k-1) ) then
248 4754772523 : cdf(isubcol,i,k) = cdf(isubcol,i,k-1)
249 : else
250 55588542677 : cdf(isubcol,i,k) = cdf(isubcol,i,k) * (1._r8 - cldf(i,k-1))
251 : end if
252 : end do
253 : end do
254 : end do
255 :
256 36590722 : do k = 1, nver
257 61580447422 : iscloudy(:,:,k) = (cdf(:,:,k) >= 1._r8 - spread(cldf(:,k), dim=1, nCopies=ngpt) )
258 : end do
259 :
260 : ! -- generate subcolumns for homogeneous clouds -----
261 : ! where there is a cloud, set the subcolumn cloud properties;
262 : ! incoming tauc should be in-cloud quantites and not grid-averaged quantities
263 36590722 : do k = 1,nver
264 581226622 : do i = 1,ncol
265 61580058159 : do isubcol = 1,ngpt
266 61543856700 : if (iscloudy(isubcol,i,k) .and. (cldf(i,k) > 0._r8) ) then
267 4805438153 : n = kdist%convert_gpt2band(isubcol)
268 4805438153 : taucmcl(isubcol,i,k) = tauc(n,i,k)
269 4805438153 : ssacmcl(isubcol,i,k) = ssac(n,i,k)
270 4805438153 : asmcmcl(isubcol,i,k) = asmc(n,i,k)
271 : else
272 56193782647 : taucmcl(isubcol,i,k) = 0._r8
273 56193782647 : ssacmcl(isubcol,i,k) = 1._r8
274 56193782647 : asmcmcl(isubcol,i,k) = 0._r8
275 : end if
276 : end do
277 : end do
278 : end do
279 :
280 389263 : call kiss_gen%finalize()
281 :
282 389263 : end subroutine mcica_subcol_sw
283 :
284 :
285 : end module mcica_subcol_gen
286 :
|