deposition fraction for a breathing model

This commit is contained in:
Luis Aleixo 2021-09-14 14:25:23 +02:00
parent f9c92329d8
commit 34b1f89d69
2 changed files with 60 additions and 40 deletions

View file

@ -39,7 +39,7 @@ print('\n')
#present_vl_er_histograms(activity='Heavy exercise', mask='No mask')
############ CDFs for comparing the QR-Values in different scenarios ############
generate_cdf_curves()
#generate_cdf_curves()
############ Deposition Fraction Graph ############
print('\n<<<<<<<<<<< Deposition Fraction for Breathing, seated >>>>>>>>>>>')

View file

@ -491,64 +491,84 @@ def calculate_cunningham_slip_factor(d:int):
return 1 + (((2*_lambda)/d) * (A1+(A2*(math.e**((-A3*d)/_lambda)))))
# Deposition factor for a breathing model
# Deposition factors for a breathing model
def calculate_deposition_factor():
exposure_mc = breathing_exposure(activity='Heavy exercise', mask='No mask')
exposure_model = exposure_mc.build_model(size=SAMPLE_SIZE)
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
############ Breathing models ############
br_seated = breathing_seated_exposure()
br_seated_model = br_seated.build_model(size=SAMPLE_SIZE)
br_light_activity = breathing_light_activity_exposure()
br_light_activity_model = br_light_activity.build_model(size=SAMPLE_SIZE)
br_heavy_exercise = breathing_heavy_exercise_exposure()
br_heavy_exercise_model = br_heavy_exercise.build_model(size=SAMPLE_SIZE)
rho_p = 1000
mu_air = 1.8*10**-5
FRC = 0.003
Vt = 0.0004
g = 9.8
BRk = exposure_model.exposed.activity.inhalation_rate
k = 1.38*10**-23
T = 300
diameters = np.linspace(0.01, 100, 200) #particle diameter (multiply later by 10**(-6))
fractions = []
for d in diameters:
d_μm = d*10**(-6)
cunningham_slip_factor = calculate_cunningham_slip_factor(d_μm)
diameters = np.linspace(0.001, 100, 200) #particle diameter (multiply later by 10**(-6))
activity_fractions = []
for scenario in (br_seated_model, br_light_activity_model, br_heavy_exercise_model):
BRk = scenario.exposed.activity.inhalation_rate
fractions = []
for d_μm in diameters:
d = d_μm*10**(-6)
cunningham_slip_factor = calculate_cunningham_slip_factor(d)
f_dep_ine = 0.08 + 0.92 / (
1 + (4.09*10**-6 * (
(((cunningham_slip_factor*rho_p*d_μm**2*(BRk/3600))/mu_air*FRC)**0.8) + (
0.01*(
((cunningham_slip_factor*g*rho_p*d_μm**2*FRC**(2/3))/(mu_air*(BRk/3600))**0.4) * (
(Vt/FRC)**0.8
)
f_dep_ine = 0.08 + 0.92 / (
1 + (4.09*10**-6 * (
(((cunningham_slip_factor*rho_p*d**2*(BRk/3600))/mu_air*FRC)**0.8) + (
0.01*(
((cunningham_slip_factor*g*rho_p*d**2*FRC**(2/3))/(mu_air*(BRk/3600))**0.4) * (
(Vt/FRC)**0.8
)
)
)
)**(-2.06)
))
f_dep_diff = 1 - (1/
((7380*(
((k * T * cunningham_slip_factor)/3*math.pi*mu_air*d_μm) *
((Vt**(1/3)) / (BRk/3600)))**0.539) * ((Vt/FRC)**0.884) + 1))
f_dep_diff = 1 - 1 / (
7380*(((k * T * cunningham_slip_factor)/(3 * math.pi * mu_air * d)*(Vt**(1/3))/(BRk/3600))**0.539 * ((Vt/FRC)**0.884)) + 1)
f_dep_sed = (0.431 * f_dep_ine) + (0.541 * f_dep_diff) + (1.060 * f_dep_ine**2) + (0.685 * f_dep_diff**2) - (1.521 * f_dep_ine * f_dep_diff)
# Eq. S.1
if d < 0.3:
f_dep = f_dep_diff
elif d >= 0.3 and d <= 1.:
f_dep = f_dep_sed
elif d > 1:
f_dep = f_dep_ine
fractions.append(f_dep)
fig = plt.figure()
ax = fig.add_subplot(1, 1, 1)
ax.plot(diameters, fractions)
f_dep_sed = (0.431 * f_dep_ine) + (0.541 * f_dep_diff) + (1.060 * f_dep_ine**2) + (0.685 * f_dep_diff**2) - (1.521 * f_dep_ine * f_dep_diff)
# Eq. S.1
if d_μm < 0.1:
f_dep = f_dep_diff
elif d_μm > 1:
f_dep = f_dep_ine
else:
f_dep = f_dep_sed
fractions.append(f_dep)
activity_fractions.append(fractions)
# Colors can be changed here
colors = ['skyblue', 'royalblue', 'navy']
breathing_rates = ['Seated', 'Light activity', 'Heavy activity']
lines = [mlines.Line2D([], [], color=color, markersize=15, label=label)
for color, label in zip(colors, breathing_rates)]
for i in range(3):
ax.plot(diameters, activity_fractions[i], c=colors[i])
ax.grid(linestyle='--')
ax.legend(handles=lines, loc='upper left')
ax.set_xscale('log')
plt.ylabel(
'Deposition fraction (f$_{dep}$)', fontsize=14)
plt.xlabel('Particle diameter (μm)', fontsize=14)
plt.ylabel('Deposition fraction (f$_{dep}$)', fontsize=14)
plt.show()
######### Auxiliar functions #########