From 5b715937d2c98090f7f9c7a46bc0131e3963dc93 Mon Sep 17 00:00:00 2001 From: jdevine Date: Thu, 21 Apr 2022 13:42:55 +0200 Subject: [PATCH 01/12] Implemented Humidity + Temperature half life function from Dabish et al. --- cara/models.py | 18 +++++++++--------- 1 file changed, 9 insertions(+), 9 deletions(-) diff --git a/cara/models.py b/cara/models.py index bfa50139..40e8a422 100644 --- a/cara/models.py +++ b/cara/models.py @@ -439,28 +439,28 @@ 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: + def decay_constant(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: # Viral inactivation per hour (h^-1) (function of humidity) - return np.log(2) / self.halflife(humidity) + 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. https://doi.org/10.1080/02786826.2020.1829536 + # with a minimum at hl = 1.1 + return max(1.1, (0.693/(0.16030 + 0.04018((inside_temp-20.615)/10.585)+0.2176((humidity-45.235)/28.665)+0.1))) Virus.types = { @@ -917,9 +917,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) + self.ventilation.air_exchange(self.room, time) ) From 52300df16b5e3c2a55faa55d6d2e055199db1a60 Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Fri, 22 Apr 2022 15:20:20 +0200 Subject: [PATCH 02/12] Changed tests - inside temp is now a property of Room --- cara/apps/calculator/model_generator.py | 5 +- cara/apps/expert.py | 28 +++++----- cara/models.py | 56 ++++++++++--------- .../apps/calculator/test_model_generator.py | 7 +-- cara/tests/conftest.py | 3 +- cara/tests/models/test_concentration_model.py | 2 +- cara/tests/test_known_quantities.py | 28 ++++------ cara/tests/test_monte_carlo.py | 4 +- cara/tests/test_monte_carlo_full_models.py | 9 +-- cara/tests/test_ventilation.py | 10 +--- 10 files changed, 66 insertions(+), 86 deletions(-) diff --git a/cara/apps/calculator/model_generator.py b/cara/apps/calculator/model_generator.py index c4593b49..ec770589 100644 --- a/cara/apps/calculator/model_generator.py +++ b/cara/apps/calculator/model_generator.py @@ -244,7 +244,7 @@ class FormData: humidity = 0.3 else: humidity = 0.5 - room = models.Room(volume=volume, humidity=humidity) + room = models.Room(volume=volume, inside_temp=models.PiecewiseConstant((0, 24), (293,)), humidity=humidity) infected_population = self.infected_population() @@ -329,13 +329,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 +342,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, 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/models.py b/cara/models.py index 40e8a422..d696067a 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,10 @@ 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) + print(room.inside_temp.transition_times) + transitions.update(room.inside_temp.transition_times) transitions.update(self.outside_temp.transition_times) return transitions @@ -304,7 +305,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,20 +440,20 @@ class Virus: #: Pre-populated examples of Viruses. types: typing.ClassVar[typing.Dict[str, "Virus"]] - def halflife(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: + def halflife(self, humidity: _VectorisedFloat, inside_temp: PiecewiseConstant, time: float) -> _VectorisedFloat: # Biological decay (inactivation of the virus in air) - virus # dependent and function of humidity raise NotImplementedError - def decay_constant(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: + def decay_constant(self, humidity: _VectorisedFloat, inside_temp: PiecewiseConstant, time: float) -> _VectorisedFloat: # Viral inactivation per hour (h^-1) (function of humidity) - return np.log(2) / self.halflife(humidity, inside_temp) + return np.log(2) / self.halflife(humidity, inside_temp, time) @dataclass(frozen=True) class SARSCoV2(Virus): - def halflife(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: + def halflife(self, humidity: _VectorisedFloat, inside_temp: PiecewiseConstant, time: float) -> _VectorisedFloat: """ Half-life changes with humidity level. Here is implemented a simple piecewise constant model (for more details see A. Henriques et al, @@ -460,7 +461,8 @@ class SARSCoV2(Virus): """ # Updated to use the formula from Dabish et al. https://doi.org/10.1080/02786826.2020.1829536 # with a minimum at hl = 1.1 - return max(1.1, (0.693/(0.16030 + 0.04018((inside_temp-20.615)/10.585)+0.2176((humidity-45.235)/28.665)+0.1))) + inside_temp = inside_temp.value(time) + return max(1.1, (0.693/(0.16030 + 0.04018*(((inside_temp-274.15)-20.615)/10.585)+0.02176*((humidity-45.235)/28.665)+0.1))) Virus.types = { @@ -919,7 +921,7 @@ class ConcentrationModel: k = (vg * 3600) / h #todo: Inside_temp needs to be exposed/added to the room; return ( - k + self.virus.decay_constant(self.room.humidity, self.room.inside_temp) + k + self.virus.decay_constant(self.room.humidity, self.room.inside_temp, time) + self.ventilation.air_exchange(self.room, time) ) @@ -950,7 +952,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 5cc9e409..6f118a8e 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=minutes_since_midnight(9 * 60)), - 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=minutes_since_midnight(9 * 60)), - 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=minutes_since_midnight(9 * 60)), - 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..3b1d9292 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, diff --git a/cara/tests/test_known_quantities.py b/cara/tests/test_known_quantities.py index 78c4541f..ed775312 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 @@ -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, 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 b97956b9..32fdd460 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, @@ -348,12 +346,11 @@ def test_report_models(mc_model, expected_pi, expected_new_cases, 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): From a85a05316f6b58f50660a47df823e677bbf0c13e Mon Sep 17 00:00:00 2001 From: jdevine Date: Mon, 25 Apr 2022 11:36:49 +0200 Subject: [PATCH 03/12] Adjusted tests, modified expected outputs for new temperature and humidity dependent halflife method --- cara/models.py | 2 +- cara/tests/models/test_exposure_model.py | 12 ++++++------ cara/tests/test_full_algorithm.py | 10 ++++++++-- cara/tests/test_known_quantities.py | 12 ++++++------ cara/tests/test_monte_carlo_full_models.py | 16 ++++++++-------- 5 files changed, 29 insertions(+), 23 deletions(-) diff --git a/cara/models.py b/cara/models.py index d696067a..1adb810e 100644 --- a/cara/models.py +++ b/cara/models.py @@ -462,7 +462,7 @@ class SARSCoV2(Virus): # Updated to use the formula from Dabish et al. https://doi.org/10.1080/02786826.2020.1829536 # with a minimum at hl = 1.1 inside_temp = inside_temp.value(time) - return max(1.1, (0.693/(0.16030 + 0.04018*(((inside_temp-274.15)-20.615)/10.585)+0.02176*((humidity-45.235)/28.665)+0.1))) + return max(1.1, (0.693/(0.16030 + 0.04018*(((inside_temp-273.15)-20.615)/10.585)+0.02176*((humidity-45.235)/28.665)+0.1))) Virus.types = { diff --git a/cara/tests/models/test_exposure_model.py b/cara/tests/models/test_exposure_model.py index 44979024..abc476d0 100644 --- a/cara/tests/models/test_exposure_model.py +++ b/cara/tests/models/test_exposure_model.py @@ -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.), 48.19316], + [(1., 1.01), 0.566368], + [(1.01, 1.02), 0.551401], + [(12., 12.01), 0.016278], + [(12., 24.), 691.21381], + [(0., 24.), 750.258043], ] ) def test_exposure_model_integral_accuracy(exposed_time_interval, diff --git a/cara/tests/test_full_algorithm.py b/cara/tests/test_full_algorithm.py index fb6ade81..71f10f05 100644 --- a/cara/tests/test_full_algorithm.py +++ b/cara/tests/test_full_algorithm.py @@ -47,6 +47,9 @@ class SimpleConcentrationModel: #: room volume (m^3) room_volume: _VectorisedFloat + #: The temperature inside the room (Kelvin). + #inside_temp: float = 293.15 + #: ventilation rate (air changes per hour) - including HEPA lambda_ventilation: _VectorisedFloat @@ -84,8 +87,11 @@ 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) ) + + return (self.lambda_ventilation + + ln2/(max(1.1, (0.693 / (0.16030 + 0.04018 * (((22) - 20.615) / 10.585) + 0.02176 * ( + (self.humidity - 45.235) / 28.665) + 0.1))))) + #6.43 if self.humidity<=0.4 else 1.1) ) @method_cache def deposition_removal_coefficient(self) -> float: diff --git a/cara/tests/test_known_quantities.py b/cara/tests/test_known_quantities.py index ed775312..4eb94d0b 100644 --- a/cara/tests/test_known_quantities.py +++ b/cara/tests/test_known_quantities.py @@ -43,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.108144e+01, 1.004741e-12, 2.108144e+01, 4.788592e-26], rtol=1e-6 ) @@ -94,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, 781.293793) def test_periodic_window(baseline_periodic_window, baseline_room): @@ -381,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', 392.994454], + ['Jun', 2127.82386], ], ) def test_exposure_hourly_dep(month,expected_deposited_exposure, baseline_sr_model): @@ -402,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', 394.000261], + ['Jun', 2239.777906], ], ) def test_exposure_hourly_dep_refined(month,expected_deposited_exposure, baseline_sr_model): diff --git a/cara/tests/test_monte_carlo_full_models.py b/cara/tests/test_monte_carlo_full_models.py index 32fdd460..58f0a57d 100644 --- a/cara/tests/test_monte_carlo_full_models.py +++ b/cara/tests/test_monte_carlo_full_models.py @@ -310,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], + ["shared_office_mc", 6.524158, 0.195725, 3.413983, 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], + ["skagit_chorale_mc",69.84901, 40.0, 121.265911, 190422], + ["bus_ride_mc", 13.311721, 8.918853, 8.6662, 5419], ["gym_mc", 0.45, 0.13, 0.208, 1145], - ["waiting_room_mc", 1.59, 0.22, 0.821, 737], + ["waiting_room_mc", 1.854373, 0.259612, 0.99173, 737], ] ) def test_report_models(mc_model, expected_pi, expected_new_cases, @@ -337,10 +337,10 @@ 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, 162], + ["No mask", "Jul", 10.966992, 12.357222, 809], + ["Type I", "Jul", 2.084648, 1.137819, 149], + ["FFP2", "Jul", 0.654566, .309637, 149], + ["Type I", "Feb", 0.612366, 0.272, 162], ], ) def test_small_shared_office_Geneva(mask_type, month, expected_pi, From 64ebb663e70d0ded43d9c80f747a210609aa8360 Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Mon, 25 Apr 2022 15:14:46 +0200 Subject: [PATCH 04/12] Handled mypy errors and increased number of decimals on test exposure model --- cara/models.py | 4 ++-- cara/tests/models/test_concentration_model.py | 2 +- cara/tests/models/test_exposure_model.py | 14 +++++++------- 3 files changed, 10 insertions(+), 10 deletions(-) diff --git a/cara/models.py b/cara/models.py index 1adb810e..65f820e9 100644 --- a/cara/models.py +++ b/cara/models.py @@ -461,8 +461,8 @@ class SARSCoV2(Virus): """ # Updated to use the formula from Dabish et al. https://doi.org/10.1080/02786826.2020.1829536 # with a minimum at hl = 1.1 - inside_temp = inside_temp.value(time) - return max(1.1, (0.693/(0.16030 + 0.04018*(((inside_temp-273.15)-20.615)/10.585)+0.02176*((humidity-45.235)/28.665)+0.1))) + temperature: _VectorisedFloat = inside_temp.value(time) + return np.maximum(1.1, (0.693/(0.16030 + 0.04018*(((temperature-273.15)-20.615)/10.585)+0.02176*((humidity-45.235)/28.665)+0.1))) Virus.types = { diff --git a/cara/tests/models/test_concentration_model.py b/cara/tests/models/test_concentration_model.py index 3b1d9292..ea11a005 100644 --- a/cara/tests/models/test_concentration_model.py +++ b/cara/tests/models/test_concentration_model.py @@ -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 abc476d0..428de12d 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.), 48.19316], - [(1., 1.01), 0.566368], - [(1.01, 1.02), 0.551401], - [(12., 12.01), 0.016278], - [(12., 24.), 691.21381], - [(0., 24.), 750.258043], + [(0., 1.), 48.193159880096644], + [(1., 1.01), 0.5663683904832492], + [(1.01, 1.02), 0.5514013220457682], + [(12., 12.01), 0.016277647772557004], + [(12., 24.), 691.2138099188268], + [(0., 24.), 750.2580429542696], ] ) def test_exposure_model_integral_accuracy(exposed_time_interval, From 71078671f99dee879bfbcf36ff033466feeb26f5 Mon Sep 17 00:00:00 2001 From: jdevine Date: Mon, 25 Apr 2022 16:20:29 +0200 Subject: [PATCH 05/12] Updated formula to reflect published correction --- cara/models.py | 7 +++++-- cara/tests/test_full_algorithm.py | 7 ++++--- 2 files changed, 9 insertions(+), 5 deletions(-) diff --git a/cara/models.py b/cara/models.py index 65f820e9..f2ac2500 100644 --- a/cara/models.py +++ b/cara/models.py @@ -459,10 +459,13 @@ class SARSCoV2(Virus): piecewise constant model (for more details see A. Henriques et al, CERN-OPEN-2021-004, DOI: 10.17181/CERN.1GDQ.5Y75) """ - # Updated to use the formula from Dabish et al. https://doi.org/10.1080/02786826.2020.1829536 + # Updated to use the formula from Dabish et al. with correction https://doi.org/10.1080/02786826.2020.1829536 # with a minimum at hl = 1.1 temperature: _VectorisedFloat = inside_temp.value(time) - return np.maximum(1.1, (0.693/(0.16030 + 0.04018*(((temperature-273.15)-20.615)/10.585)+0.02176*((humidity-45.235)/28.665)+0.1))) + return np.maximum(1.1, (0.693/((0.16030 + 0.04018*(((temperature-273.15)-20.615)/10.585) + +0.02176*((humidity-45.235)/28.665) + -0.14369 + -0.2636*((temperature-273.15)-20.615)/10.585)))) Virus.types = { diff --git a/cara/tests/test_full_algorithm.py b/cara/tests/test_full_algorithm.py index 71f10f05..8a07321c 100644 --- a/cara/tests/test_full_algorithm.py +++ b/cara/tests/test_full_algorithm.py @@ -89,9 +89,10 @@ class SimpleConcentrationModel: """ return (self.lambda_ventilation - + ln2/(max(1.1, (0.693 / (0.16030 + 0.04018 * (((22) - 20.615) / 10.585) + 0.02176 * ( - (self.humidity - 45.235) / 28.665) + 0.1))))) - #6.43 if self.humidity<=0.4 else 1.1) ) + + ln2/(max(1.1, (0.693 / ((0.16030 + 0.04018 * (((22) - 20.615) / 10.585) + + 0.02176 * ((self.humidity - 45.235) / 28.665) + - 0.14369 + - 0.2636((22-20.615)/10.585))))))) @method_cache def deposition_removal_coefficient(self) -> float: From c665f31df98ff4089eb0b0c03974cf2bc9d87b7d Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Mon, 25 Apr 2022 17:24:03 +0200 Subject: [PATCH 06/12] adjusted tests for the formula from Dabish et al --- cara/models.py | 2 +- cara/tests/models/test_exposure_model.py | 12 ++++++------ cara/tests/test_full_algorithm.py | 11 ++++------- cara/tests/test_known_quantities.py | 12 ++++++------ cara/tests/test_monte_carlo_full_models.py | 20 ++++++++++---------- 5 files changed, 27 insertions(+), 30 deletions(-) diff --git a/cara/models.py b/cara/models.py index f2ac2500..0a94a2dd 100644 --- a/cara/models.py +++ b/cara/models.py @@ -465,7 +465,7 @@ class SARSCoV2(Virus): return np.maximum(1.1, (0.693/((0.16030 + 0.04018*(((temperature-273.15)-20.615)/10.585) +0.02176*((humidity-45.235)/28.665) -0.14369 - -0.2636*((temperature-273.15)-20.615)/10.585)))) + -0.02636*((temperature-273.15)-20.615)/10.585)))) Virus.types = { diff --git a/cara/tests/models/test_exposure_model.py b/cara/tests/models/test_exposure_model.py index 428de12d..5a6a5449 100644 --- a/cara/tests/models/test_exposure_model.py +++ b/cara/tests/models/test_exposure_model.py @@ -179,12 +179,12 @@ def sr_model(): @pytest.mark.parametrize( ["exposed_time_interval", "expected_deposited_exposure"], [ - [(0., 1.), 48.193159880096644], - [(1., 1.01), 0.5663683904832492], - [(1.01, 1.02), 0.5514013220457682], - [(12., 12.01), 0.016277647772557004], - [(12., 24.), 691.2138099188268], - [(0., 24.), 750.2580429542696], + [(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], ] ) def test_exposure_model_integral_accuracy(exposed_time_interval, diff --git a/cara/tests/test_full_algorithm.py b/cara/tests/test_full_algorithm.py index 8a07321c..3f701e5e 100644 --- a/cara/tests/test_full_algorithm.py +++ b/cara/tests/test_full_algorithm.py @@ -47,9 +47,6 @@ class SimpleConcentrationModel: #: room volume (m^3) room_volume: _VectorisedFloat - #: The temperature inside the room (Kelvin). - #inside_temp: float = 293.15 - #: ventilation rate (air changes per hour) - including HEPA lambda_ventilation: _VectorisedFloat @@ -89,10 +86,10 @@ class SimpleConcentrationModel: """ return (self.lambda_ventilation - + ln2/(max(1.1, (0.693 / ((0.16030 + 0.04018 * (((22) - 20.615) / 10.585) - + 0.02176 * ((self.humidity - 45.235) / 28.665) + + ln2/(np.maximum(1.1, (0.693 / ((0.16030 + 0.04018 * (((21) - 20.615) / 10.585) + + 0.02176*((self.humidity - 45.235) / 28.665) - 0.14369 - - 0.2636((22-20.615)/10.585))))))) + - 0.02636*((21-20.615)/10.585))))))) @method_cache def deposition_removal_coefficient(self) -> float: @@ -460,7 +457,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 4eb94d0b..ed775312 100644 --- a/cara/tests/test_known_quantities.py +++ b/cara/tests/test_known_quantities.py @@ -43,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, 2.108144e+01, 1.004741e-12, 2.108144e+01, 4.788592e-26], + [0.000000e+00, 20.805628, 6.602814e-13, 20.805628, 2.09545e-26], rtol=1e-6 ) @@ -94,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, 781.293793) + npt.assert_allclose(r0, 776.941990) def test_periodic_window(baseline_periodic_window, baseline_room): @@ -381,8 +381,8 @@ def build_exposure_model(concentration_model, short_range_model): @pytest.mark.parametrize( "month, expected_deposited_exposure", [ - ['Jan', 392.994454], - ['Jun', 2127.82386], + ['Jan', 377.440565819], + ['Jun', 1721.03336729], ], ) def test_exposure_hourly_dep(month,expected_deposited_exposure, baseline_sr_model): @@ -402,8 +402,8 @@ def test_exposure_hourly_dep(month,expected_deposited_exposure, baseline_sr_mode @pytest.mark.parametrize( "month, expected_deposited_exposure", [ - ['Jan', 394.000261], - ['Jun', 2239.777906], + ['Jan', 383.339206111], + ['Jun', 1799.17597184], ], ) def test_exposure_hourly_dep_refined(month,expected_deposited_exposure, baseline_sr_model): diff --git a/cara/tests/test_monte_carlo_full_models.py b/cara/tests/test_monte_carlo_full_models.py index 58f0a57d..48682758 100644 --- a/cara/tests/test_monte_carlo_full_models.py +++ b/cara/tests/test_monte_carlo_full_models.py @@ -310,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.524158, 0.195725, 3.413983, 809], - ["classroom_mc", 9.5, 1.85, 9.478, 5624], - ["ski_cabin_mc", 16.0, 0.5, 17.315, 7966], - ["skagit_chorale_mc",69.84901, 40.0, 121.265911, 190422], - ["bus_ride_mc", 13.311721, 8.918853, 8.6662, 5419], + ["shared_office_mc", 6.03, 0.18, 3.198, 809], + ["classroom_mc", 8.63, 1.64, 8.082, 5624], + ["ski_cabin_mc", 16.0, 0.47, 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.854373, 0.259612, 0.99173, 737], + ["waiting_room_mc", 1.59, 0.22, 0.821, 737], ] ) def test_report_models(mc_model, expected_pi, expected_new_cases, @@ -337,10 +337,10 @@ 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", 10.966992, 12.357222, 809], - ["Type I", "Jul", 2.084648, 1.137819, 149], - ["FFP2", "Jul", 0.654566, .309637, 149], - ["Type I", "Feb", 0.612366, 0.272, 162], + ["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, 162], ], ) def test_small_shared_office_Geneva(mask_type, month, expected_pi, From 71f36e0aa41d5acab314ab5d2a2f03fb62e5845b Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Wed, 27 Apr 2022 14:09:24 +0200 Subject: [PATCH 07/12] Changed inside_temp to be passed as vectorised_float in the decay_constant method call --- cara/models.py | 16 +++++++--------- 1 file changed, 7 insertions(+), 9 deletions(-) diff --git a/cara/models.py b/cara/models.py index 0a94a2dd..cb389be0 100644 --- a/cara/models.py +++ b/cara/models.py @@ -294,7 +294,6 @@ class WindowOpening(Ventilation): def transition_times(self, room: Room) -> typing.Set[float]: transitions = super().transition_times(room) - print(room.inside_temp.transition_times) transitions.update(room.inside_temp.transition_times) transitions.update(self.outside_temp.transition_times) return transitions @@ -440,20 +439,20 @@ class Virus: #: Pre-populated examples of Viruses. types: typing.ClassVar[typing.Dict[str, "Virus"]] - def halflife(self, humidity: _VectorisedFloat, inside_temp: PiecewiseConstant, time: float) -> _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, inside_temp: PiecewiseConstant, time: float) -> _VectorisedFloat: + def decay_constant(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: # Viral inactivation per hour (h^-1) (function of humidity) - return np.log(2) / self.halflife(humidity, inside_temp, time) + return np.log(2) / self.halflife(humidity, inside_temp) @dataclass(frozen=True) class SARSCoV2(Virus): - def halflife(self, humidity: _VectorisedFloat, inside_temp: PiecewiseConstant, time: float) -> _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, @@ -461,11 +460,10 @@ class SARSCoV2(Virus): """ # Updated to use the formula from Dabish et al. with correction https://doi.org/10.1080/02786826.2020.1829536 # with a minimum at hl = 1.1 - temperature: _VectorisedFloat = inside_temp.value(time) - return np.maximum(1.1, (0.693/((0.16030 + 0.04018*(((temperature-273.15)-20.615)/10.585) + return np.maximum(1.1, (0.693/((0.16030 + 0.04018*(((inside_temp-273.15)-20.615)/10.585) +0.02176*((humidity-45.235)/28.665) -0.14369 - -0.02636*((temperature-273.15)-20.615)/10.585)))) + -0.02636*((inside_temp-273.15)-20.615)/10.585)))) Virus.types = { @@ -924,7 +922,7 @@ class ConcentrationModel: k = (vg * 3600) / h #todo: Inside_temp needs to be exposed/added to the room; return ( - k + self.virus.decay_constant(self.room.humidity, self.room.inside_temp, time) + k + self.virus.decay_constant(self.room.humidity, self.room.inside_temp.value(time)) + self.ventilation.air_exchange(self.room, time) ) From 5dcb698622666c293692e7c986b9e4c74c6d6d4f Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Wed, 27 Apr 2022 15:23:46 +0200 Subject: [PATCH 08/12] Updated humidity value to be in percentage --- cara/models.py | 6 +++--- cara/tests/models/test_exposure_model.py | 12 +++++------ cara/tests/models/test_virus.py | 24 ++++++++++++++++++++++ cara/tests/test_full_algorithm.py | 2 +- cara/tests/test_known_quantities.py | 12 +++++------ cara/tests/test_monte_carlo_full_models.py | 20 +++++++++--------- 6 files changed, 50 insertions(+), 26 deletions(-) create mode 100644 cara/tests/models/test_virus.py diff --git a/cara/models.py b/cara/models.py index cb389be0..333f84bf 100644 --- a/cara/models.py +++ b/cara/models.py @@ -445,7 +445,7 @@ class Virus: raise NotImplementedError def decay_constant(self, humidity: _VectorisedFloat, inside_temp: _VectorisedFloat) -> _VectorisedFloat: - # Viral inactivation per hour (h^-1) (function of humidity) + # Viral inactivation per hour (h^-1) (function of humidity and inside temperature) return np.log(2) / self.halflife(humidity, inside_temp) @@ -459,9 +459,9 @@ class SARSCoV2(Virus): CERN-OPEN-2021-004, DOI: 10.17181/CERN.1GDQ.5Y75) """ # Updated to use the formula from Dabish et al. with correction https://doi.org/10.1080/02786826.2020.1829536 - # with a minimum at hl = 1.1 + # with a minimum at hl = 1.1. Note that humidity is in percentage and inside_temp in °C. return np.maximum(1.1, (0.693/((0.16030 + 0.04018*(((inside_temp-273.15)-20.615)/10.585) - +0.02176*((humidity-45.235)/28.665) + +0.02176*(((humidity*100)-45.235)/28.665) -0.14369 -0.02636*((inside_temp-273.15)-20.615)/10.585)))) diff --git a/cara/tests/models/test_exposure_model.py b/cara/tests/models/test_exposure_model.py index 5a6a5449..b1959ca8 100644 --- a/cara/tests/models/test_exposure_model.py +++ b/cara/tests/models/test_exposure_model.py @@ -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.), 49.60355486823376], + [(1., 1.01), 0.5876509666617122], + [(1.01, 1.02), 0.5726560233646302], + [(12., 12.01), 0.016288631412456556], + [(12., 24.), 716.6229828782525], + [(0., 24.), 777.8717785392312], ] ) 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..0f875fa8 --- /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, 35.67710693238622, 0.01942834607844098], + [272.15, 0.7, 96.40459058258793, 0.007189981061805762], + [300.15, 1., 10.418034697539541, 0.0665333914393324], + [300.15, 0., 1.1, 0.6301338005090411], + [np.array([272.15, 300.15]), np.array([0.7, 0.]), + np.array([96.40459058258793, 1.1]), np.array([0.007189981061805762, 0.6301338005090411])], + [np.array([293.15, 300.15]), np.array([0.5, 1.]), + np.array([35.67710693238622, 10.418034697539541]), np.array([0.01942834607844098, 0.0665333914393324])] + ], +) +def test_decay_constant(inside_temp, humidity, expected_halflife, expected_decay_constant): + npt.assert_equal(models.Virus.types['SARS_CoV_2'].halflife(humidity, inside_temp), + expected_halflife) + npt.assert_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 3f701e5e..1fac1c8b 100644 --- a/cara/tests/test_full_algorithm.py +++ b/cara/tests/test_full_algorithm.py @@ -87,7 +87,7 @@ class SimpleConcentrationModel: return (self.lambda_ventilation + ln2/(np.maximum(1.1, (0.693 / ((0.16030 + 0.04018 * (((21) - 20.615) / 10.585) - + 0.02176*((self.humidity - 45.235) / 28.665) + + 0.02176*(((self.humidity * 100) - 45.235) / 28.665) - 0.14369 - 0.02636*((21-20.615)/10.585))))))) diff --git a/cara/tests/test_known_quantities.py b/cara/tests/test_known_quantities.py index ed775312..98a90acb 100644 --- a/cara/tests/test_known_quantities.py +++ b/cara/tests/test_known_quantities.py @@ -43,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.122276e+01, 1.240684e-12, 2.122276e+01, 7.253047e-26], rtol=1e-6 ) @@ -94,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, 783.490035) def test_periodic_window(baseline_periodic_window, baseline_room): @@ -381,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', 401.300989], + ['Jun', 2420.383151], ], ) def test_exposure_hourly_dep(month,expected_deposited_exposure, baseline_sr_model): @@ -402,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', 402.348745], + ['Jun', 2558.632473], ], ) def test_exposure_hourly_dep_refined(month,expected_deposited_exposure, baseline_sr_model): diff --git a/cara/tests/test_monte_carlo_full_models.py b/cara/tests/test_monte_carlo_full_models.py index 48682758..d60dbcaf 100644 --- a/cara/tests/test_monte_carlo_full_models.py +++ b/cara/tests/test_monte_carlo_full_models.py @@ -310,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", 8.63, 1.64, 8.082, 5624], + ["shared_office_mc", 6.75, 0.20, 3.594, 809], + ["classroom_mc", 9.58, 1.82, 9.664, 5624], ["ski_cabin_mc", 16.0, 0.47, 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], + ["skagit_chorale_mc",72.20, 43.32, 134.234, 190422], + ["bus_ride_mc", 14.08, 9.43, 9.26, 5419], + ["gym_mc", 0.49, 0.14, 0.226, 1145], + ["waiting_room_mc", 2.02, 0.28, 1.104, 737], ] ) def test_report_models(mc_model, expected_pi, expected_new_cases, @@ -337,10 +337,10 @@ 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, 162], + ["No mask", "Jul", 11.81, 14.137, 809], + ["Type I", "Jul", 2.33, 1.305, 149], + ["FFP2", "Jul", 0.73, 0.351, 149], + ["Type I", "Feb", 0.62, 0.291, 162], ], ) def test_small_shared_office_Geneva(mask_type, month, expected_pi, From f518a893d506058968c1cd722a6ea5e7a04fbc05 Mon Sep 17 00:00:00 2001 From: jdevine Date: Mon, 9 May 2022 14:04:27 +0200 Subject: [PATCH 09/12] Corrected handling of min-1/hr-1 half life Changed half life 'limit' value to 6.43 (as paper model, used when formula gives negative values). --- cara/models.py | 9 ++++++--- 1 file changed, 6 insertions(+), 3 deletions(-) diff --git a/cara/models.py b/cara/models.py index 06f45d25..faee2e3f 100644 --- a/cara/models.py +++ b/cara/models.py @@ -459,11 +459,14 @@ class SARSCoV2(Virus): CERN-OPEN-2021-004, DOI: 10.17181/CERN.1GDQ.5Y75) """ # Updated to use the formula from Dabish et al. with correction https://doi.org/10.1080/02786826.2020.1829536 - # with a minimum at hl = 1.1. Note that humidity is in percentage and inside_temp in °C. - return np.maximum(1.1, (0.693/((0.16030 + 0.04018*(((inside_temp-273.15)-20.615)/10.585) + # with a maximum at hl = 6.43. Note that humidity is in percentage and inside_temp in °C. + hl_calc = ((0.693/((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)))) + -0.02636*((inside_temp-273.15)-20.615)/10.585)))/60) + if (hl_calc <= 0): + hl_calc = 6.43 + return hl_calc Virus.types = { From 1e689af19c8a34acb1aedd5d56165ee8ece62a96 Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Mon, 9 May 2022 16:29:52 +0200 Subject: [PATCH 10/12] updated half-life formula and respective tests --- cara/models.py | 11 ++++++----- cara/tests/models/test_exposure_model.py | 12 ++++++------ cara/tests/models/test_virus.py | 20 ++++++++++---------- cara/tests/test_full_algorithm.py | 9 +++++---- cara/tests/test_known_quantities.py | 12 ++++++------ cara/tests/test_monte_carlo_full_models.py | 20 ++++++++++---------- 6 files changed, 43 insertions(+), 41 deletions(-) diff --git a/cara/models.py b/cara/models.py index faee2e3f..34c0d666 100644 --- a/cara/models.py +++ b/cara/models.py @@ -459,14 +459,15 @@ class SARSCoV2(Virus): CERN-OPEN-2021-004, DOI: 10.17181/CERN.1GDQ.5Y75) """ # 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. Note that humidity is in percentage and inside_temp in °C. - hl_calc = ((0.693/((0.16030 + 0.04018*(((inside_temp-273.15)-20.615)/10.585) + # 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) - if (hl_calc <= 0): - hl_calc = 6.43 - return hl_calc + + return np.where(hl_calc <= 0, 6.43, np.minimum(6.43, hl_calc)) Virus.types = { diff --git a/cara/tests/models/test_exposure_model.py b/cara/tests/models/test_exposure_model.py index b1959ca8..896c1c6b 100644 --- a/cara/tests/models/test_exposure_model.py +++ b/cara/tests/models/test_exposure_model.py @@ -179,12 +179,12 @@ def sr_model(): @pytest.mark.parametrize( ["exposed_time_interval", "expected_deposited_exposure"], [ - [(0., 1.), 49.60355486823376], - [(1., 1.01), 0.5876509666617122], - [(1.01, 1.02), 0.5726560233646302], - [(12., 12.01), 0.016288631412456556], - [(12., 24.), 716.6229828782525], - [(0., 24.), 777.8717785392312], + [(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 index 0f875fa8..83f4d1f0 100644 --- a/cara/tests/models/test_virus.py +++ b/cara/tests/models/test_virus.py @@ -7,18 +7,18 @@ from cara import models @pytest.mark.parametrize( "inside_temp, humidity, expected_halflife, expected_decay_constant", [ - [293.15, 0.5, 35.67710693238622, 0.01942834607844098], - [272.15, 0.7, 96.40459058258793, 0.007189981061805762], - [300.15, 1., 10.418034697539541, 0.0665333914393324], - [300.15, 0., 1.1, 0.6301338005090411], + [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([96.40459058258793, 1.1]), np.array([0.007189981061805762, 0.6301338005090411])], + 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([35.67710693238622, 10.418034697539541]), np.array([0.01942834607844098, 0.0665333914393324])] + 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_equal(models.Virus.types['SARS_CoV_2'].halflife(humidity, inside_temp), - expected_halflife) - npt.assert_equal(models.Virus.types['SARS_CoV_2'].decay_constant(humidity, inside_temp), - expected_decay_constant) \ No newline at end of file + 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 7afbd1ce..53a4e99e 100644 --- a/cara/tests/test_full_algorithm.py +++ b/cara/tests/test_full_algorithm.py @@ -84,12 +84,13 @@ class SimpleConcentrationModel: """ removal rate lambda in h^-1, excluding the deposition rate. """ + 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.maximum(1.1, (0.693 / ((0.16030 + 0.04018 * (((21) - 20.615) / 10.585) - + 0.02176*(((self.humidity * 100) - 45.235) / 28.665) - - 0.14369 - - 0.02636*((21-20.615)/10.585))))))) + + ln2/(np.where(hl_calc <= 0, 6.43, np.minimum(6.43, hl_calc)))) @method_cache def deposition_removal_coefficient(self) -> float: diff --git a/cara/tests/test_known_quantities.py b/cara/tests/test_known_quantities.py index 98a90acb..055c587b 100644 --- a/cara/tests/test_known_quantities.py +++ b/cara/tests/test_known_quantities.py @@ -43,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, 2.122276e+01, 1.240684e-12, 2.122276e+01, 7.253047e-26], + [0.000000e+00, 2.046096e+01, 3.846725e-13, 2.046096e+01, 7.231966e-27], rtol=1e-6 ) @@ -94,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, 783.490035) + npt.assert_allclose(r0, 771.380385) def test_periodic_window(baseline_periodic_window, baseline_room): @@ -381,8 +381,8 @@ def build_exposure_model(concentration_model, short_range_model): @pytest.mark.parametrize( "month, expected_deposited_exposure", [ - ['Jan', 401.300989], - ['Jun', 2420.383151], + ['Jan', 359.140499], + ['Jun', 1385.917562], ], ) def test_exposure_hourly_dep(month,expected_deposited_exposure, baseline_sr_model): @@ -402,8 +402,8 @@ def test_exposure_hourly_dep(month,expected_deposited_exposure, baseline_sr_mode @pytest.mark.parametrize( "month, expected_deposited_exposure", [ - ['Jan', 402.348745], - ['Jun', 2558.632473], + ['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_full_models.py b/cara/tests/test_monte_carlo_full_models.py index 7d75b2d3..3bce07ae 100644 --- a/cara/tests/test_monte_carlo_full_models.py +++ b/cara/tests/test_monte_carlo_full_models.py @@ -310,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.75, 0.20, 3.594, 809], - ["classroom_mc", 9.58, 1.82, 9.664, 5624], + ["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",72.20, 43.32, 134.234, 190422], - ["bus_ride_mc", 14.08, 9.43, 9.26, 5419], - ["gym_mc", 0.49, 0.14, 0.226, 1145], - ["waiting_room_mc", 2.02, 0.28, 1.104, 737], + ["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, @@ -337,10 +337,10 @@ 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", 11.81, 14.137, 809], - ["Type I", "Jul", 2.33, 1.305, 149], - ["FFP2", "Jul", 0.73, 0.351, 149], - ["Type I", "Feb", 0.62, 0.291, 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, From 0b9582ce310a698338a631077ce1e47c8c8a46f6 Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Tue, 10 May 2022 11:41:19 +0200 Subject: [PATCH 11/12] exposure of the inside_temp and humidity to the API --- cara/apps/calculator/model_generator.py | 15 +++++++++++---- cara/apps/calculator/static/js/form.js | 3 +++ cara/apps/templates/base/calculator.form.html.j2 | 2 ++ 3 files changed, 16 insertions(+), 4 deletions(-) diff --git a/cara/apps/calculator/model_generator.py b/cara/apps/calculator/model_generator.py index 08cbc247..dd7f1e37 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, inside_temp=models.PiecewiseConstant((0, 24), (293,)), 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() 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/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 @@ + +
From 96af7827370f342b11757992e33b0f4381f613a3 Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Tue, 10 May 2022 11:56:49 +0200 Subject: [PATCH 12/12] Added humidity and inside_temp to the default values for the baseline form data --- cara/apps/calculator/model_generator.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/cara/apps/calculator/model_generator.py b/cara/apps/calculator/model_generator.py index dd7f1e37..c7ebd410 100644 --- a/cara/apps/calculator/model_generator.py +++ b/cara/apps/calculator/model_generator.py @@ -693,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', @@ -702,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',