From 4722a5511c139599265a7fcc7a59a7bba1614f56 Mon Sep 17 00:00:00 2001 From: Phil Elson Date: Mon, 26 Oct 2020 20:12:54 +0100 Subject: [PATCH 1/2] Simplify the concentration calculation. --- cara/models.py | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/cara/models.py b/cara/models.py index 1e2ca4bc..dd8262bb 100644 --- a/cara/models.py +++ b/cara/models.py @@ -266,6 +266,13 @@ class Model: return k + self.virus.decay_constant + self.ventilation.air_exchange(self.room, time) + def collect_time_state_changes(self): + """ + All time dependent entities on this model must provide information about + the times at which their state changes. + + """ + @functools.lru_cache() def concentration(self, time: float) -> float: t = time @@ -286,4 +293,5 @@ class Model: else: # Concentration while infected not present. end_concentration = self.concentration(t1) - return (end_concentration + ((np.exp(IVRR * t1) - 1) * end_concentration)) * np.exp(-IVRR * t) + fac = np.exp(IVRR * (t1 - t)) + return end_concentration * fac From 9e9f1a28bc8eef815520bde63ef5eb58086a68b9 Mon Sep 17 00:00:00 2001 From: Phil Elson Date: Mon, 26 Oct 2020 20:21:15 +0100 Subject: [PATCH 2/2] Fix a bug in the calculation of the concentration. --- cara/models.py | 5 ++--- cara/tests/test_known_quantities.py | 2 +- 2 files changed, 3 insertions(+), 4 deletions(-) diff --git a/cara/models.py b/cara/models.py index dd8262bb..69f642a6 100644 --- a/cara/models.py +++ b/cara/models.py @@ -287,9 +287,8 @@ class Model: elif t0 < t <= t1: # Concentration while infected present. init_concentration = self.concentration(t0) - time_present = t - t0 - fac = np.exp(-IVRR * time_present) - return ((ER + Ni) / (IVRR * V)) * (1 - fac) + init_concentration * fac + fac = np.exp(IVRR * (t0 - t)) + return ((ER * Ni) / (IVRR * V)) * (1 - fac) + init_concentration * fac else: # Concentration while infected not present. end_concentration = self.concentration(t1) diff --git a/cara/tests/test_known_quantities.py b/cara/tests/test_known_quantities.py index 8976a558..b87a318b 100644 --- a/cara/tests/test_known_quantities.py +++ b/cara/tests/test_known_quantities.py @@ -64,7 +64,7 @@ def test_r0(baseline_model): concentrations = [baseline_model.concentration(t) for t in ts] npt.assert_allclose( concentrations, - [0.000000e+00, 2.909211e-01, 1.273836e-04, 2.909210e-01, 5.577662e-08], + [0.000000e+00, 2.891970e-01, 1.266287e-04, 2.891969e-01, 5.544607e-08], rtol=1e-5 )