Corona simulation challenge

[hellomike]

Just a thought.

Washington Post Article

Is there anyone smart enough to code the simulation using BB4W and (or) BBC BASIC for SDL 2.0 as is shown in the mentioned site?
I myself would probably tackle it somehow but it would take a very long time and certainly not in the most effective way as I have no experience using the graphics library capabilities.

Mike

[KenDown]

Unfortunately one needs to be a subscriber to use the link you have provided.

[RichardRussell]

Is there anyone smart enough to code the simulation using BB4W and (or) BBC BASIC for SDL 2.0 as is shown in the mentioned site?

I don't think coding it would be a challenge, it's just random motion and collision detection. To make it look as nice would mean using anti-aliased graphics for sure, and checking collisions between 200 objects is not going to be quick, so whether it could be made to run in 'real time' in BB4W or BBCSDL is questionable.

I just wish the real virus had the characteristic of those simulations: everybody recovers and nobody dies! Perhaps a BBC BASIC version could be made more realistic (and depressing).

[RichardRussell]

Here's my attempt, in BBCSDL (sorry, it's not BB4W compatible although it wouldn't be that difficult to convert), which reproduces the simplest of the simulations. I've concentrated on the graphical aspects (so everything drawn antialiased for quality) and it needs a fairly fast PC to run at full rate:

Code: Select all

      INSTALL @lib$ + "aagfxlib"
      VDU 23,22,800;460;8,20,16,128

      PEOPLE = 200
      SIZE = 10
      COLL = SIZE^2
      SPEED = 5
      HEALTHY = 0
      SICK = 1
      DURATION = 300
      RECOVERED = 2
      DEAD = 3

      LEDGE = SIZE
      REDGE = 1600 - SIZE
      BEDGE = SIZE
      TEDGE =  800 - SIZE

      LINE 0, 800, 1600, 800
      VDU 28, 0, 22, 99, 3, 5, 30

      OSCLI "font """ + @lib$ + "DejaVuSans"", 11"
      PRINT "  Recovered" ' "  Healthy" ' "  Sick" ;
      OSCLI "font """ + @lib$ + "DejaVuSans"", 11, B"

      DIM x(PEOPLE), y(PEOPLE), u(PEOPLE), v(PEOPLE), f%(PEOPLE), t%(PEOPLE)
      DIM col%(3), num%(3)
      col%() = &FFCAC6AA, &FF1D64BB, &FFC08ACB, &00000000
      COLOUR 1,&AA,&C6,&CA : COLOUR 2,&BB,&64,&1D : COLOUR 3,&CB,&8A,&C0

      ON ERROR OSCLI "REFRESH ON" : CLS : REPORT : END

      REM Initialise to random positions and velocities:
      FOR I% = 1 TO PEOPLE
        direction = 2 * PI * RND(1)
        x(I%) = (1598 - 2*SIZE) * RND(1) + SIZE
        y(I%) = ( 798 - 2*SIZE) * RND(1) + SIZE
        u(I%) = SPEED * COS(direction)
        v(I%) = SPEED * SIN(direction)
      NEXT

      REM Infect one person:
      f%(RND(PEOPLE)) = SICK

      REM Animate
      *REFRESH OFF
      graph = 300
      REPEAT
        CLS
        num%() = 0

        REM Update positions:
        FOR I% = 1 TO PEOPLE
          x = x(I%) : y = y(I%)
          x += u(I%) : IF x > LEDGE IF x < REDGE x(I%) = x ELSE u(I%) *= -1
          y += v(I%) : IF y > BEDGE IF y < TEDGE y(I%) = y ELSE v(I%) *= -1
        NEXT

        REM Check for collisions:
        FOR I% = 1 TO PEOPLE - 1
          x = x(I%) : y = y(I%)
          FOR J% = I% + 1 TO PEOPLE
            IF (x-x(J%))^2 + (y-y(J%))^2 < COLL THEN
              u(I%) *= -1 : v(I%) *= -1 : u(J%) *= -1 : v(J%) *= -1
              x(I%) += u(I%) : y(I%) += v(I%) : x(J%) += u(J%) : y(J%) += v(J%)
              IF f%(I%) = SICK IF f%(J%) = HEALTHY f%(J%) = SICK
              IF f%(J%) = SICK IF f%(I%) = HEALTHY f%(I%) = SICK
            ENDIF
          NEXT
        NEXT I%

        REM Display and count:
        FOR I% = 1 TO PEOPLE
          PROC_aasector(x(I%), y(I%), SIZE, SIZE, 0, 360, col%(f%(I%)))
          IF f%(I%) = SICK t%(I%) += 1 : IF t%(I%) > DURATION f%(I%) = RECOVERED
          num%(f%(I%)) += 1
        NEXT

        REM Report and graph:
        GCOL 15 : RECTANGLE FILL 200, 820, 80, 100
        GCOL 3 : MOVE 200,918 : PRINT ;num%(2);
        GCOL 1 : MOVE 200,882 : PRINT ;num%(0);
        GCOL 2 : MOVE 200,846 : PRINT ;num%(1);

        y0 = 810
        y1 = y0 + num%(1) / 2
        y2 = y1 + num%(0) / 2
        y3 = y2 + num%(2) / 2
        IF y1 <> y0 PROC_aaline(graph, y0, graph, y1, 2, col%(1), 0)
        IF y2 <> y1 PROC_aaline(graph, y1, graph, y2, 2, col%(0), 0)
        IF y3 <> y2 PROC_aaline(graph, y2, graph, y3, 2, col%(2), 0)
        graph += 1

        *REFRESH
      UNTIL FALSE
      END

