Introduction

This is a very brief post, with a simulation written in JavaScript. It’s a very simple model of how a disease can spread into a population. You may get the code, adjust the parameters and watch different outcomes.

The idea of this project came from this: Why outbreaks like coronavirus spread exponentially, and how to ‘flatten the curve’. I highly recommend reading it.

Theory

I recommend watching this video, for some theory:

It’s a good start, there is a lot more information about this on the internet, so I won’t bother with formulae here.

Simulation

The model resembles a lot the one from the link pointed above, with an addition: this model also simulates death, not all ‘people’ recover.

Briefly, it’s just molecular dynamics in 2D, in a ‘box’, with elastic collisions. All ‘particles’ have the same mass and they can be either not infected, infected, cured and dead. The dead ones disappear from the simulation.
Once an infected ‘particle’ hits one that wasn’t infected yet, it will infect it. A ‘particle’ once having the illness, has a chance of dying (but only in the last two thirds of the infected period of time). After being sick for a while, if it doesn’t die in the meantime, it will recover and have immunity, so it cannot be infected again.

There is also a chart showing the number of infections over time. You’ll notice that it will level out under 100%, the difference being the ones that did not get an infection, being lucky, and the dead ones.

The code has a feature that it’s usually a bug: because of the too big time step, some ‘particles’ can stick together. Basically when they collide, they overlap a little (more, if the time step is big) and if they cannot go apart in one time step, remaining overlapped, they will be considered as colliding again and so on. There are various ways of solving this, but for this model I consider it a feature: people sometimes do something analogous.

The blue ‘particles’ are not infected, the red ones are infected and the green ones are cured. The ones that die disappear from the canvas.

The computation will be a little slow on a phone, for this, you should better use a desktop.

Related to this, for a better implementation on such molecular dynamics where there is short range interaction only, visit Event Driven Molecular Dynamics. It’s also 3D and with nicer graphics, using OpenGL.

The Code

There are some comments in the code, hopefully it’s clear enough. I wrote it very fast, so it’s far from perfect but it seems to work as intended.

