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/async
.
@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.