[RichardRussell]

Here's a slightly modified version which includes a death rate of 3%.

Code: Select all

      INSTALL @lib$ + "aagfxlib"
      VDU 23,22,800;480;8,20,16,128

      PEOPLE = 200
      SIZE = 10
      COLL = 4*SIZE^2
      SPEED = 5
      FATAL = 0.03
      DURATION = 300
      HEALTHY = 0
      SICK = 1
      RECOVERED = 2
      DEAD = 3

      LEDGE = SIZE
      REDGE = 1600 - SIZE
      BEDGE = SIZE
      TEDGE =  800 - SIZE

      LINE 0, 800, 1600, 800
      VDU 28, 0, 23, 99, 4, 5, 30

      OSCLI "font """ + @lib$ + "DejaVuSans"", 11"
      PRINT "  Dead" ' "  Recovered"' "  Healthy" ' "  Sick" ;
      OSCLI "font """ + @lib$ + "DejaVuSans"", 11, B"

      DIM x(PEOPLE), y(PEOPLE), u(PEOPLE), v(PEOPLE), f%(PEOPLE), t%(PEOPLE)
      DIM col%(3), num%(3)
      col%() = &FFCAC6AA, &FF1D64BB, &FFC08ACB, &00000000
      COLOUR 1,&AA,&C6,&CA : COLOUR 2,&BB,&64,&1D : COLOUR 3,&CB,&8A,&C0

      ON ERROR OSCLI "REFRESH ON" : CLS : REPORT : END

      REM Initialise to random positions and velocities:
      FOR I% = 1 TO PEOPLE
        direction = 2 * PI * RND(1)
        x(I%) = (1598 - 2*SIZE) * RND(1) + SIZE
        y(I%) = ( 798 - 2*SIZE) * RND(1) + SIZE
        u(I%) = SPEED * COS(direction)
        v(I%) = SPEED * SIN(direction)
      NEXT

      REM Infect one person:
      f%(RND(PEOPLE)) = SICK

      REM Animate
      *REFRESH OFF
      graph = 300
      REPEAT
        CLS
        num%() = 0

        REM Update positions, display and count:
        FOR I% = 1 TO PEOPLE
          x = x(I%) : y = y(I%)
          x += u(I%) : IF x > LEDGE IF x < REDGE x(I%) = x ELSE u(I%) *= -1
          y += v(I%) : IF y > BEDGE IF y < TEDGE y(I%) = y ELSE v(I%) *= -1
          PROC_aasector(x(I%), y(I%), SIZE, SIZE, 0, 360, col%(f%(I%)))
          IF f%(I%) = SICK THEN
            t%(I%) += 1
            IF t%(I%) > DURATION f%(I%) = RECOVERED : IF RND(1) < FATAL f%(I%) = DEAD
          ENDIF
          num%(f%(I%)) += 1
        NEXT

        REM Check for collisions:
        FOR I% = 1 TO PEOPLE - 1
          x = x(I%) : y = y(I%)
          FOR J% = I% + 1 TO PEOPLE
            d = (x-x(J%))^2 + (y-y(J%))^2
            IF d < COLL IF f%(I%) <> DEAD IF f%(J%) <> DEAD THEN
              IF d < (x-u(I%)-x(J%)+u(J%))^2 + (y-v(I%)-y(J%)+v(J%))^2 THEN
                u(I%) *= -1 : v(I%) *= -1
                u(J%) *= -1 : v(J%) *= -1
              ENDIF
              IF f%(I%) = SICK IF f%(J%) = HEALTHY f%(J%) = SICK
              IF f%(J%) = SICK IF f%(I%) = HEALTHY f%(I%) = SICK
            ENDIF
          NEXT
        NEXT I%

        REM Report and graph:
        GCOL 15 : RECTANGLE FILL 200, 820, 80, 140
        GCOL 0 : MOVE 200,958 : PRINT ;num%(3);
        GCOL 3 : MOVE 200,922 : PRINT ;num%(2);
        GCOL 1 : MOVE 200,886 : PRINT ;num%(0);
        GCOL 2 : MOVE 200,846 : PRINT ;num%(1);

        y0 = 830
        y1 = y0 + num%(1) / 2
        y2 = y1 + num%(0) / 2
        y3 = y2 + num%(2) / 2
        y4 = y3 + num%(3) / 2
        IF y1 <> y0 PROC_aaline(graph, y0, graph, y1, 2, col%(1), 0)
        IF y2 <> y1 PROC_aaline(graph, y1, graph, y2, 2, col%(0), 0)
        IF y3 <> y2 PROC_aaline(graph, y2, graph, y3, 2, col%(2), 0)
        IF y4 <> y3 PROC_aaline(graph, y3, graph, y4, 2, &FF000000, 0)
        graph += 1

        *REFRESH
      UNTIL FALSE
      END

