37 subroutine mosaic_init(env_state, aero_data, del_t, do_optical)
40 use module_data_mosaic_aero
, only: alpha_astem, rtol_eqb_astem, &
41 ptol_mol_astem, mgas_aer_xfer, mdynamic_solver
43 use module_data_mosaic_main
, only: tbeg_sec, dt_sec, rlon, rlat, &
44 zalt_m, rh, te, pr_atm, cair_mlc, cair_molm3, ppb, avogad, &
45 mmode, mgas, maer, mcld, maeroptic, mshellcore, &
46 msolar, mphoto, lun_aeroptic, naerbin
54 real(kind=dp),
intent(in) :: del_t
56 logical,
intent(in) :: do_optical
61 subroutine loadperoxyparameters()
62 end subroutine loadperoxyparameters
63 subroutine init_data_modules()
64 end subroutine init_data_modules
65 subroutine allocatememory()
66 end subroutine allocatememory
69 call init_data_modules
92 rtol_eqb_astem = 0.01d0
93 ptol_mol_astem = 0.01d0
97 tbeg_sec = env_state%start_day*24*3600 + &
98 nint(env_state%start_time)
101 rlon =
deg2rad(env_state%longitude)
102 rlat =
deg2rad(env_state%latitude)
103 zalt_m = env_state%altitude
106 rh = env_state%rel_humid * 100.d0
108 pr_atm = env_state%pressure / const%air_std_press
109 cair_mlc = avogad*pr_atm/(82.056d0*te)
110 cair_molm3 = 1d6*pr_atm/(82.056d0*te)
113 call loadperoxyparameters
116 if (lun_aeroptic <= 0 ) lun_aeroptic =
get_unit()
119 call
assert_msg(111041803, aero_data%i_water > 0, &
120 "MOSAIC requires H2O as an aerosol species")
131 #ifdef PMC_USE_MOSAIC
134 subroutine deallocatememory()
135 end subroutine deallocatememory
138 call deallocatememory()
147 aero_state, gas_data, gas_state)
149 #ifdef PMC_USE_MOSAIC
150 use module_data_mosaic_aero
, only: nbin_a, aer, num_a, jhyst_leg, &
153 use module_data_mosaic_main
, only: tbeg_sec, tcur_sec, tmid_sec, &
154 dt_sec, dt_min, dt_aeroptic_min, rh, te, pr_atm, cnn, cair_mlc, &
155 cair_molm3, ppb, avogad, msolar, naerbin
169 #ifdef PMC_USE_MOSAIC
171 real(kind=dp) :: time_utc
172 real(kind=dp) :: tmar21_sec
173 real(kind=dp) :: conv_fac(aero_data%n_spec), dum_var
174 integer :: i_part, i_spec, i_spec_mosaic
176 real(kind=dp) :: num_conc
180 subroutine allocatememory()
181 end subroutine allocatememory
182 subroutine deallocatememory()
183 end subroutine deallocatememory
187 tmar21_sec =
real((79*24 + 12)*3600, kind=dp)
188 tcur_sec =
real(tbeg_sec, kind=dp) + env_state%elapsed_time
191 time_utc = env_state%start_time/3600d0
192 time_utc = time_utc + dt_sec/3600d0
194 do while (time_utc >= 24d0)
195 time_utc = time_utc - 24d0
198 tmid_sec = tcur_sec + 0.5d0*dt_sec
199 if(tmid_sec .ge. tmar21_sec)
then
200 tmid_sec = tmid_sec - tmar21_sec
202 tmid_sec = tmid_sec &
203 + dble(((365-79)*24 - 12)*3600)
209 dt_aeroptic_min = 0d0
212 do i_spec = 1,aero_data%n_spec
214 conv_fac(i_spec) = 1.d9 * aero_data%density(i_spec) &
215 / aero_data%molec_weight(i_spec)
219 rh = env_state%rel_humid * 100.d0
221 pr_atm = env_state%pressure / const%air_std_press
222 cair_mlc = avogad*pr_atm/(82.056d0*te)
223 cair_molm3 = 1d6*pr_atm/(82.056d0*te)
228 if (nbin_a > naerbin)
then
229 call deallocatememory()
231 call allocatememory()
236 do i_part = aero_state%apa%n_part,1,-1
237 particle => aero_state%apa%particle(i_part)
239 do i_spec = 1,aero_data%n_spec
240 i_spec_mosaic = aero_data%mosaic_index(i_spec)
241 if (i_spec_mosaic > 0)
then
243 aer(i_spec_mosaic, 3, i_part) &
244 = particle%vol(i_spec) * conv_fac(i_spec) * num_conc
249 water_a(i_part) = particle%vol(aero_data%i_water) &
250 * aero_data%density(aero_data%i_water) * num_conc
251 num_a(i_part) = 1d-6 * num_conc
252 jhyst_leg(i_part) = particle%water_hyst_leg
257 do i_spec = 1,gas_data%n_spec
258 i_spec_mosaic = gas_data%mosaic_index(i_spec)
259 if (i_spec_mosaic > 0)
then
261 cnn(i_spec_mosaic) = gas_state%mix_rat(i_spec) * cair_mlc / ppb
274 #ifdef PMC_USE_MOSAIC
275 use module_data_mosaic_aero
, only: nbin_a, aer, num_a, jhyst_leg, &
278 use module_data_mosaic_main
, only: tbeg_sec, tcur_sec, tmid_sec, &
279 dt_sec, dt_min, dt_aeroptic_min, rh, te, pr_atm, cnn, cair_mlc, &
280 cair_molm3, ppb, avogad, msolar, cos_sza
294 #ifdef PMC_USE_MOSAIC
296 real(kind=dp) :: conv_fac(aero_data%n_spec), dum_var, num_conc
297 integer :: i_part, i_spec, i_spec_mosaic
299 real(kind=dp) :: reweight_num_conc(aero_state%apa%n_part)
302 do i_spec = 1,aero_data%n_spec
304 conv_fac(i_spec) = 1d9 * aero_data%density(i_spec) &
305 / aero_data%molec_weight(i_spec)
309 env_state%rel_humid = rh / 100d0
311 env_state%pressure = pr_atm * const%air_std_press
312 if (msolar == 1)
then
313 env_state%solar_zenith_angle = acos(cos_sza)
315 cair_mlc = avogad*pr_atm/(82.056d0*te)
316 cair_molm3 = 1d6*pr_atm/(82.056d0*te)
320 aero_state%valid_sort = .false.
324 do i_part = 1,aero_state%apa%n_part,1
325 particle => aero_state%apa%particle(i_part)
327 do i_spec = 1,aero_data%n_spec
328 i_spec_mosaic = aero_data%mosaic_index(i_spec)
329 if (i_spec_mosaic > 0)
then
330 particle%vol(i_spec) = &
332 aer(i_spec_mosaic, 3, i_part) &
333 / (conv_fac(i_spec) * num_conc)
336 particle%water_hyst_leg = jhyst_leg(i_part)
339 particle%vol(aero_data%i_water) = water_a(i_part) &
340 / aero_data%density(aero_data%i_water) / num_conc
346 do i_spec = 1,gas_data%n_spec
347 i_spec_mosaic = gas_data%mosaic_index(i_spec)
348 if (i_spec_mosaic > 0)
then
350 gas_state%mix_rat(i_spec) = cnn(i_spec_mosaic) / cair_mlc * ppb
367 gas_state, do_optical)
369 #ifdef PMC_USE_MOSAIC
370 use module_data_mosaic_main
, only: msolar
384 logical,
intent(in) :: do_optical
386 #ifdef PMC_USE_MOSAIC
389 subroutine solarzenithangle()
390 end subroutine solarzenithangle
391 subroutine integratechemistry()
392 end subroutine integratechemistry
393 subroutine aerosol_optical()
394 end subroutine aerosol_optical
401 if (msolar == 1)
then
402 call solarzenithangle
405 call integratechemistry
413 aero_state, gas_data, gas_state)
431 aero_state, gas_data, gas_state)
433 #ifdef PMC_USE_MOSAIC
434 use module_data_mosaic_aero
, only: ri_shell_a, ri_core_a, &
435 ext_cross, scat_cross, asym_particle
449 #ifdef PMC_USE_MOSAIC
452 subroutine aerosol_optical()
453 end subroutine aerosol_optical
468 do i_part = aero_state%apa%n_part,1,-1
469 particle => aero_state%apa%particle(i_part)
470 particle%absorb_cross_sect = (ext_cross(i_part) &
471 - scat_cross(i_part)) / 1d4
472 particle%scatter_cross_sect = scat_cross(i_part) / 1d4
473 particle%asymmetry = asym_particle(i_part)
474 particle%refract_shell = cmplx(ri_shell_a(i_part), kind=dc)
475 particle%refract_core = cmplx(ri_core_a(i_part), kind=dc)
478 particle%core_vol = 0d0
The aero_data_t structure and associated subroutines.
The env_state_t structure and associated subroutines.
subroutine aero_state_num_conc_for_reweight(aero_state, reweight_num_conc)
Save the correct number concentrations for later use by aero_state_reweight().
subroutine assert_msg(code, condition_ok, error_msg)
Errors unless condition_ok is true.
subroutine mosaic_to_partmc(env_state, aero_data, aero_state, gas_data, gas_state)
Map all data MOSAIC -> PartMC.
real(kind=dp) function aero_weight_array_num_conc(aero_weight_array, aero_particle)
Compute the number concentration for a particle (m^{-3}).
The gas_data_t structure and associated subroutines.
Interface to the MOSAIC aerosol and gas phase chemistry code.
Common utility subroutines.
Current environment state.
logical function mosaic_support()
Whether MOSAIC support is compiled in.
The aero_state_t structure and assocated subroutines.
subroutine mosaic_timestep(env_state, aero_data, aero_state, gas_data, gas_state, do_optical)
Do one timestep with MOSAIC.
real(kind=dp) function deg2rad(deg)
Convert degrees to radians.
The gas_state_t structure and associated subroutines.
integer function aero_state_total_particles(aero_state, i_group, i_class)
Returns the total number of particles in an aerosol distribution.
subroutine aero_state_reweight(aero_state, reweight_num_conc)
Reweight all particles after their constituent volumes have been altered.
The current collection of aerosol particles.
Single aerosol particle data structure.
subroutine mosaic_aero_optical(env_state, aero_data, aero_state, gas_data, gas_state)
Compute the optical properties of each aerosol particle. FIXME: currently disabled.
integer function get_unit()
Returns an available unit number. This should be freed by free_unit().
subroutine mosaic_from_partmc(env_state, aero_data, aero_state, gas_data, gas_state)
Map all data PartMC -> MOSAIC.
subroutine mosaic_init(env_state, aero_data, del_t, do_optical)
Initialize all MOSAIC data-structures.
Current state of the gas mixing ratios in the system.
subroutine mosaic_cleanup()
Clean-up after running MOSAIC, deallocating memory.
Aerosol material properties and associated data.