Line data Source code
1 : module hygrocoreshell_aerosol_optics_mod
2 : use shr_kind_mod, only: r8 => shr_kind_r8
3 : use aerosol_optics_mod, only: aerosol_optics
4 : use aerosol_state_mod, only: aerosol_state
5 : use aerosol_properties_mod, only: aerosol_properties
6 : use table_interp_mod, only: table_interp, table_interp_wghts, table_interp_calcwghts
7 :
8 : implicit none
9 :
10 : private
11 : public :: hygrocoreshell_aerosol_optics
12 :
13 : !> hygrocoreshell_aerosol_optics
14 : !! Table look up implementation of aerosol_optics to parameterize aerosol
15 : !! radiative properties in terms of core mass fraction, black carbon/dust fraction,
16 : !! kappa and relative humidity
17 : type, extends(aerosol_optics) :: hygrocoreshell_aerosol_optics
18 :
19 : real(r8), allocatable :: totalmmr(:,:) ! total mmr of the aerosol
20 : real(r8), allocatable :: corefrac(:,:) ! mass fraction that is core
21 : real(r8), allocatable :: bcdust(:,:) ! mass fraction of bc vs (bc + dust)
22 : real(r8), allocatable :: kappa(:,:) ! hygroscopicity
23 : real(r8), allocatable :: relh(:,:) ! relative humidity
24 :
25 : real(r8), pointer :: sw_hygro_coreshell_ext(:,:,:,:,:) => null() ! short wave extinction table
26 : real(r8), pointer :: sw_hygro_coreshell_ssa(:,:,:,:,:) => null() ! short wave single-scatter albedo table
27 : real(r8), pointer :: sw_hygro_coreshell_asm(:,:,:,:,:) => null() ! short wave asymmetry table
28 : real(r8), pointer :: lw_hygro_coreshell_abs(:,:,:,:,:) => null() ! long wave absorption table
29 :
30 : real(r8), pointer :: tbl_corefrac(:) => null() ! core fraction dimension values
31 : real(r8), pointer :: tbl_bcdust(:) => null() ! bc/(bc + dust) fraction dimension values
32 : real(r8), pointer :: tbl_kap(:) => null() ! hygroscopicity dimension values
33 : real(r8), pointer :: tbl_relh(:) => null() ! relative humidity dimension values
34 :
35 : integer :: nfrac = -1 ! core fraction dimension size
36 : integer :: nbcdust = -1 ! bc/(bc + dust) fraction dimension size
37 : integer :: nkap = -1 ! hygroscopicity dimension size
38 : integer :: nrelh = -1 ! relative humidity dimension size
39 :
40 : contains
41 :
42 : procedure :: sw_props
43 : procedure :: lw_props
44 :
45 : final :: destructor
46 :
47 : end type hygrocoreshell_aerosol_optics
48 :
49 : interface hygrocoreshell_aerosol_optics
50 : procedure :: constructor
51 : end interface hygrocoreshell_aerosol_optics
52 :
53 : contains
54 :
55 : !------------------------------------------------------------------------------
56 : !------------------------------------------------------------------------------
57 0 : function constructor(aero_props, aero_state, ilist, ibin, ncol, nlev, relhum) result(newobj)
58 :
59 : class(aerosol_properties),intent(in) :: aero_props ! aerosol_properties object
60 : class(aerosol_state),intent(in) :: aero_state ! aerosol_state object
61 : integer, intent(in) :: ilist ! climate or a diagnostic list number
62 : integer, intent(in) :: ibin ! bin number
63 : integer, intent(in) :: ncol ! number of columns
64 : integer, intent(in) :: nlev ! number of levels
65 : real(r8),intent(in) :: relhum(ncol,nlev) ! relative humidity
66 :
67 : type(hygrocoreshell_aerosol_optics), pointer :: newobj
68 :
69 : integer :: ierr, nspec
70 : integer :: ilev, ispec, icol
71 :
72 0 : real(r8), pointer :: specmmr(:,:) ! species mass mixing ratio
73 :
74 0 : real(r8) :: coremmr(ncol,nlev)
75 0 : real(r8) :: coredustmmr(ncol,nlev)
76 0 : real(r8) :: corebcmmr(ncol,nlev)
77 0 : real(r8) :: shellmmr(ncol,nlev)
78 0 : real(r8) :: bcdustmmr(ncol,nlev)
79 :
80 : character(len=32) :: spectype ! species type
81 : character(len=32) :: specmorph
82 : real(r8) :: specdens ! species density (kg/m3)
83 :
84 0 : allocate(newobj, stat=ierr)
85 0 : if (ierr/=0) then
86 0 : nullify(newobj)
87 : return
88 : end if
89 :
90 0 : allocate(newobj%totalmmr(ncol,nlev),stat=ierr)
91 0 : if (ierr/=0) then
92 0 : nullify(newobj)
93 0 : return
94 : end if
95 :
96 0 : allocate(newobj%corefrac(ncol,nlev),stat=ierr)
97 0 : if (ierr/=0) then
98 0 : nullify(newobj)
99 0 : return
100 : end if
101 :
102 0 : allocate(newobj%bcdust(ncol,nlev),stat=ierr)
103 0 : if (ierr/=0) then
104 0 : nullify(newobj)
105 0 : return
106 : end if
107 :
108 0 : allocate(newobj%kappa(ncol,nlev),stat=ierr)
109 0 : if (ierr/=0) then
110 0 : nullify(newobj)
111 0 : return
112 : end if
113 :
114 0 : allocate(newobj%relh(ncol,nlev),stat=ierr)
115 0 : if (ierr/=0) then
116 0 : nullify(newobj)
117 0 : return
118 : end if
119 :
120 0 : nspec = aero_props%nspecies(ilist,ibin)
121 :
122 0 : coremmr(:,:) = 0._r8
123 0 : coredustmmr(:,:) = 0._r8
124 0 : corebcmmr(:,:) = 0._r8
125 0 : shellmmr(:,:) = 0._r8
126 :
127 0 : do ispec = 1,nspec
128 :
129 0 : call aero_state%get_ambient_mmr(ilist,ispec,ibin,specmmr)
130 :
131 : call aero_props%get(ibin, ispec, list_ndx=ilist, density=specdens, &
132 0 : spectype=spectype, specmorph=specmorph)
133 :
134 0 : if (trim(specmorph) == 'core') then
135 0 : if (trim(spectype) == 'dust') then
136 0 : coredustmmr(:ncol,:nlev) = coredustmmr(:ncol,:nlev) + specmmr(:ncol,:nlev)
137 : end if
138 0 : if (trim(spectype) == 'black-c') then
139 0 : corebcmmr(:ncol,:nlev) = corebcmmr(:ncol,:nlev) + specmmr(:ncol,:nlev)
140 : end if
141 0 : coremmr(:ncol,:nlev) = coremmr(:ncol,:nlev) + specmmr(:ncol,:nlev)
142 0 : else if (trim(specmorph) == 'shell') then
143 0 : shellmmr(:ncol,:nlev) = shellmmr(:ncol,:nlev) + specmmr(:ncol,:nlev)
144 : else
145 0 : nullify(newobj)
146 0 : return
147 : end if
148 :
149 : end do
150 :
151 0 : newobj%totalmmr(:,:) = coremmr(:,:) + shellmmr(:,:)
152 0 : bcdustmmr(:,:) = corebcmmr(:,:) + coredustmmr(:,:)
153 :
154 0 : do ilev = 1, nlev
155 0 : do icol = 1, ncol
156 :
157 0 : if (newobj%totalmmr(icol,ilev) > 0._r8) then
158 0 : newobj%corefrac(icol,ilev) = coremmr(icol,ilev) / newobj%totalmmr(icol,ilev)
159 : else
160 0 : newobj%corefrac(icol,ilev) = 0._r8
161 : end if
162 0 : newobj%corefrac(icol,ilev) = max(0._r8, min(1.0_r8, newobj%corefrac(icol,ilev)))
163 :
164 0 : if (bcdustmmr(icol,ilev) > 0._r8) then
165 0 : newobj%bcdust(icol,ilev) = corebcmmr(icol,ilev) / bcdustmmr(icol,ilev)
166 : else
167 0 : newobj%bcdust(icol,ilev) = 0._r8
168 : end if
169 0 : newobj%bcdust(icol,ilev) = max(0._r8, min(1.0_r8, newobj%bcdust(icol,ilev)))
170 :
171 : end do
172 : end do
173 :
174 0 : call aero_state%hygroscopicity(ilist, ibin, newobj%kappa)
175 :
176 : call aero_props%optics_params(ilist, ibin, &
177 : corefrac=newobj%tbl_corefrac, kap=newobj%tbl_kap, &
178 : bcdust=newobj%tbl_bcdust, relh=newobj%tbl_relh, &
179 : nfrac=newobj%nfrac, nbcdust=newobj%nbcdust, &
180 0 : nkap=newobj%nkap, nrelh=newobj%nrelh)
181 :
182 0 : newobj%relh(:ncol,:) = relhum(:ncol,:)
183 :
184 : ! long wave optical properties table
185 : call aero_props%optics_params(ilist, ibin, &
186 : sw_hygro_coreshell_ext=newobj%sw_hygro_coreshell_ext, &
187 : sw_hygro_coreshell_ssa=newobj%sw_hygro_coreshell_ssa, &
188 : sw_hygro_coreshell_asm=newobj%sw_hygro_coreshell_asm, &
189 0 : lw_hygro_coreshell_ext=newobj%lw_hygro_coreshell_abs)
190 :
191 0 : end function constructor
192 :
193 : !------------------------------------------------------------------------------
194 : ! returns short wave aerosol optics properties
195 : !------------------------------------------------------------------------------
196 0 : subroutine sw_props(self, ncol, ilev, iwav, pext, pabs, palb, pasm)
197 :
198 : class(hygrocoreshell_aerosol_optics), intent(in) :: self
199 : integer, intent(in) :: ncol ! number of columns
200 : integer, intent(in) :: ilev ! vertical level index
201 : integer, intent(in) :: iwav ! wave length index
202 : real(r8),intent(out) :: pext(ncol) ! parameterized specific extinction (m2/kg)
203 : real(r8),intent(out) :: pabs(ncol) ! parameterized specific absorption (m2/kg)
204 : real(r8),intent(out) :: palb(ncol) ! parameterized single scattering albedo
205 : real(r8),intent(out) :: pasm(ncol) ! parameterized asymmetry factor
206 :
207 : integer :: icol
208 :
209 0 : type(table_interp_wghts) :: rhwghts(ncol)
210 0 : type(table_interp_wghts) :: cfwghts(ncol)
211 0 : type(table_interp_wghts) :: bcwghts(ncol)
212 0 : type(table_interp_wghts) :: kpwghts(ncol)
213 :
214 0 : rhwghts = table_interp_calcwghts( self%nrelh, self%tbl_relh, ncol, self%relh(:ncol,ilev) )
215 0 : cfwghts = table_interp_calcwghts( self%nfrac, self%tbl_corefrac, ncol, self%corefrac(:ncol,ilev) )
216 0 : bcwghts = table_interp_calcwghts( self%nbcdust, self%tbl_bcdust, ncol, self%bcdust(:ncol,ilev) )
217 0 : kpwghts = table_interp_calcwghts( self%nkap, self%tbl_kap, ncol, self%kappa(:ncol,ilev) )
218 :
219 0 : pext = table_interp( ncol, self%nrelh,self%nfrac,self%nbcdust,self%nkap, rhwghts,cfwghts,bcwghts,kpwghts, self%sw_hygro_coreshell_ext(:,iwav,:,:,:))
220 0 : pabs = (1._r8-table_interp( ncol, self%nrelh,self%nfrac,self%nbcdust,self%nkap, rhwghts,cfwghts,bcwghts,kpwghts, self%sw_hygro_coreshell_ssa(:,iwav,:,:,:)))*pext
221 0 : pasm = table_interp( ncol, self%nrelh,self%nfrac,self%nbcdust,self%nkap, rhwghts,cfwghts,bcwghts,kpwghts, self%sw_hygro_coreshell_asm(:,iwav,:,:,:))
222 :
223 0 : do icol = 1, ncol
224 :
225 0 : pext(icol) = pext(icol)*self%totalmmr(icol,ilev)
226 0 : pabs(icol) = pabs(icol)*self%totalmmr(icol,ilev)
227 0 : pabs(icol) = max(0._r8,pabs(icol))
228 0 : pabs(icol) = min(pext(icol),pabs(icol))
229 :
230 0 : palb(icol) = 1._r8-pabs(icol)/max(pext(icol),1.e-40_r8)
231 :
232 : end do
233 :
234 0 : end subroutine sw_props
235 :
236 : !------------------------------------------------------------------------------
237 : ! returns long wave aerosol optics properties
238 : !------------------------------------------------------------------------------
239 0 : subroutine lw_props(self, ncol, ilev, iwav, pabs)
240 :
241 : class(hygrocoreshell_aerosol_optics), intent(in) :: self
242 : integer, intent(in) :: ncol ! number of columns
243 : integer, intent(in) :: ilev ! vertical level index
244 : integer, intent(in) :: iwav ! wave length index
245 : real(r8),intent(out) :: pabs(ncol) ! parameterized specific absorption (m2/kg)
246 :
247 : integer :: icol
248 :
249 0 : type(table_interp_wghts) :: rhwghts(ncol)
250 0 : type(table_interp_wghts) :: cfwghts(ncol)
251 0 : type(table_interp_wghts) :: bcwghts(ncol)
252 0 : type(table_interp_wghts) :: kpwghts(ncol)
253 :
254 0 : rhwghts = table_interp_calcwghts( self%nrelh, self%tbl_relh, ncol, self%relh(:ncol,ilev) )
255 0 : cfwghts = table_interp_calcwghts( self%nfrac, self%tbl_corefrac, ncol, self%corefrac(:ncol,ilev) )
256 0 : bcwghts = table_interp_calcwghts( self%nbcdust, self%tbl_bcdust, ncol, self%bcdust(:ncol,ilev) )
257 0 : kpwghts = table_interp_calcwghts( self%nkap, self%tbl_kap, ncol, self%kappa(:ncol,ilev) )
258 :
259 0 : pabs = table_interp( ncol, self%nrelh,self%nfrac,self%nbcdust,self%nkap, rhwghts,cfwghts,bcwghts,kpwghts, self%lw_hygro_coreshell_abs(:,iwav,:,:,:))
260 :
261 0 : do icol = 1, ncol
262 0 : pabs(icol) = pabs(icol)*self%totalmmr(icol,ilev)
263 0 : pabs(icol) = max(0._r8,pabs(icol))
264 : end do
265 :
266 0 : end subroutine lw_props
267 :
268 : !------------------------------------------------------------------------------
269 : !------------------------------------------------------------------------------
270 0 : subroutine destructor(self)
271 :
272 : type(hygrocoreshell_aerosol_optics), intent(inout) :: self
273 :
274 0 : deallocate(self%totalmmr)
275 0 : deallocate(self%corefrac)
276 0 : deallocate(self%bcdust)
277 0 : deallocate(self%kappa)
278 0 : deallocate(self%relh)
279 :
280 0 : nullify(self%tbl_corefrac)
281 0 : nullify(self%tbl_bcdust)
282 0 : nullify(self%tbl_kap)
283 0 : nullify(self%tbl_relh)
284 0 : nullify(self%sw_hygro_coreshell_ext)
285 0 : nullify(self%sw_hygro_coreshell_ssa)
286 0 : nullify(self%sw_hygro_coreshell_asm)
287 0 : nullify(self%lw_hygro_coreshell_abs)
288 :
289 0 : end subroutine destructor
290 :
291 0 : end module hygrocoreshell_aerosol_optics_mod
|