diff --git a/cara/apps/calculator/model_generator.py b/cara/apps/calculator/model_generator.py index 196487bc..c7ebd410 100644 --- a/cara/apps/calculator/model_generator.py +++ b/cara/apps/calculator/model_generator.py @@ -45,6 +45,7 @@ class FormData: floor_area: float hepa_amount: float hepa_option: bool + humidity: str infected_coffee_break_option: str #Used if infected_dont_have_breaks_with_exposed infected_coffee_duration: int #Used if infected_dont_have_breaks_with_exposed infected_dont_have_breaks_with_exposed: bool @@ -54,6 +55,7 @@ class FormData: infected_lunch_start: minutes_since_midnight #Used if infected_dont_have_breaks_with_exposed infected_people: int infected_start: minutes_since_midnight + inside_temp: float location_name: str location_latitude: float location_longitude: float @@ -100,6 +102,7 @@ class FormData: 'floor_area': 0., 'hepa_amount': 0., 'hepa_option': False, + 'humidity': '', 'infected_coffee_break_option': 'coffee_break_0', 'infected_coffee_duration': 5, 'infected_dont_have_breaks_with_exposed': False, @@ -109,6 +112,7 @@ class FormData: 'infected_lunch_start': '12:30', 'infected_people': _NO_DEFAULT, 'infected_start': '08:30', + 'inside_temp': 293., 'location_latitude': _NO_DEFAULT, 'location_longitude': _NO_DEFAULT, 'location_name': _NO_DEFAULT, @@ -240,11 +244,14 @@ class FormData: volume = self.room_volume else: volume = self.floor_area * self.ceiling_height - if self.room_heating_option: - humidity = 0.3 + if self.humidity == '': + if self.room_heating_option: + humidity = 0.3 + else: + humidity = 0.5 else: - humidity = 0.5 - room = models.Room(volume=volume, humidity=humidity) + humidity = float(self.humidity) + room = models.Room(volume=volume, inside_temp=models.PiecewiseConstant((0, 24), (self.inside_temp,)), humidity=humidity) infected_population = self.infected_population() @@ -329,13 +336,11 @@ class FormData: window_interval = always_on outside_temp = self.outside_temp() - inside_temp = models.PiecewiseConstant((0, 24), (293,)) ventilation: models.Ventilation if self.window_type == 'window_sliding': ventilation = models.SlidingWindow( active=window_interval, - inside_temp=inside_temp, outside_temp=outside_temp, window_height=self.window_height, opening_length=self.opening_distance, @@ -344,7 +349,6 @@ class FormData: elif self.window_type == 'window_hinged': ventilation = models.HingedWindow( active=window_interval, - inside_temp=inside_temp, outside_temp=outside_temp, window_height=self.window_height, window_width=self.window_width, @@ -689,6 +693,7 @@ def baseline_raw_form_data() -> typing.Dict[str, typing.Union[str, float]]: 'floor_area': '', 'hepa_amount': '250', 'hepa_option': '0', + 'humidity': '', 'infected_coffee_break_option': 'coffee_break_4', 'infected_coffee_duration': '10', 'infected_dont_have_breaks_with_exposed': '1', @@ -698,6 +703,7 @@ def baseline_raw_form_data() -> typing.Dict[str, typing.Union[str, float]]: 'infected_lunch_start': '12:30', 'infected_people': '1', 'infected_start': '09:00', + 'inside_temp': 293., 'location_latitude': 46.20833, 'location_longitude': 6.14275, 'location_name': 'Geneva', diff --git a/cara/apps/calculator/static/js/form.js b/cara/apps/calculator/static/js/form.js index fd4f3c73..e4597acc 100644 --- a/cara/apps/calculator/static/js/form.js +++ b/cara/apps/calculator/static/js/form.js @@ -783,6 +783,9 @@ $(document).ready(function () { templateSelection: formatLocationSelection }); + // Logic for the API requests. Always set humity input as the empty string so that we can profit from the "room_heating_option default" values for humidity. + $("[name='humidity']").val(""); + function formatlocation(suggestedLocation) { // Function is called for each location from the geocoding API. diff --git a/cara/apps/expert.py b/cara/apps/expert.py index 84beed27..84690d3c 100644 --- a/cara/apps/expert.py +++ b/cara/apps/expert.py @@ -243,19 +243,29 @@ class ModelWidgets(View): def _build_room(self, node): room_volume = widgets.IntSlider(value=node.volume, min=10, max=150) + inside_temp = widgets.IntSlider(value=node.inside_temp.values[0]-273.15, min=15., max=25.) - def on_value_change(change): + def on_volume_change(change): node.volume = change['new'] + def on_insidetemp_change(change): + node.inside_temp.values = (change['new']+273.15,) + # TODO: Link the state back to the widget, not just the other way around. - room_volume.observe(on_value_change, names=['value']) + room_volume.observe(on_volume_change, names=['value']) + inside_temp.observe(on_insidetemp_change, names=['value']) + def on_state_change(): room_volume.value = node.volume + inside_temp.value = node.inside_temp + node.dcs_observe(on_state_change) widget = collapsible( [widget_group( - [[widgets.Label('Room volume (m³)'), room_volume]] + [[widgets.Label('Room volume (m³)'), room_volume], + [widgets.Label('Inside temperature (℃)'), inside_temp], + ] )], title='Specification of workplace', ) @@ -281,7 +291,6 @@ class ModelWidgets(View): def _build_window(self, node) -> WidgetGroup: period = widgets.IntSlider(value=node.active.period, min=0, max=240) interval = widgets.IntSlider(value=node.active.duration, min=0, max=240) - inside_temp = widgets.IntSlider(value=node.inside_temp.values[0]-273.15, min=15., max=25.) def on_period_change(change): node.active.period = change['new'] @@ -289,13 +298,9 @@ class ModelWidgets(View): def on_interval_change(change): node.active.duration = change['new'] - def insidetemp_change(change): - node.inside_temp.values = (change['new']+273.15,) - # TODO: Link the state back to the widget, not just the other way around. period.observe(on_period_change, names=['value']) interval.observe(on_interval_change, names=['value']) - inside_temp.observe(insidetemp_change, names=['value']) outsidetemp_widgets = { 'Fixed': self._build_outsidetemp(node.outside_temp), @@ -327,10 +332,6 @@ class ModelWidgets(View): widgets.Label('Duration of opening (minutes)', layout=auto_width), interval, ), - ( - widgets.Label('Inside temperature (℃)', layout=auto_width), - inside_temp, - ), ( widgets.Label('Outside temperature scheme', layout=auto_width), outsidetemp_w, @@ -485,10 +486,9 @@ class ModelWidgets(View): baseline_model = models.ExposureModel( concentration_model=models.ConcentrationModel( - room=models.Room(volume=75), + room=models.Room(volume=75, inside_temp=models.PiecewiseConstant((0., 24.), (293.15,))), ventilation=models.SlidingWindow( active=models.PeriodicInterval(period=120, duration=15), - inside_temp=models.PiecewiseConstant((0., 24.), (293.15,)), outside_temp=models.PiecewiseConstant((0., 24.), (283.15,)), window_height=1.6, opening_length=0.6, ), diff --git a/cara/apps/templates/base/calculator.form.html.j2 b/cara/apps/templates/base/calculator.form.html.j2 index 499fc761..e1b66074 100644 --- a/cara/apps/templates/base/calculator.form.html.j2 +++ b/cara/apps/templates/base/calculator.form.html.j2 @@ -137,6 +137,8 @@ + +
diff --git a/cara/models.py b/cara/models.py index c6df3464..34c0d666 100644 --- a/cara/models.py +++ b/cara/models.py @@ -57,15 +57,6 @@ _VectorisedFloat = typing.Union[float, np.ndarray] _VectorisedInt = typing.Union[int, np.ndarray] -@dataclass(frozen=True) -class Room: - #: The total volume of the room - volume: _VectorisedFloat - - #: The humidity in the room (from 0 to 1 - e.g. 0.5 is 50% humidity) - humidity: _VectorisedFloat = 0.5 - - Time_t = typing.TypeVar('Time_t', float, int) BoundaryPair_t = typing.Tuple[Time_t, Time_t] BoundarySequence_t = typing.Union[typing.Tuple[BoundaryPair_t, ...], typing.Tuple] @@ -195,6 +186,18 @@ class PiecewiseConstant: ) +@dataclass(frozen=True) +class Room: + #: The total volume of the room + volume: _VectorisedFloat + + #: The temperature inside the room (Kelvin). + inside_temp: PiecewiseConstant = PiecewiseConstant((0, 24), (293,)) + + #: The humidity in the room (from 0 to 1 - e.g. 0.5 is 50% humidity) + humidity: _VectorisedFloat = 0.5 + + @dataclass(frozen=True) class _VentilationBase: """ @@ -207,7 +210,7 @@ class _VentilationBase: mechanical air exchange through a filter. """ - def transition_times(self) -> typing.Set[float]: + def transition_times(self, room: Room) -> typing.Set[float]: raise NotImplementedError("Subclass must implement") def air_exchange(self, room: Room, time: float) -> _VectorisedFloat: @@ -228,7 +231,7 @@ class Ventilation(_VentilationBase): #: The interval in which the ventilation is active. active: Interval - def transition_times(self) -> typing.Set[float]: + def transition_times(self, room: Room) -> typing.Set[float]: return self.active.transition_times() @@ -243,10 +246,10 @@ class MultipleVentilation(_VentilationBase): """ ventilations: typing.Tuple[_VentilationBase, ...] - def transition_times(self) -> typing.Set[float]: + def transition_times(self, room: Room) -> typing.Set[float]: transitions = set() for ventilation in self.ventilations: - transitions.update(ventilation.transition_times()) + transitions.update(ventilation.transition_times(room)) return transitions def air_exchange(self, room: Room, time: float) -> _VectorisedFloat: @@ -265,9 +268,6 @@ class WindowOpening(Ventilation): #: The interval in which the window is open. active: Interval - #: The temperature inside the room (Kelvin). - inside_temp: PiecewiseConstant - #: The temperature outside of the window (Kelvin). outside_temp: PiecewiseConstant @@ -292,9 +292,9 @@ class WindowOpening(Ventilation): """ raise NotImplementedError("Unknown discharge coefficient") - def transition_times(self) -> typing.Set[float]: - transitions = super().transition_times() - transitions.update(self.inside_temp.transition_times) + def transition_times(self, room: Room) -> typing.Set[float]: + transitions = super().transition_times(room) + transitions.update(room.inside_temp.transition_times) transitions.update(self.outside_temp.transition_times) return transitions @@ -304,7 +304,7 @@ class WindowOpening(Ventilation): return 0. # Reminder, no dependence on time in the resulting calculation. - inside_temp: _VectorisedFloat = self.inside_temp.value(time) + inside_temp: _VectorisedFloat = room.inside_temp.value(time) outside_temp: _VectorisedFloat = self.outside_temp.value(time) # The inside_temperature is forced to be always at least min_deltaT degree @@ -439,28 +439,35 @@ class Virus: #: Pre-populated examples of Viruses. types: typing.ClassVar[typing.Dict[str, "Virus"]] - def halflife(self, humidity: _VectorisedFloat) -> _VectorisedFloat: + def halflife(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: # Biological decay (inactivation of the virus in air) - virus # dependent and function of humidity raise NotImplementedError - def decay_constant(self, humidity: _VectorisedFloat) -> _VectorisedFloat: - # Viral inactivation per hour (h^-1) (function of humidity) - return np.log(2) / self.halflife(humidity) + def decay_constant(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: + # Viral inactivation per hour (h^-1) (function of humidity and inside temperature) + return np.log(2) / self.halflife(humidity, inside_temp) @dataclass(frozen=True) class SARSCoV2(Virus): - def halflife(self, humidity: _VectorisedFloat) -> _VectorisedFloat: + def halflife(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: """ Half-life changes with humidity level. Here is implemented a simple piecewise constant model (for more details see A. Henriques et al, CERN-OPEN-2021-004, DOI: 10.17181/CERN.1GDQ.5Y75) """ - # Taken from Morris et al (https://doi.org/10.7554/eLife.65902) data at T = 22°C and RH = 40 %, - # and from Doremalen et al (https://www.nejm.org/doi/10.1056/NEJMc2004973). - return np.piecewise(humidity, [humidity <= 0.4, humidity > 0.4], [6.43, 1.1]) + # Updated to use the formula from Dabish et al. with correction https://doi.org/10.1080/02786826.2020.1829536 + # with a maximum at hl = 6.43 (compensate for the negative decay values in the paper). + # Note that humidity is in percentage and inside_temp in °C. + # factor np.log(2) -> decay rate to half-life; factor 60 -> minutes to hours + hl_calc = ((np.log(2)/((0.16030 + 0.04018*(((inside_temp-273.15)-20.615)/10.585) + +0.02176*(((humidity*100)-45.235)/28.665) + -0.14369 + -0.02636*((inside_temp-273.15)-20.615)/10.585)))/60) + + return np.where(hl_calc <= 0, 6.43, np.minimum(6.43, hl_calc)) Virus.types = { @@ -917,9 +924,9 @@ class ConcentrationModel: h = 1.5 # Deposition rate (h^-1) k = (vg * 3600) / h - + #todo: Inside_temp needs to be exposed/added to the room; return ( - k + self.virus.decay_constant(self.room.humidity) + k + self.virus.decay_constant(self.room.humidity, self.room.inside_temp.value(time)) + self.ventilation.air_exchange(self.room, time) ) @@ -950,7 +957,7 @@ class ConcentrationModel: """ state_change_times = {0.} state_change_times.update(self.infected.presence.transition_times()) - state_change_times.update(self.ventilation.transition_times()) + state_change_times.update(self.ventilation.transition_times(self.room)) return sorted(state_change_times) @method_cache diff --git a/cara/tests/apps/calculator/test_model_generator.py b/cara/tests/apps/calculator/test_model_generator.py index 3bfb08f5..6298fd5c 100644 --- a/cara/tests/apps/calculator/test_model_generator.py +++ b/cara/tests/apps/calculator/test_model_generator.py @@ -60,7 +60,6 @@ def test_ventilation_slidingwindow(baseline_form: model_generator.FormData): window = models.SlidingWindow( active=models.PeriodicInterval(period=120, duration=10, start=9), - inside_temp=models.PiecewiseConstant((0, 24), (293,)), outside_temp=baseline_window.outside_temp, window_height=1.6, opening_length=0.6, ) @@ -92,7 +91,6 @@ def test_ventilation_hingedwindow(baseline_form: model_generator.FormData): window = models.HingedWindow( active=models.PeriodicInterval(period=120, duration=10, start=9), - inside_temp=models.PiecewiseConstant((0, 24), (293,)), outside_temp=baseline_window.outside_temp, window_height=1.6, window_width=1., opening_length=0.6, ) @@ -106,7 +104,7 @@ def test_ventilation_hingedwindow(baseline_form: model_generator.FormData): def test_ventilation_mechanical(baseline_form: model_generator.FormData): - room = models.Room(75) + room = models.Room(volume=75, inside_temp=models.PiecewiseConstant((0, 24), (293,))) mech = models.HVACMechanical( active=models.PeriodicInterval(period=120, duration=120), q_air_mech=500., @@ -121,7 +119,7 @@ def test_ventilation_mechanical(baseline_form: model_generator.FormData): def test_ventilation_airchanges(baseline_form: model_generator.FormData): - room = models.Room(75) + room = models.Room(75, inside_temp=models.PiecewiseConstant((0, 24), (293,))) airchange = models.AirChange( active=models.PeriodicInterval(period=120, duration=120), air_exch=3., @@ -153,7 +151,6 @@ def test_ventilation_window_hepa(baseline_form: model_generator.FormData): # Now build the equivalent ventilation instance directly, and compare. window = models.SlidingWindow( active=models.PeriodicInterval(period=120, duration=10, start=9), - inside_temp=models.PiecewiseConstant((0, 24), (293,)), outside_temp=baseline_window.outside_temp, window_height=1.6, opening_length=0.6, ) diff --git a/cara/tests/conftest.py b/cara/tests/conftest.py index 3cce3eb3..ef4bff1d 100644 --- a/cara/tests/conftest.py +++ b/cara/tests/conftest.py @@ -8,7 +8,7 @@ import pytest @pytest.fixture def baseline_concentration_model(): model = models.ConcentrationModel( - room=models.Room(volume=75), + room=models.Room(volume=75, inside_temp=models.PiecewiseConstant((0., 24.), (293,))), ventilation=models.AirChange( active=models.SpecificInterval(((0., 24.), )), air_exch=30., @@ -55,7 +55,6 @@ def exposure_model_w_outside_temp_changes(baseline_exposure_model: models.Exposu baseline_exposure_model, { 'concentration_model.ventilation': models.SlidingWindow( active=models.PeriodicInterval(2.2 * 60, 1.8 * 60), - inside_temp=models.PiecewiseConstant((0., 24.), (293,)), outside_temp=cara.data.GenevaTemperatures['Jan'], window_height=1.6, opening_length=0.6, diff --git a/cara/tests/models/test_concentration_model.py b/cara/tests/models/test_concentration_model.py index 432b4d80..ea11a005 100644 --- a/cara/tests/models/test_concentration_model.py +++ b/cara/tests/models/test_concentration_model.py @@ -26,7 +26,7 @@ def test_concentration_model_vectorisation(override_params): always = models.PeriodicInterval(240, 240) # TODO: This should be a thing on an interval. c_model = models.ConcentrationModel( - models.Room(defaults['volume'], defaults['humidity']), + models.Room(defaults['volume'], models.PiecewiseConstant((0., 24.), (293,)), defaults['humidity']), models.AirChange(always, defaults['air_change']), models.InfectedPopulation( number=1, @@ -59,7 +59,7 @@ def test_concentration_model_vectorisation(override_params): def simple_conc_model(): interesting_times = models.SpecificInterval(([0.5, 1.], [1.1, 2], [2., 3.]), ) return models.ConcentrationModel( - models.Room(75), + models.Room(75, models.PiecewiseConstant((0., 24.), (293,))), models.AirChange(interesting_times, 100), models.InfectedPopulation( number=1, diff --git a/cara/tests/models/test_exposure_model.py b/cara/tests/models/test_exposure_model.py index 44979024..896c1c6b 100644 --- a/cara/tests/models/test_exposure_model.py +++ b/cara/tests/models/test_exposure_model.py @@ -152,7 +152,7 @@ def conc_model(): ) always = models.SpecificInterval(((0., 24.), )) return models.ConcentrationModel( - models.Room(25), + models.Room(25, models.PiecewiseConstant((0., 24.), (293,))), models.AirChange(always, 5), models.EmittingPopulation( number=1, @@ -179,12 +179,12 @@ def sr_model(): @pytest.mark.parametrize( ["exposed_time_interval", "expected_deposited_exposure"], [ - [(0., 1.), 45.6008710], - [(1., 1.01), 0.5280401], - [(1.01, 1.02), 0.51314096385], - [(12., 12.01), 0.016255813185], - [(12., 24.), 645.63619275], - [(0., 24.), 700.7322474], + [(0., 1.), 42.63222033436878], + [(1., 1.01), 0.485377549596179], + [(1.01, 1.02), 0.47058239520823814], + [(12., 12.01), 0.01622776617499709], + [(12., 24.), 595.1115223695439], + [(0., 24.), 645.8401125684933], ] ) def test_exposure_model_integral_accuracy(exposed_time_interval, diff --git a/cara/tests/models/test_virus.py b/cara/tests/models/test_virus.py new file mode 100644 index 00000000..83f4d1f0 --- /dev/null +++ b/cara/tests/models/test_virus.py @@ -0,0 +1,24 @@ +import numpy as np +import numpy.testing as npt +import pytest + +from cara import models + +@pytest.mark.parametrize( + "inside_temp, humidity, expected_halflife, expected_decay_constant", + [ + [293.15, 0.5, 0.5947447349860315, 1.1654532436949188], + [272.15, 0.7, 1.6070844193207476, 0.4313072619127947], + [300.15, 1., 0.17367078830147223, 3.9911558376571805], + [300.15, 0., 6.43, 0.10779893943389507], + [np.array([272.15, 300.15]), np.array([0.7, 0.]), + np.array([1.60708442, 6.43]), np.array([0.43130726, 0.10779894])], + [np.array([293.15, 300.15]), np.array([0.5, 1.]), + np.array([0.59474473, 0.17367079]), np.array([1.16545324, 3.99115584])] + ], +) +def test_decay_constant(inside_temp, humidity, expected_halflife, expected_decay_constant): + npt.assert_almost_equal(models.Virus.types['SARS_CoV_2'].halflife(humidity, inside_temp), + expected_halflife) + npt.assert_almost_equal(models.Virus.types['SARS_CoV_2'].decay_constant(humidity, inside_temp), + expected_decay_constant) \ No newline at end of file diff --git a/cara/tests/test_full_algorithm.py b/cara/tests/test_full_algorithm.py index 372ca8b1..53a4e99e 100644 --- a/cara/tests/test_full_algorithm.py +++ b/cara/tests/test_full_algorithm.py @@ -84,8 +84,13 @@ class SimpleConcentrationModel: """ removal rate lambda in h^-1, excluding the deposition rate. """ - return (self.lambda_ventilation - + ln2/(6.43 if self.humidity<=0.4 else 1.1) ) + hl_calc = ((ln2/((0.16030 + 0.04018*(((293-273.15)-20.615)/10.585) + +0.02176*(((self.humidity*100)-45.235)/28.665) + -0.14369 + -0.02636*((293-273.15)-20.615)/10.585)))/60) + + return (self.lambda_ventilation + + ln2/(np.where(hl_calc <= 0, 6.43, np.minimum(6.43, hl_calc)))) @method_cache def deposition_removal_coefficient(self) -> float: @@ -461,7 +466,7 @@ interaction_intervals = (models.SpecificInterval(present_times=((10.5, 11.0),)), @pytest.fixture def c_model() -> mc.ConcentrationModel: return mc.ConcentrationModel( - room=models.Room(volume=50, humidity=0.3), + room=models.Room(volume=50, inside_temp=models.PiecewiseConstant((0., 24.), (293,)), humidity=0.3), ventilation=models.AirChange(active=models.PeriodicInterval(period=120, duration=120), air_exch=1.), infected=mc.InfectedPopulation( number=1, diff --git a/cara/tests/test_known_quantities.py b/cara/tests/test_known_quantities.py index 78c4541f..055c587b 100644 --- a/cara/tests/test_known_quantities.py +++ b/cara/tests/test_known_quantities.py @@ -19,7 +19,6 @@ def test_no_mask_superspeading_emission_rate(baseline_concentration_model): def baseline_periodic_window(): return models.SlidingWindow( active=models.PeriodicInterval(period=120, duration=15), - inside_temp=models.PiecewiseConstant((0., 24.), (293,)), outside_temp=models.PiecewiseConstant((0., 24.), (283,)), window_height=1.6, opening_length=0.6, ) @@ -27,7 +26,7 @@ def baseline_periodic_window(): @pytest.fixture def baseline_room(): - return models.Room(volume=75) + return models.Room(volume=75, inside_temp=models.PiecewiseConstant((0., 24.), (293,))) @pytest.fixture @@ -44,7 +43,7 @@ def test_concentrations(baseline_concentration_model): concentrations = [baseline_concentration_model.concentration(float(t)) for t in ts] npt.assert_allclose( concentrations, - [0.000000e+00, 20.805628, 6.602814e-13, 20.805628, 2.09545e-26], + [0.000000e+00, 2.046096e+01, 3.846725e-13, 2.046096e+01, 7.231966e-27], rtol=1e-6 ) @@ -95,7 +94,7 @@ def test_r0(baseline_exposure_model): # expected r0 was computed with a trapezoidal integration, using # a mesh of 100'000 pts per exposed presence interval. r0 = baseline_exposure_model.reproduction_number() - npt.assert_allclose(r0, 776.941990) + npt.assert_allclose(r0, 771.380385) def test_periodic_window(baseline_periodic_window, baseline_room): @@ -131,11 +130,10 @@ def test_periodic_hepa(baseline_periodic_hepa, baseline_room): ], ) def test_multiple_ventilation_HEPA_window(baseline_periodic_hepa, time, expected_value): - room = models.Room(volume=68.) + room = models.Room(volume=68., inside_temp=models.PiecewiseConstant((0., 24.),(293.15,))) tempOutside = models.PiecewiseConstant((0., 1., 2.5),(273.15, 283.15)) - tempInside = models.PiecewiseConstant((0., 24.),(293.15,)) window = models.SlidingWindow(active=models.SpecificInterval([(1 / 60, 24.)]), - inside_temp=tempInside,outside_temp=tempOutside, + outside_temp=tempOutside, window_height=1.,opening_length=0.6) vent = models.MultipleVentilation([window, baseline_periodic_hepa]) npt.assert_allclose(vent.air_exchange(room,time), expected_value, rtol=1e-5) @@ -143,12 +141,12 @@ def test_multiple_ventilation_HEPA_window(baseline_periodic_hepa, time, expected def test_multiple_ventilation_HEPA_window_transitions(baseline_periodic_hepa): tempOutside = models.PiecewiseConstant((0., 1., 2.5),(273.15, 283.15)) - tempInside = models.PiecewiseConstant((0., 24.),(293.15,)) + room = models.Room(68, models.PiecewiseConstant((0., 24.),(293.15,))) window = models.SlidingWindow(active=models.SpecificInterval([(1 / 60, 24.)]), - inside_temp=tempInside,outside_temp=tempOutside, + outside_temp=tempOutside, window_height=1.,opening_length=0.6) vent = models.MultipleVentilation([window, baseline_periodic_hepa]) - assert set(vent.transition_times()) == set([0.0, 1/60, 0.25, 1.0, 2.0, 2.25, + assert set(vent.transition_times(room)) == set([0.0, 1/60, 0.25, 1.0, 2.0, 2.25, 2.5, 4.0, 4.25, 6.0, 6.25, 8.0, 8.25, 10.0, 10.25, 12.0, 12.25, 14.0, 14.25, 16.0, 16.25, 18.0, 18.25, 20.0, 20.25, 22.0, 22.25, 24.]) @@ -188,14 +186,13 @@ def test_multiple_ventilation_HEPA_HVAC_AirChange(volume, expected_value): ) def test_windowopening(time, expected_value): tempOutside = models.PiecewiseConstant((0., 10., 24.),(273.15, 283.15)) - tempInside = models.PiecewiseConstant((0., 24.), (293.15,)) w = models.SlidingWindow( active=models.SpecificInterval([(0., 24.)]), - inside_temp=tempInside,outside_temp=tempOutside, + outside_temp=tempOutside, window_height=1., opening_length=0.6, ) npt.assert_allclose( - w.air_exchange(models.Room(volume=68), time), expected_value, rtol=1e-5 + w.air_exchange(models.Room(volume=68, inside_temp=models.PiecewiseConstant((0., 24.), (293.15, ))), time), expected_value, rtol=1e-5 ) @@ -223,10 +220,9 @@ def build_hourly_dependent_model( outside_temp = temperatures[month] model = models.ConcentrationModel( - room=models.Room(volume=75), + room=models.Room(volume=75, inside_temp=models.PiecewiseConstant((0., 24.), (293, ))), ventilation=models.SlidingWindow( active=models.SpecificInterval(intervals_open), - inside_temp=models.PiecewiseConstant((0., 24.), (293, )), outside_temp=outside_temp, window_height=1.6, opening_length=0.6, ), @@ -246,10 +242,9 @@ def build_hourly_dependent_model( def build_constant_temp_model(outside_temp, intervals_open=((7.5, 8.5),)): model = models.ConcentrationModel( - room=models.Room(volume=75), + room=models.Room(volume=75, inside_temp=models.PiecewiseConstant((0., 24.), (293,))), ventilation=models.SlidingWindow( active=models.SpecificInterval(intervals_open), - inside_temp=models.PiecewiseConstant((0., 24.), (293,)), outside_temp=models.PiecewiseConstant((0., 24.), (outside_temp,)), window_height=1.6, opening_length=0.6, ), @@ -271,7 +266,6 @@ def build_hourly_dependent_model_multipleventilation(month, intervals_open=((7.5 vent = models.MultipleVentilation(( models.SlidingWindow( active=models.SpecificInterval(intervals_open), - inside_temp=models.PiecewiseConstant((0., 24.), (293,)), outside_temp=data.GenevaTemperatures[month], window_height=1.6, opening_length=0.6, ), @@ -281,7 +275,7 @@ def build_hourly_dependent_model_multipleventilation(month, intervals_open=((7.5 ), )) model = models.ConcentrationModel( - room=models.Room(volume=75), + room=models.Room(volume=75, inside_temp=models.PiecewiseConstant((0., 24.), (293,))), ventilation=vent, infected=models.EmittingPopulation( number=1, @@ -387,8 +381,8 @@ def build_exposure_model(concentration_model, short_range_model): @pytest.mark.parametrize( "month, expected_deposited_exposure", [ - ['Jan', 377.440565819], - ['Jun', 1721.03336729], + ['Jan', 359.140499], + ['Jun', 1385.917562], ], ) def test_exposure_hourly_dep(month,expected_deposited_exposure, baseline_sr_model): @@ -408,8 +402,8 @@ def test_exposure_hourly_dep(month,expected_deposited_exposure, baseline_sr_mode @pytest.mark.parametrize( "month, expected_deposited_exposure", [ - ['Jan', 383.339206111], - ['Jun', 1799.17597184], + ['Jan', 359.983716], + ['Jun', 1439.267381], ], ) def test_exposure_hourly_dep_refined(month,expected_deposited_exposure, baseline_sr_model): diff --git a/cara/tests/test_monte_carlo.py b/cara/tests/test_monte_carlo.py index a398e716..b948f4b9 100644 --- a/cara/tests/test_monte_carlo.py +++ b/cara/tests/test_monte_carlo.py @@ -40,10 +40,10 @@ def test_type_annotations(): @pytest.fixture def baseline_mc_concentration_model() -> cara.monte_carlo.ConcentrationModel: mc_model = cara.monte_carlo.ConcentrationModel( - room=cara.monte_carlo.Room(volume=cara.monte_carlo.sampleable.Normal(75, 20)), + room=cara.monte_carlo.Room(volume=cara.monte_carlo.sampleable.Normal(75, 20), + inside_temp=cara.models.PiecewiseConstant((0., 24.), (293,))), ventilation=cara.monte_carlo.SlidingWindow( active=cara.models.PeriodicInterval(period=120, duration=120), - inside_temp=cara.models.PiecewiseConstant((0., 24.), (293,)), outside_temp=cara.models.PiecewiseConstant((0., 24.), (283,)), window_height=1.6, opening_length=0.6, ), diff --git a/cara/tests/test_monte_carlo_full_models.py b/cara/tests/test_monte_carlo_full_models.py index caed9c52..3bce07ae 100644 --- a/cara/tests/test_monte_carlo_full_models.py +++ b/cara/tests/test_monte_carlo_full_models.py @@ -45,12 +45,11 @@ def shared_office_mc(): Corresponds to the 1st line of Table 4 in https://doi.org/10.1101/2021.10.14.21264988 """ concentration_mc = mc.ConcentrationModel( - room=models.Room(volume=50, humidity=0.5), + room=models.Room(volume=50, inside_temp=models.PiecewiseConstant((0., 24.), (298,)), humidity=0.5), ventilation=models.MultipleVentilation( ventilations=( models.SlidingWindow( active=models.PeriodicInterval(period=120, duration=120), - inside_temp=models.PiecewiseConstant((0., 24.), (298,)), outside_temp=data.GenevaTemperatures['Jun'], window_height=1.6, opening_length=0.2, @@ -88,12 +87,11 @@ def classroom_mc(): Corresponds to the 2nd line of Table 4 in https://doi.org/10.1101/2021.10.14.21264988 """ concentration_mc = mc.ConcentrationModel( - room=models.Room(volume=160, humidity=0.3), + room=models.Room(volume=160, inside_temp=models.PiecewiseConstant((0., 24.), (293,)), humidity=0.3), ventilation=models.MultipleVentilation( ventilations=( models.SlidingWindow( active=models.PeriodicInterval(period=120, duration=120), - inside_temp=models.PiecewiseConstant((0., 24.), (293,)), outside_temp=TorontoTemperatures['Dec'], window_height=1.6, opening_length=0.2, @@ -312,13 +310,13 @@ def waiting_room_mc(): @pytest.mark.parametrize( "mc_model, expected_pi, expected_new_cases, expected_dose, expected_ER", [ - ["shared_office_mc", 6.03, 0.18, 3.198, 809], - ["classroom_mc", 9.5, 1.85, 9.478, 5624], - ["ski_cabin_mc", 16.0, 0.5, 17.315, 7966], - ["skagit_chorale_mc",65.7, 40.0, 102.213, 190422], - ["bus_ride_mc", 12.0, 8.0, 7.65, 5419], - ["gym_mc", 0.45, 0.13, 0.208, 1145], - ["waiting_room_mc", 1.59, 0.22, 0.821, 737], + ["shared_office_mc", 5.55, 0.17, 2.699, 809], + ["classroom_mc", 9.58, 1.82, 9.034, 5624], + ["ski_cabin_mc", 16.0, 0.47, 17.315, 7966], + ["skagit_chorale_mc",61.01, 36.53, 84.730, 190422], + ["bus_ride_mc", 10.59, 7.06, 6.65, 5419], + ["gym_mc", 0.43, 0.12, 0.197, 1145], + ["waiting_room_mc", 1.34, 0.18, 0.670, 737], ] ) def test_report_models(mc_model, expected_pi, expected_new_cases, @@ -339,21 +337,20 @@ def test_report_models(mc_model, expected_pi, expected_new_cases, @pytest.mark.parametrize( "mask_type, month, expected_pi, expected_dose, expected_ER", [ - ["No mask", "Jul", 9.52, 9.920, 809], - ["Type I", "Jul", 1.7, 0.913, 149], - ["FFP2", "Jul", 0.51, 0.239, 149], - ["Type I", "Feb", 0.57, 0.272, 149], + ["No mask", "Jul", 8.46, 8.113, 809], + ["Type I", "Jul", 1.44, 0.727, 149], + ["FFP2", "Jul", 0.43, 0.197, 149], + ["Type I", "Feb", 0.54, 0.253, 149], ], ) def test_small_shared_office_Geneva(mask_type, month, expected_pi, expected_dose, expected_ER): concentration_mc = mc.ConcentrationModel( - room=models.Room(volume=33, humidity=0.5), + room=models.Room(volume=33, inside_temp=models.PiecewiseConstant((0., 24.), (293,)), humidity=0.5), ventilation=models.MultipleVentilation( ( models.SlidingWindow( active=models.SpecificInterval(((0., 24.),)), - inside_temp=models.PiecewiseConstant((0., 24.), (293,)), outside_temp=data.GenevaTemperatures[month], window_height=1.5, opening_length=0.2, ), diff --git a/cara/tests/test_ventilation.py b/cara/tests/test_ventilation.py index 134acf57..201da2c7 100644 --- a/cara/tests/test_ventilation.py +++ b/cara/tests/test_ventilation.py @@ -11,7 +11,6 @@ from cara import models def baseline_slidingwindow(): return models.SlidingWindow( active=models.SpecificInterval(((0, 4), (5, 9))), - inside_temp=models.PiecewiseConstant((0, 24), (293,)), outside_temp=models.PiecewiseConstant((0, 24), (283,)), window_height=1.6, opening_length=0.6, ) @@ -21,14 +20,13 @@ def baseline_slidingwindow(): def baseline_hingedwindow(): return models.HingedWindow( active=models.SpecificInterval(((0, 4), (5, 9))), - inside_temp=models.PiecewiseConstant((0, 24), (293,)), outside_temp=models.PiecewiseConstant((0, 24), (283,)), window_height=1.6, opening_length=0.6, window_width=1., ) def test_number_of_windows(baseline_slidingwindow): - room = models.Room(75) + room = models.Room(volume=75, inside_temp=models.PiecewiseConstant((0, 24), (293,))) two_windows = dataclasses.replace(baseline_slidingwindow, number_of_windows=2) one_window_exchange = baseline_slidingwindow.air_exchange(room, 1) @@ -63,9 +61,6 @@ def test_hinged_window(baseline_hingedwindow, window_width, {'outside_temp': models.PiecewiseConstant( (0, 2, 3), (np.array([20, 30, 28]), np.array([25, 30, 27])) )}, - {'inside_temp': models.PiecewiseConstant( - (0, 20), (np.array([20, 30, 25]), ) - )}, ] ) def test_hinged_window_vectorisation(override_params): @@ -73,11 +68,10 @@ def test_hinged_window_vectorisation(override_params): 'window_height': 0.15, 'window_width': 0.15, 'opening_length': 0.15, - 'inside_temp': models.PiecewiseConstant((0, 2, 3), (20, 25)), 'outside_temp': models.PiecewiseConstant((0, 2, 3), (10, 15)), } defaults.update(override_params) - room = models.Room(volume=75) + room = models.Room(volume=75, inside_temp=models.PiecewiseConstant((0, 2, 3), (20, 25))) t = 0.5 window = models.HingedWindow(models.PeriodicInterval(60, 30), **defaults) if {'window_height', 'opening_length', 'window_width'}.intersection(override_params):