https://www.youtube.com/watch?v=I8XTXDDJ1lM

[hellomike]

Impressive attempt Richard and super fast.

For myself I added a NOMOVE value for people that don't move around and changed code as follows:

Code: Select all

        IF RND(PEOPLE / NOMOVE) > 1 THEN
          u(I%) = SPEED * COS(direction)
          v(I%) = SPEED * SIN(direction)
        ENDIF

(which fails when patient zero happens to not move around )

Yes, it's a pretty crude simulation but must have been a fun exercise for you.
Thanks again.

Mike

[RichardRussell]

Before somebody else points it out, the line:

Code: Select all

              u(I%) *= -1 : v(I%) *= -1 : u(J%) *= -1 : v(J%) *= -1

should be:

Code: Select all

              SWAP u(I%), u(J%) : SWAP v(I%), v(J%)

but the simulation of the spread of the virus isn't affected.

[RichardRussell]

For myself I added a NOMOVE value for people that don't move around

A nice modification, but unfortunately it can seriously slow down the program if significant numbers of 'non moving' people are closer together than the collision threshold (they seem to collide every frame!). The collision maths has to be changed too, otherwise a collision between a moving person and a non-moving person isn't physically realistic.

It sounds as though it ought to be a simple change, but it's not!

[hellomike]

Yes Richard, you are right. I see what you mean.
I wonder if the collision math can be more like:

Code: Select all

          newx = x(I%) + u(I%) : newy = y(I%) + v(I%)
          IF POINT(newx + SIZE, newy + SIZE) <> 15 THEN
            REM I'm about to bump into someone at the new position....
            ....
          ENDIF
          ....

Did some testing in that direction but without succes.

[RichardRussell]

I wonder if the collision math can be more like:

POINT() is far too slow (in BBCSDL) to be used in that kind of situation, it has to be used very sparingly or not at all. If you want to reproduce the Washington Post simulation accurately, which presumably is the objective, you will need to:

1. Ensure that 'patient zero' is moving (that at least is easy).
   
   
2. Ignore 'collisions' (overlaps) between two non-moving people because they just waste time. That will need to be done with care so as not to slow it down too much.
   
   
3. Calculate the rebound angle differently for collisions with non-moving people. Hitting 'head on' will cause a rebound along the original trajectory, hitting slightly to the left will cause a rebound in that direction and hitting to the right similarly.

Good luck with that!