I am working on an epidemiological agent-based model. Suppose I have an array of Agents
mutable struct Human isinfected end agents = [Human(i) for i = 1:10000]
When the simulation starts, I introduce a random infected person into the population, i.e.
setup_random_infected(), after which the main discrete-time simulation starts:
## main simulation loop. for t=1:SIM_TIME contact_dynamics() infection_dynamics() data_collection() end
Problem: Sometimes the initial infected person recovers without making anyone else sick, due to stochascity or contact patterns. When this happens, I would like to kill the simulation (i.e. return from the function) without running through the entire for loop.
I can have an
if statement, but i fear that this is too costly. Sure, when the initial person does recover and the basic reproduction number R0 is 0, this kills the loop. However, in the simulations where the epidemic takes off, I fear that constantly evaluating
if statement is costly. This is because the
array of humans (a “human” is a large
mutable struct) is heap-allocated, and checking this person’s property,
human.isinfected constantly is not going to be optimized.
## main simulation loop. for t=1:SIM_TIME if initial_person == recover && R0 == 0 return contact_dynamics() infection_dynamics() data_collection() end
I was wondering if there is a way to have the for loop run without an if statement, but still watch for that infected person and kill the loop if the epidemic dies.
My Idea*: I was thinking something alone the lines of using channels plus
@sync begin @async start_condition_check() @async start_main_time_loop() end
Real problem: Although what I’ve described above is barebones, my real problem is that I am trying to calculate R0 on the fly. Out of, say 2000 Monte Carlo simulations, I want to keep only those simulations running that the R0 is some threshold for the initial infected. If R0 dosn’t reach the threshold, I would like to kill the simulation.
Ideas and tips would be appreciated.