-
Notifications
You must be signed in to change notification settings - Fork 9
/
carmaelement_mod.F90
323 lines (274 loc) · 14 KB
/
carmaelement_mod.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
!! The CARMAELEMENT module contains configuration information about a particle
!! element used by CARMA.
!!
!! @version March-2010
!! @author Chuck Bardeen
module CARMAELEMENT_mod
use carma_precision_mod
use carma_enums_mod
use carma_constants_mod
use carma_types_mod
! CARMA explicitly declares all variables.
implicit none
! All CARMA variables and procedures are private except those explicitly declared to be public.
private
! Declare the public methods.
public CARMAELEMENT_Create
public CARMAELEMENT_Destroy
public CARMAELEMENT_Get
public CARMAELEMENT_Print
contains
!! Defines a gas used by CARMA for nucleation and growth of cloud and
!! aerosol particles.
!!
!! NOTE: The element density can be specifeid per bin using rhobin; however,
!! if only the bulk density is provided (rho) then the same value will be used
!! for all bins. The bulk density allows for backward compatability and ease of
!! configuration. If rhobin is provided, then rho is ignored.
!!
!! @author Chuck Bardeen
!! @version March-2010
!!
!! @see CARMA_AddGas
!! @see CARMAELEMENT_Destroy
subroutine CARMAELEMENT_Create(carma, ielement, igroup, name, rho, itype, icomposition, rc, &
shortname, isolute, rhobin, arat, kappa, refidx, isShell)
type(carma_type), intent(inout) :: carma !! the carma object
integer, intent(in) :: ielement !! the element index
integer, intent(in) :: igroup !! Group to which the element belongs
character(*), intent(in) :: name !! the element name, maximum of 255 characters
real(kind=f), intent(in) :: rho !! bulk mass density of particle element [g/cm^3]
integer, intent(in) :: itype !! Particle type specification
integer, intent(in) :: icomposition !! Particle compound specification
integer, intent(out) :: rc !! return code, negative indicates failure
character(*), optional, intent(in) :: shortname !! the element shortname, maximum of 6 characters
integer, optional, intent(in) :: isolute !! Index of solute for the particle element
real(kind=f), optional, intent(in) :: rhobin(carma%f_NBIN)!! mass density per bin of particle element [g/cm^3]
real(kind=f), optional, intent(in) :: arat(carma%f_NBIN) !! projected area ratio
real(kind=f), optional, intent(in) :: kappa !! hygroscopicity parameter for the particle element [units?]
complex(kind=f), optional, intent(in) :: refidx(carma%f_NWAVE, carma%f_NREFIDX) !! refractive indices
logical, optional, intent(in) :: isShell !! Is this element part of the shell or the core?
! Local variables
integer :: ier
! Assume success.
rc = RC_OK
! Make sure there are enough elements allocated.
if (ielement > carma%f_NELEM) then
if (carma%f_do_print) write(carma%f_LUNOPRT, *) "CARMAELEMENT_Create:: ERROR - The specifed element (", &
ielement, ") is larger than the number of elements (", carma%f_NELEM, ")."
rc = RC_ERROR
return
end if
! Make sure there are enough groups allocated.
if (igroup > carma%f_NGROUP) then
if (carma%f_do_print) write(carma%f_LUNOPRT, *) "CARMAELEMENT_Create:: ERROR - The specifed group (", &
igroup, ") is larger than the number of groups (", carma%f_NGROUP, ")."
rc = RC_ERROR
return
end if
allocate( &
carma%f_element(ielement)%f_rho(carma%f_NBIN), &
stat=ier)
if(ier /= 0) then
if (carma%f_do_print) write(carma%f_LUNOPRT, *) "CARMAELEMENT_Add: ERROR allocating, status=", ier
rc = RC_ERROR
return
end if
if ((carma%f_NWAVE > 0) .and. (carma%f_NREFIDX > 0)) then
allocate( &
carma%f_element(ielement)%f_refidx(carma%f_NWAVE, carma%f_NREFIDX), &
stat=ier)
if(ier /= 0) then
if (carma%f_do_print) write(carma%f_LUNOPRT, *) "CARMAELEMENT_Add: ERROR allocating, status=", ier
rc = RC_ERROR
return
end if
carma%f_element(ielement)%f_refidx(:,:) = CMPLX(0._f, 0._f, kind=f)
end if
! Save off the settings.
carma%f_element(ielement)%f_igroup = igroup
carma%f_element(ielement)%f_name = name
carma%f_element(ielement)%f_rho(:) = rho
carma%f_element(ielement)%f_itype = itype
carma%f_element(ielement)%f_icomposition = icomposition
! Defaults for optional parameters
carma%f_element(ielement)%f_shortname = ""
carma%f_element(ielement)%f_isolute = 0
carma%f_element(ielement)%f_kappa = 0.0_f
carma%f_element(ielement)%f_isShell = .true.
! Set optional parameters.
if (present(shortname)) carma%f_element(ielement)%f_shortname = shortname
if (present(kappa)) carma%f_element(ielement)%f_kappa = kappa
if (present(refidx)) carma%f_element(ielement)%f_refidx(:,:) = refidx(:,:)
if (present(isShell)) carma%f_element(ielement)%f_isShell = isShell
if (present(isolute)) then
! Make sure there are enough solutes allocated.
if (isolute > carma%f_NSOLUTE) then
if (carma%f_do_print) write(carma%f_LUNOPRT, *) "CARMAELEMENT_Create:: ERROR - The specifed solute (", &
isolute, ") is larger than the number of solutes (", carma%f_NSOLUTE, ")."
rc = RC_ERROR
return
end if
carma%f_element(ielement)%f_isolute = isolute
end if
if (present(rhobin)) carma%f_element(ielement)%f_rho(:) = rhobin(:)
! If the area ratio is specfied (usually along with rhobin), then set this
! for the group.
if (present(arat)) carma%f_group(igroup)%f_arat(:) = arat(:)
! Keep track of the fact that another element has been added to the group.
carma%f_group(igroup)%f_nelem = carma%f_group(igroup)%f_nelem + 1
return
end subroutine CARMAELEMENT_Create
!! Deallocates the memory associated with a CARMAELEMENT object.
!!
!! @author Chuck Bardeen
!! @version March-2010
!!
!! @see CARMAELEMENT_Create
subroutine CARMAELEMENT_Destroy(carma, ielement, rc)
type(carma_type), intent(inout) :: carma !! the carma object
integer, intent(in) :: ielement !! the element index
integer, intent(out) :: rc !! return code, negative indicates failure
! Local variables
integer :: ier
! Assume success.
rc = RC_OK
! Make sure there are enough elements allocated.
if (ielement > carma%f_NELEM) then
if (carma%f_do_print) write(carma%f_LUNOPRT, *) "CARMAELEMENT_Destroy:: ERROR - The specifed element (", &
ielement, ") is larger than the number of elements (", carma%f_NELEM, ")."
rc = RC_ERROR
return
end if
if (allocated(carma%f_element(ielement)%f_rho)) then
deallocate(carma%f_element(ielement)%f_rho,stat=ier)
if(ier /= 0) then
if (carma%f_do_print) then
write(carma%f_LUNOPRT, *) "CARMAELEMENT_Destroy: ERROR deallocating f_rho, status=", ier
endif
rc = RC_ERROR
return
endif
endif
if (allocated(carma%f_element(ielement)%f_refidx)) then
deallocate(carma%f_element(ielement)%f_refidx, stat=ier)
if(ier /= 0) then
if (carma%f_do_print) then
write(carma%f_LUNOPRT, *) "CARMAELEMENT_Destroy: ERROR deallocating f_refidx, status=", ier
endif
rc = RC_ERROR
return
endif
endif
return
end subroutine CARMAELEMENT_Destroy
!! Gets information about a particle element.
!!
!! The group name and other properties are available after a call to
!! CARMAELEMENT_Create().
!!
!! @author Chuck Bardeen
!! @version March-2010
!!
!! @see CARMAELEMENT_Create
!! @see CARMA_GetElement
subroutine CARMAELEMENT_Get(carma, ielement, rc, igroup, name, shortname, rho, itype, icomposition, isolute, kappa, refidx, isShell)
type(carma_type), intent(in) :: carma !! the carma object
integer, intent(in) :: ielement !! the element index
integer, intent(out) :: rc !! return code, negative indicates failure
integer, optional, intent(out) :: igroup !! Group to which the element belongs
character(len=*), optional, intent(out) :: name !! the element name
character(len=*), optional, intent(out) :: shortname !! the element short name
real(kind=f), optional, intent(out) :: rho(carma%f_NBIN) !! Mass density of particle element [g/cm^3]
integer, optional, intent(out) :: itype !! Particle type specification
integer, optional, intent(out) :: icomposition !! Particle compound specification
integer, optional, intent(out) :: isolute !! Index of solute for the particle element
real(kind=f), optional, intent(out) :: kappa !! hygroscopicity parameter for the particle element [units?]
complex(kind=f), optional, intent(out) :: refidx(carma%f_NWAVE, carma%f_NREFIDX) !! Refractive indices
logical, optional, intent(out) :: isShell !! Is this element part of the shell or the core?
! Assume success.
rc = RC_OK
! Make sure there are enough elements allocated.
if (ielement > carma%f_NELEM) then
if (carma%f_do_print) write(carma%f_LUNOPRT, *) "CARMAELEMENT_Get:: ERROR - The specifed element (", &
ielement, ") is larger than the number of elements (", carma%f_NELEM, ")."
rc = RC_ERROR
return
end if
! Return any requested properties of the group.
if (present(igroup)) igroup = carma%f_element(ielement)%f_igroup
if (present(name)) name = carma%f_element(ielement)%f_name
if (present(shortname)) shortname = carma%f_element(ielement)%f_shortname
if (present(rho)) rho(:) = carma%f_element(ielement)%f_rho(:)
if (present(itype)) itype = carma%f_element(ielement)%f_itype
if (present(icomposition)) icomposition = carma%f_element(ielement)%f_icomposition
if (present(isolute)) isolute = carma%f_element(ielement)%f_isolute
if (present(kappa)) kappa = carma%f_element(ielement)%f_kappa
if (present(isShell)) isShell = carma%f_element(ielement)%f_isShell
if ((carma%f_NWAVE == 0) .or. (carma%f_NREFIDX == 0)) then
if (present(refidx)) then
if (carma%f_do_print) write(carma%f_LUNOPRT, *) "CARMAGROUP_Get: ERROR no refidx defined."
rc = RC_ERROR
return
end if
else
if (present(refidx)) refidx(:,:) = carma%f_element(ielement)%f_refidx(:,:)
end if
return
end subroutine CARMAELEMENT_Get
!! Prints information about an element.
!!
!! @author Chuck Bardeen
!! @version March-2010
!!
!! @see CARMAELEMENT_Get
subroutine CARMAELEMENT_Print(carma, ielement, rc)
type(carma_type), intent(in) :: carma !! the carma object
integer, intent(in) :: ielement !! the element index
integer, intent(out) :: rc !! return code, negative indicates failure
! Local variables
character(len=CARMA_NAME_LEN) :: name ! name
character(len=CARMA_SHORT_NAME_LEN) :: shortname ! shortname
real(kind=f) :: rho(carma%f_NBIN)! density (g/cm3)
integer :: igroup ! Group to which the element belongs
integer :: itype ! Particle type specification
integer :: icomposition ! Particle compound specification
integer :: isolute ! Index of solute for the particle element
real(kind=f) :: kappa ! hygroscopicity factor
complex(kind=f) :: refidx(carma%f_NWAVE, carma%f_NREFIDX) ! Refractive indices
logical :: isShell ! Is this element part of the shell or the core?
! Assume success.
rc = RC_OK
! Test out the Get method.
if (carma%f_do_print) then
call CARMAELEMENT_Get(carma, ielement, rc, name=name, shortname=shortname, igroup=igroup, &
itype=itype, icomposition=icomposition, rho=rho, isolute=isolute, &
kappa=kappa, refidx=refidx, isShell=isShell)
if (rc < 0) return
write(carma%f_LUNOPRT,*) " name : ", trim(name)
write(carma%f_LUNOPRT,*) " igroup : ", igroup
write(carma%f_LUNOPRT,*) " shortname : ", trim(shortname)
write(carma%f_LUNOPRT,*) " rho : ", rho, " (g/cm3)"
select case(itype)
case (I_INVOLATILE)
write(carma%f_LUNOPRT,*) " itype : involatile"
case (I_VOLATILE)
write(carma%f_LUNOPRT,*) " itype : volatile"
case (I_COREMASS)
write(carma%f_LUNOPRT,*) " itype : core mass"
case (I_VOLCORE)
write(carma%f_LUNOPRT,*) " itype : volatile core"
case (I_CORE2MOM)
write(carma%f_LUNOPRT,*) " itype : core mass - second moment"
case default
write(carma%f_LUNOPRT,*) " itype : unknown, ", itype
end select
write(carma%f_LUNOPRT,*) " icomposition : ", icomposition
write(carma%f_LUNOPRT,*) " isolute : ", isolute
write(carma%f_LUNOPRT,*) " kappa : ", kappa
write(carma%f_LUNOPRT,*) " isShell : ", isShell
write(carma%f_LUNOPRT,*) " ref. index : ", refidx
end if
return
end subroutine CARMAELEMENT_Print
end module CARMAELEMENT_mod