37 #ifdef PMC_USE_SUNDIALS
43 logical,
parameter :: CONDENSE_DO_TEST_JAC_SOLVE = .false.
46 logical,
parameter :: CONDENSE_DO_TEST_COUNTS = .false.
49 integer,
parameter :: PMC_CONDENSE_SOLVER_SUCCESS = 0
51 integer,
parameter :: PMC_CONDENSE_SOLVER_INIT_Y = 1
53 integer,
parameter :: PMC_CONDENSE_SOLVER_INIT_ABSTOL = 2
55 integer,
parameter :: PMC_CONDENSE_SOLVER_INIT_CVODE_MEM = 3
57 integer,
parameter :: PMC_CONDENSE_SOLVER_INIT_CVODE = 4
59 integer,
parameter :: PMC_CONDENSE_SOLVER_SVTOL = 5
61 integer,
parameter :: PMC_CONDENSE_SOLVER_SET_MAX_STEPS = 6
63 integer,
parameter :: PMC_CONDENSE_SOLVER_FAIL = 7
79 real(kind=dp) :: V_comp
83 real(kind=dp) :: D_dry
85 real(kind=dp) :: kappa
94 real(kind=dp) :: Hdot_i
96 real(kind=dp) :: Hdot_env
98 real(kind=dp) :: dDdot_dD
100 real(kind=dp) :: dDdot_dH
102 real(kind=dp) :: dHdoti_dD
104 real(kind=dp) :: dHdoti_dH
106 real(kind=dp) :: dHdotenv_dD
108 real(kind=dp) :: dHdotenv_dH
119 real(kind=dp) :: condense_saved_Tdot
122 real(kind=dp) :: condense_saved_pdot
125 real(kind=dp),
allocatable :: condense_saved_kappa(:)
128 real(kind=dp),
allocatable :: condense_saved_D_dry(:)
131 real(kind=dp),
allocatable :: condense_saved_num_conc(:)
135 integer,
save :: condense_count_vf
138 integer,
save :: condense_count_solve
148 env_state_final, del_t)
159 type(env_state_t),
intent(inout) :: env_state_final
161 real(kind=dp),
intent(in) :: del_t
163 integer :: i_part, n_eqn, i_eqn
165 real(kind=dp) :: state(aero_state%apa%n_part + 1), init_time, final_time
166 real(kind=dp) :: abs_tol_vector(aero_state%apa%n_part + 1)
167 real(kind=dp) :: num_conc
168 real(kind=dp) :: reweight_num_conc(aero_state%apa%n_part)
169 real(kind=dp) :: water_vol_conc_initial, water_vol_conc_final
170 real(kind=dp) :: vapor_vol_conc_initial, vapor_vol_conc_final
171 real(kind=dp) :: d_water_vol_conc, d_vapor_vol_conc
172 real(kind=dp) :: v_comp_ratio, water_rel_error
173 #ifdef PMC_USE_SUNDIALS
174 real(kind=c_double),
target :: state_f(aero_state%apa%n_part + 1)
175 real(kind=c_double),
target :: abstol_f(aero_state%apa%n_part + 1)
176 type(c_ptr
) :: state_f_p, abstol_f_p
177 integer(kind=c_int) :: n_eqn_f, solver_stat
178 real(kind=c_double) :: reltol_f, t_initial_f, t_final_f
181 #ifdef PMC_USE_SUNDIALS
182 #ifndef DOXYGEN_SKIP_DOC
185 reltol_f, t_initial_f, t_final_f) bind(c)
187 integer(kind=c_int), value :: neq
188 type(c_ptr
), value :: x_f
189 type(c_ptr
), value :: abstol_f
190 real(kind=c_double), value :: reltol_f
191 real(kind=c_double), value :: t_initial_f
192 real(kind=c_double), value :: t_final_f
199 water_vol_conc_initial = 0d0
200 do i_part = 1,aero_state%apa%n_part
201 aero_particle => aero_state%apa%particle(i_part)
203 water_vol_conc_initial = water_vol_conc_initial &
204 + aero_particle%vol(aero_data%i_water) * num_conc
211 call
env_state_copy(env_state_initial, condense_saved_env_state_initial)
212 condense_saved_tdot &
213 = (env_state_final%temp - env_state_initial%temp) / del_t
214 condense_saved_pdot &
215 = (env_state_final%pressure - env_state_initial%pressure) / del_t
218 allocate(condense_saved_kappa(aero_state%apa%n_part))
219 allocate(condense_saved_d_dry(aero_state%apa%n_part))
220 allocate(condense_saved_num_conc(aero_state%apa%n_part))
224 do i_part = aero_state%apa%n_part,1,-1
225 aero_particle => aero_state%apa%particle(i_part)
226 condense_saved_kappa(i_part) &
228 condense_saved_d_dry(i_part) =
vol2diam(&
230 condense_saved_num_conc(i_part) &
233 abs_tol_vector(i_part) = max(1d-30, &
234 1d-8 * (state(i_part) - condense_saved_d_dry(i_part)))
236 state(aero_state%apa%n_part + 1) = env_state_initial%rel_humid
237 abs_tol_vector(aero_state%apa%n_part + 1) = 1d-10
239 #ifdef PMC_USE_SUNDIALS
241 n_eqn = aero_state%apa%n_part + 1
242 n_eqn_f = int(n_eqn, kind=c_int)
243 reltol_f =
real(1d-8, kind=c_double)
244 t_initial_f =
real(0, kind=c_double)
245 t_final_f =
real(del_t, kind=c_double)
247 state_f(i_eqn) =
real(state(i_eqn), kind=c_double)
248 abstol_f(i_eqn) =
real(abs_tol_vector(i_eqn), kind=c_double)
250 state_f_p = c_loc(state_f)
251 abstol_f_p = c_loc(abstol_f)
252 condense_count_vf = 0
253 condense_count_solve = 0
254 solver_stat =
condense_solver(n_eqn_f, state_f_p, abstol_f_p, reltol_f, &
255 t_initial_f, t_final_f)
256 call condense_check_solve(solver_stat)
257 if (condense_do_test_counts)
then
258 write(0,*)
'condense_count_vf ', condense_count_vf
259 write(0,*)
'condense_count_solve ', condense_count_solve
262 state(i_eqn) =
real(state_f(i_eqn), kind=dp)
267 env_state_final%rel_humid = state(aero_state%apa%n_part + 1)
272 water_vol_conc_final = 0d0
274 do i_part = 1,aero_state%apa%n_part
275 aero_particle => aero_state%apa%particle(i_part)
279 aero_particle%vol(aero_data%i_water) =
diam2vol(state(i_part)) &
283 aero_particle%vol(aero_data%i_water) = max(0d0, &
284 aero_particle%vol(aero_data%i_water))
287 water_vol_conc_final = water_vol_conc_final &
288 + aero_particle%vol(aero_data%i_water) * num_conc
297 v_comp_ratio = env_state_final%temp * env_state_initial%pressure &
298 / (env_state_initial%temp * env_state_final%pressure)
299 vapor_vol_conc_initial = aero_data%molec_weight(aero_data%i_water) &
300 / (const%univ_gas_const * env_state_initial%temp) &
302 * env_state_initial%rel_humid &
304 vapor_vol_conc_final = aero_data%molec_weight(aero_data%i_water) &
305 / (const%univ_gas_const * env_state_final%temp) &
307 * env_state_final%rel_humid &
309 d_vapor_vol_conc = vapor_vol_conc_final - vapor_vol_conc_initial
310 d_water_vol_conc = water_vol_conc_final - water_vol_conc_initial
311 water_rel_error = (d_vapor_vol_conc + d_water_vol_conc) &
312 / (vapor_vol_conc_final + water_vol_conc_final)
314 "condensation water imbalance too high: " &
317 deallocate(condense_saved_kappa)
318 deallocate(condense_saved_d_dry)
319 deallocate(condense_saved_num_conc)
327 #ifdef PMC_USE_SUNDIALS
329 subroutine condense_check_solve(value)
332 integer(kind=c_int),
intent(in) :: value
334 if (value == pmc_condense_solver_success)
then
336 elseif (value == pmc_condense_solver_init_y)
then
337 call
die_msg(123749472,
"condense_solver: " &
338 //
"failed to allocate y vector")
339 elseif (value == pmc_condense_solver_init_abstol)
then
340 call
die_msg(563665949,
"condense_solver: " &
341 //
"failed to allocate abstol vector")
342 elseif (value == pmc_condense_solver_init_cvode_mem)
then
343 call
die_msg(700541443,
"condense_solver: " &
344 //
"failed to create the solver")
345 elseif (value == pmc_condense_solver_init_cvode)
then
346 call
die_msg(297559183,
"condense_solver: " &
347 //
"failure to initialize the solver")
348 elseif (value == pmc_condense_solver_svtol)
then
349 call
die_msg(848342417,
"condense_solver: " &
350 //
"failed to set tolerances")
351 elseif (value == pmc_condense_solver_set_max_steps)
then
352 call
die_msg(275591501,
"condense_solver: " &
353 //
"failed to set maximum steps")
354 elseif (value == pmc_condense_solver_fail)
then
355 call
die_msg(862254233,
"condense_solver: solver failed")
357 call
die_msg(635697577,
"condense_solver: unknown return code: " &
361 end subroutine condense_check_solve
376 real(kind=dp) :: rho_w, m_w, p_0, dp0_dt_div_p0, rho_air, k_a, d_v, u
377 real(kind=dp) :: v, w, x, y, z, k_ap, dkap_dd, d_vp, ddvp_dd
378 real(kind=dp) :: a_w, daw_dd, delta_star, h, dh_ddelta, dh_dd
379 real(kind=dp) :: dh_dh, ddeltastar_dd, ddeltastar_dh
380 integer :: newton_step
382 rho_w = const%water_density
383 m_w = const%water_molec_weight
384 p_0 = const%water_eq_vap_press &
385 * 10d0**(7.45d0 * (inputs%T - const%water_freeze_temp) &
387 dp0_dt_div_p0 = 7.45d0 * log(10d0) * (const%water_freeze_temp - 38d0) &
388 / (inputs%T - 38d0)**2
389 rho_air = const%air_molec_weight * inputs%p &
390 / (const%univ_gas_const * inputs%T)
392 k_a = 1d-3 * (4.39d0 + 0.071d0 * inputs%T)
393 d_v = 0.211d-4 / (inputs%p / const%air_std_press) &
394 * (inputs%T / 273d0)**1.94d0
395 u = const%water_latent_heat * rho_w / (4d0 * inputs%T)
396 v = 4d0 * m_w * p_0 / (rho_w * const%univ_gas_const * inputs%T)
397 w = const%water_latent_heat * m_w / (const%univ_gas_const * inputs%T)
398 x = 4d0 * m_w * const%water_surf_eng &
399 / (const%univ_gas_const * inputs%T * rho_w)
400 y = 2d0 * k_a / (const%accom_coeff * rho_air &
401 * const%air_spec_heat) &
402 * sqrt(2d0 * const%pi * const%air_molec_weight &
403 / (const%univ_gas_const * inputs%T))
404 z = 2d0 * d_v / const%accom_coeff * sqrt(2d0 * const%pi * m_w &
405 / (const%univ_gas_const * inputs%T))
407 outputs%Hdot_env = - dp0_dt_div_p0 * inputs%Tdot * inputs%H &
408 + inputs%H * inputs%pdot / inputs%p
409 outputs%dHdotenv_dD = 0d0
410 outputs%dHdotenv_dH = - dp0_dt_div_p0 * inputs%Tdot &
411 + inputs%pdot / inputs%p
413 if (inputs%D <= inputs%D_dry)
then
414 k_ap = k_a / (1d0 + y / inputs%D_dry)
416 d_vp = d_v / (1d0 + z / inputs%D_dry)
421 delta_star = u * v * d_vp * inputs%H / k_ap
423 outputs%Ddot = k_ap * delta_star / (u * inputs%D_dry)
424 outputs%Hdot_i = - 2d0 * const%pi / (v * inputs%V_comp) &
425 * inputs%D_dry**2 * outputs%Ddot
429 dh_dh = - u * v * d_vp
431 ddeltastar_dd = - dh_dd / dh_ddelta
432 ddeltastar_dh = - dh_dh / dh_ddelta
434 outputs%dDdot_dD = 0d0
435 outputs%dDdot_dH = k_ap / (u * inputs%D_dry) * ddeltastar_dh
436 outputs%dHdoti_dD = - 2d0 * const%pi / (v * inputs%V_comp) &
437 * inputs%D_dry**2 * outputs%dDdot_dD
438 outputs%dHdoti_dH = - 2d0 * const%pi / (v * inputs%V_comp) &
439 * inputs%D_dry**2 * outputs%dDdot_dH
444 k_ap = k_a / (1d0 + y / inputs%D)
445 dkap_dd = k_a * y / (inputs%D + y)**2
446 d_vp = d_v / (1d0 + z / inputs%D)
447 ddvp_dd = d_v * z / (inputs%D + z)**2
448 a_w = (inputs%D**3 - inputs%D_dry**3) &
449 / (inputs%D**3 + (inputs%kappa - 1d0) * inputs%D_dry**3)
450 daw_dd = 3d0 * inputs%D**2 * inputs%kappa * inputs%D_dry**3 &
451 / (inputs%D**3 + (inputs%kappa - 1d0) * inputs%D_dry**3)**2
459 delta_star = delta_star - h / dh_ddelta
460 h = k_ap * delta_star - u * v * d_vp &
461 * (inputs%H - a_w / (1d0 + delta_star) &
462 * exp(w * delta_star / (1d0 + delta_star) &
463 + (x / inputs%D) / (1d0 + delta_star)))
465 k_ap - u * v * d_vp * a_w / (1d0 + delta_star)**2 &
466 * (1d0 - w / (1d0 + delta_star) &
467 + (x / inputs%D) / (1d0 + delta_star)) &
468 * exp(w * delta_star / (1d0 + delta_star) &
469 + (x / inputs%D) / (1d0 + delta_star))
472 abs(h) < 1d3 * epsilon(1d0) * abs(u * v * d_vp * inputs%H), &
473 "condensation newton loop did not satisfy convergence tolerance")
475 outputs%Ddot = k_ap * delta_star / (u * inputs%D)
476 outputs%Hdot_i = - 2d0 * const%pi / (v * inputs%V_comp) &
477 * inputs%D**2 * outputs%Ddot
479 dh_dd = dkap_dd * delta_star &
480 - u * v * ddvp_dd * inputs%H + u * v &
481 * (a_w * ddvp_dd + d_vp * daw_dd &
482 - d_vp * a_w * (x / inputs%D**2) / (1d0 + delta_star)) &
483 * (1d0 / (1d0 + delta_star)) &
484 * exp((w * delta_star) / (1d0 + delta_star) &
485 + (x / inputs%D) / (1d0 + delta_star))
486 dh_dh = - u * v * d_vp
488 ddeltastar_dd = - dh_dd / dh_ddelta
489 ddeltastar_dh = - dh_dh / dh_ddelta
491 outputs%dDdot_dD = dkap_dd * delta_star / (u * inputs%D) &
492 + k_ap * ddeltastar_dd / (u * inputs%D) &
493 - k_ap * delta_star / (u * inputs%D**2)
494 outputs%dDdot_dH = k_ap / (u * inputs%D) * ddeltastar_dh
495 outputs%dHdoti_dD = - 2d0 * const%pi / (v * inputs%V_comp) &
496 * (2d0 * inputs%D * outputs%Ddot + inputs%D**2 * outputs%dDdot_dD)
497 outputs%dHdoti_dH = - 2d0 * const%pi / (v * inputs%V_comp) &
498 * inputs%D**2 * outputs%dDdot_dH
504 #ifdef PMC_USE_SUNDIALS
507 subroutine condense_vf_f(n_eqn, time, state_p, state_dot_p) bind(c)
510 integer(kind=c_int), value,
intent(in) :: n_eqn
512 real(kind=c_double), value,
intent(in) :: time
514 type(c_ptr
), value,
intent(in) :: state_p
516 type(c_ptr
), value,
intent(in) :: state_dot_p
518 real(kind=c_double),
pointer :: state(:)
519 real(kind=c_double),
pointer :: state_dot(:)
520 real(kind=dp) :: hdot
525 condense_count_vf = condense_count_vf + 1
527 call c_f_pointer(state_p, state, (/ n_eqn /))
528 call c_f_pointer(state_dot_p, state_dot, (/ n_eqn /))
530 inputs%T = condense_saved_env_state_initial%temp &
531 + time * condense_saved_tdot
532 inputs%p = condense_saved_env_state_initial%pressure &
533 + time * condense_saved_pdot
534 inputs%Tdot = condense_saved_tdot
535 inputs%pdot = condense_saved_pdot
536 inputs%H = state(n_eqn)
539 do i_part = 1,(n_eqn - 1)
540 inputs%D = state(i_part)
541 inputs%D_dry = condense_saved_d_dry(i_part)
542 inputs%V_comp = (inputs%T &
543 * condense_saved_env_state_initial%pressure) &
544 / (condense_saved_env_state_initial%temp * inputs%p) &
545 / condense_saved_num_conc(i_part)
546 inputs%kappa = condense_saved_kappa(i_part)
548 state_dot(i_part) = outputs%Ddot
549 hdot = hdot + outputs%Hdot_i
551 hdot = hdot + outputs%Hdot_env
553 state_dot(n_eqn) = hdot
555 end subroutine condense_vf_f
560 #ifdef PMC_USE_SUNDIALS
564 subroutine condense_jac(n_eqn, time, state_p, dDdot_dD, dDdot_dH, &
568 integer(kind=c_int),
intent(in) :: n_eqn
570 real(kind=c_double),
intent(in) :: time
572 type(c_ptr
),
intent(in) :: state_p
574 real(kind=dp),
intent(out) :: dddot_dd(n_eqn - 1)
576 real(kind=dp),
intent(out) :: dddot_dh(n_eqn - 1)
578 real(kind=dp),
intent(out) :: dhdot_dd(n_eqn - 1)
580 real(kind=dp),
intent(out) :: dhdot_dh
582 real(kind=c_double),
pointer :: state(:)
587 call c_f_pointer(state_p, state, (/ n_eqn /))
589 inputs%T = condense_saved_env_state_initial%temp &
590 + time * condense_saved_tdot
591 inputs%p = condense_saved_env_state_initial%pressure &
592 + time * condense_saved_pdot
593 inputs%Tdot = condense_saved_tdot
594 inputs%pdot = condense_saved_pdot
595 inputs%H = state(n_eqn)
598 do i_part = 1,(n_eqn - 1)
599 inputs%D = state(i_part)
600 inputs%D_dry = condense_saved_d_dry(i_part)
601 inputs%V_comp = (inputs%T &
602 * condense_saved_env_state_initial%pressure) &
603 / (condense_saved_env_state_initial%temp * inputs%p) &
604 / condense_saved_num_conc(i_part)
605 inputs%kappa = condense_saved_kappa(i_part)
607 dddot_dd(i_part) = outputs%dDdot_dD
608 dddot_dh(i_part) = outputs%dDdot_dH
609 dhdot_dd(i_part) = outputs%dHdoti_dD + outputs%dHdotenv_dD
610 dhdot_dh = dhdot_dh + outputs%dHdoti_dH
612 dhdot_dh = dhdot_dh + outputs%dHdotenv_dH
614 end subroutine condense_jac
619 #ifdef PMC_USE_SUNDIALS
623 subroutine condense_jac_solve_f(n_eqn, time, state_p, state_dot_p, &
624 rhs_p, gamma) bind(c)
627 integer(kind=c_int), value,
intent(in) :: n_eqn
629 real(kind=c_double), value,
intent(in) :: time
631 type(c_ptr
), value,
intent(in) :: state_p
633 type(c_ptr
), value,
intent(in) :: state_dot_p
635 type(c_ptr
), value,
intent(in) :: rhs_p
637 real(kind=c_double), value,
intent(in) :: gamma
639 real(kind=c_double),
pointer :: state(:), state_dot(:), rhs(:)
640 real(kind=c_double) :: soln(n_eqn)
641 real(kind=dp) :: dddot_dd(n_eqn - 1), dddot_dh(n_eqn - 1)
642 real(kind=dp) :: dhdot_dd(n_eqn - 1), dhdot_dh
643 real(kind=dp) :: lhs_n, rhs_n
644 real(kind=c_double) :: residual(n_eqn)
645 real(kind=dp) :: rhs_norm, soln_norm, residual_norm
648 condense_count_solve = condense_count_solve + 1
650 call condense_jac(n_eqn, time, state_p, dddot_dd, dddot_dh, &
653 call c_f_pointer(state_p, state, (/ n_eqn /))
654 call c_f_pointer(state_dot_p, state_dot, (/ n_eqn /))
655 call c_f_pointer(rhs_p, rhs, (/ n_eqn /))
658 lhs_n = 1d0 - gamma * dhdot_dh
660 do i_part = 1,(n_eqn - 1)
661 lhs_n = lhs_n - (- gamma * dddot_dh(i_part)) &
662 * (- gamma * dhdot_dd(i_part)) / (1d0 - gamma * dddot_dd(i_part))
663 rhs_n = rhs_n - (- gamma * dhdot_dd(i_part)) * rhs(i_part) &
664 / (1d0 - gamma * dddot_dd(i_part))
666 soln(n_eqn) = rhs_n / lhs_n
668 do i_part = 1,(n_eqn - 1)
669 soln(i_part) = (rhs(i_part) &
670 - (- gamma * dddot_dh(i_part)) * soln(n_eqn)) &
671 / (1d0 - gamma * dddot_dd(i_part))
674 if (condense_do_test_jac_solve)
then
678 residual(n_eqn) = sum(dhdot_dd * soln(1:(n_eqn-1))) &
679 + dhdot_dh * soln(n_eqn)
680 residual(1:(n_eqn-1)) = dddot_dd * soln(1:(n_eqn-1)) &
681 + dddot_dh * soln(n_eqn)
683 residual = rhs - (soln - gamma * residual)
684 rhs_norm = sqrt(sum(rhs**2))
685 soln_norm = sqrt(sum(soln**2))
686 residual_norm = sqrt(sum(residual**2))
687 write(0,*)
'rhs, soln, residual, residual/rhs = ', &
688 rhs_norm, soln_norm, residual_norm, residual_norm / rhs_norm
693 end subroutine condense_jac_solve_f
709 real(kind=dp) :: x, kappa, d_dry, d, g, dg_dd, a_w, daw_dd
710 integer :: newton_step
712 x = 4d0 * const%water_molec_weight * const%water_surf_eng &
713 / (const%univ_gas_const * env_state%temp &
714 * const%water_density)
721 do newton_step = 1,20
723 a_w = (d**3 - d_dry**3) / (d**3 + (kappa - 1d0) * d_dry**3)
724 daw_dd = 3d0 * d**2 * kappa * d_dry**3 &
725 / (d**3 + (kappa - 1d0) * d_dry**3)**2
726 g = env_state%rel_humid - a_w * exp(x / d)
727 dg_dd = - daw_dd * exp(x / d) + a_w * exp(x / d) * (x / d**2)
730 "convergence problem in equilibriation")
751 real(kind=dp) :: reweight_num_conc(aero_state%apa%n_part)
754 aero_state%valid_sort = .false.
757 do i_part = aero_state%apa%n_part,1,-1
759 aero_state%apa%particle(i_part))
real(kind=dp) elemental function vol2diam(v)
Convert volume (m^3) to diameter (m).
The aero_data_t structure and associated subroutines.
subroutine die_msg(code, error_msg)
Error immediately.
real(kind=dp) function aero_particle_solute_kappa(aero_particle, aero_data)
Returns the average of the solute kappas (1).
real(kind=dp) elemental function diam2vol(d)
Convert diameter (m) to volume (m^3).
Water condensation onto aerosol particles.
The aero_particle_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().
elemental real(kind=dp) function aero_particle_diameter(aero_particle)
Total diameter of the particle (m).
subroutine condense_equilib_particles(env_state, aero_data, aero_state)
Call condense_equilib_particle() on each particle in the aerosol to ensure that every particle has it...
real(kind=dp) function aero_weight_array_num_conc(aero_weight_array, aero_particle)
Compute the number concentration for a particle (m^{-3}).
real(kind=dp) function aero_particle_water_density(aero_data)
Returns the water density (kg/m^3).
subroutine warn_assert_msg(code, condition_ok, warning_msg)
Prints a warning message if condition_ok is false.
subroutine aero_data_deallocate(aero_data)
Frees all storage.
subroutine condense_equilib_particle(env_state, aero_data, aero_particle)
Determine the water equilibrium state of a single particle.
Common utility subroutines.
Current environment state.
subroutine env_state_deallocate(env_state)
Free all storage.
int condense_solver(int neq, double *x_f, double *abstol_f, double reltol_f, double t_initial_f, double t_final_f)
Call the ODE solver.
The aero_state_t structure and assocated subroutines.
subroutine condense_rates(inputs, outputs)
Compute the rate of change of particle diameter and relative humidity for a single particle...
character(len=pmc_util_convert_string_len) function integer_to_string(val)
Convert an integer to a string format.
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 env_state_copy(env_from, env_to)
env_to = env_from
subroutine condense_particles(aero_state, aero_data, env_state_initial, env_state_final, del_t)
Do condensation to all the particles for a given time interval, including updating the environment to...
subroutine aero_data_copy(aero_data_from, aero_data_to)
Copy structure.
real(kind=dp) function aero_particle_solute_volume(aero_particle, aero_data)
Returns the total solute volume (m^3).
Internal-use structure for storing the outputs from the rate-calculation function.
subroutine aero_data_allocate(aero_data)
Allocate storage for aero_data.
subroutine env_state_allocate(env_state)
Allocate an empty environment.
real(kind=dp) function env_state_sat_vapor_pressure(env_state)
Computes the current saturation vapor pressure (Pa).
character(len=pmc_util_convert_string_len) function real_to_string(val)
Convert a real to a string format.
Aerosol material properties and associated data.