(function () {
    var canvas = document.getElementById("Canvas");
    var chart = document.getElementById("Chart");

    canvas.width = Math.min(window.innerWidth - 50, window.innerHeight - 50);
    if (canvas.width < 300) canvas.width = 300;
    else if (canvas.width > 800) canvas.width = 800;
    canvas.height = canvas.width;

    chart.width = canvas.width;
    chart.height = chart.width / 3;

    var ctx = canvas.getContext("2d");
    var canvasText = document.getElementById("canvasText");

    var ctxChart = chart.getContext("2d");
    ctxChart.strokeStyle = "#FF0000";
    ctxChart.lineWidth = 4;

    var people = [];

    // some parameters
    var nrPeople = 500;
    var speed = canvas.width / 50.;
    var radius = canvas.width / 100.;
    var deltat = 0.01;
    var cureTime = 5.;
    var deathProb = 0.05;
    var infectProb = 1.;

    // chart
    var chartPosX = 0;
    var chartValY = 0;
    var chartXScale = 0.1;
    var chartYScale = chart.height / nrPeople;

    // statistics
    var deaths = 0;
    var infections = 1;
    var cures = 0;

    function Init() {
        deaths = 0;
        infections = 1;
        cures = 0;
        people = [];
        chartPosX = 0;
        chartValY = 0;

        for (i = 0; i < nrPeople; ++i) {
            var person =
            {
                posX: 0,
                posY: 0,
                velX: 0,
                velY: 0,
                dead: false,
                infected: false,
                cured: false,
                infectedTime: 0.0,

                NormalizeVelocity: function () {
                    var len = Math.sqrt(this.velX * this.velX + this.velY * this.velY);
                    this.velX /= len;
                    this.velY /= len;
                },

                Distance: function (other) {
                    var distX = this.posX - other.posX;
                    var distY = this.posY - other.posY;

                    return Math.sqrt(distX * distX + distY * distY);
                },

                Collision: function (other) {
                    var dist = this.Distance(other);

                    return dist <= 2. * radius;
                },

                Collide: function (other) {
                    var velXdif = this.velX - other.velX;
                    var velYdif = this.velY - other.velY;

                    var posXdif = this.posX - other.posX;
                    var posYdif = this.posY - other.posY;

                    var dist2 = posXdif * posXdif + posYdif * posYdif;
                    var dotProd = velXdif * posXdif + velYdif * posYdif;
                    dotProd /= dist2;

                    this.velX -= dotProd * posXdif;
                    this.velY -= dotProd * posYdif;

                    other.velX += dotProd * posXdif;
                    other.velY += dotProd * posYdif;
                }
            };

            for (; ;) {
                var X = Math.floor(Math.random() * (canvas.width - 2. * radius)) + radius;
                var Y = Math.floor(Math.random() * (canvas.height - 2 * radius)) + radius;

                person.posX = X;
                person.posY = Y;

                overlap = false;
                for (j = 0; j < i; ++j) {
                    var person2 = people[j];
                    if (person2.Collision(person)) {
                        overlap = true;
                        break;
                    }
                }
                if (!overlap) break;
            }

            person.velX = Math.random() - 0.5;
            person.velY = Math.random() - 0.5;
            person.NormalizeVelocity();
            person.velX *= speed;
            person.velY *= speed;

            if (i == 0) person.infected = true;

            people.push(person);
        }
    }

    Init();

    function Advance() {
        // for each from the population
        for (i = 0; i < nrPeople; ++i) {
            var person = people[i];
            if (person.dead) continue;

            // move

            person.posX += person.velX * deltat;
            person.posY += person.velY * deltat;

            // collide / infect / cure

            // first, walls collisions
            if (person.posX <= radius) {
                person.velX *= -1;
                person.posX = radius;
            }
            else if (person.posX >= canvas.width - radius) {
                person.velX *= -1;
                person.posX = canvas.width - radius;
            }

            if (person.posY <= radius) {
                person.velY *= -1;
                person.posY = radius;
            }
            else if (person.posY >= canvas.height - radius) {
                person.velY *= -1;
                person.posY = canvas.height - radius;
            }

            // keep track of how long the infection lasts
            if (person.infected)
                person.infectedTime += deltat;

            // collisions and infections between them
            for (j = 0; j < i; ++j) {
                var person2 = people[j];
                if (person2.dead) continue;

                if (person.Collision(person2)) {
                    person.Collide(person2);

                    if (person.infected && !person2.infected && !person2.cured) {
                        if (Math.random() < infectProb)
                        {
                           person2.infected = true;
                           ++infections;
                        }
                    }
                    else if (person2.infected && !person.infected && !person.cured) {
                        if (Math.random() < infectProb)
                        {
                           person.infected = true;
                           ++infections;
                        }
                    }
                }
            }

            // cure
            if (person.infected && person.infectedTime > cureTime) {
                person.infected = false;
                person.cured = true;
                ++cures;
            }

            // kill
            if (person.infected && person.infectedTime > cureTime / 3.) {
                if (Math.random() < deathProb * deltat / (cureTime * 2. / 3.)) {
                    person.dead = true;
                    ++deaths;
                }
            }
        }

        // display
        ctx.clearRect(0, 0, canvas.width, canvas.height);
        for (i = 0; i < nrPeople; ++i) {
            var person = people[i];
            if (person.dead) continue;
            ctx.beginPath();
            ctx.arc(person.posX, person.posY, radius, 0, 2 * Math.PI, false);
            ctx.fillStyle = person.infected ? 'red' : (person.cured ? 'green' : 'blue');
            ctx.fill();
            ctx.stroke();
        }

        canvasText.innerHTML = "Total: " + nrPeople + " Infected: " + infections + 
        " Deaths: " + deaths + " Cured: " + cures + " Sick: " + (infections - cures - deaths);

        ctxChart.beginPath();
        ctxChart.moveTo(chartPosX * chartXScale, chart.height - chartValY * chartYScale);
        ++chartPosX;  
        chartValY = infections;      
        ctxChart.lineTo(chartPosX * chartXScale, chart.height - chartValY * chartYScale);
        ctxChart.stroke();
       
        if (infections - deaths == cures) // there is no more left to cure
        {
            ctxChart.clearRect(0, 0, chart.width, chart.height);
            Init();
        }
    }

    setInterval(Advance, 10);
})();

Conclusion

As usual, if you find mistakes or have suggestions, please point them out.

I’ll end up with a quote from Will the pandemic kill you? by Nobel Laureate Peter Doherty:

“What influenza also tells us is that, while viruses that spread via the respiratory route are the most likely cause of some future pandemic, only the most draconian and immediate restrictions on human travel are likely to limit the spread of infection, and then only briefly.”

“Such limitations are likely to be applied quickly if we are faced with a situation in which, for example, more than 30% of those affected develop severe or even fatal consequences. The more dangerous situation may be when mortality rates are in the 1-3% range, causing (ultimately) 70 million to 210 million deaths globally. Such an infection could “get away” before we realised what was happening.”

Originally appeared on: Epidemics | Computational Physics (go.ro). Follow this link to see the code in action.

Epidemic - Simple Model of How a Disease Can Spread into a Population