#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#define N 128          /* taille du reseau */
#define M 4000         /* nombre de marcheurs */
#define L 20           /* duree de survie des malades */
#define NITERMAX 5000  /* nombre maximal d'iterations temporelles */
int main(void)
{
/* Declarations ============================================================ */
  int x[M], y[M], infect[M], survie[M] ;
  int ninfect, nmort, iter, j, k, jj, xnew, ynew, d2 ;
/* Executable ============================================================== */
  for (j=0 ; j<M ; j++) {                       /* M marcheurs sur le reseau */
    x[j]     =floor(1.*N*rand()/RAND_MAX) ;
    y[j]     =floor(1.*N*rand()/RAND_MAX) ;
    infect[j]=0 ; survie[j]=0 ;
  }
  jj        =1.*M*rand()/RAND_MAX ;             /* infection d'un marcheur */
  infect[jj]=1 ; survie[jj]=L ; ninfect=1 ; nmort=0 ; iter=0 ;
  while ( ninfect > 0 && iter < NITERMAX ) {    /* Iteration temporelle */
    for (j=0 ; j<M ; j++ ) {                    /* boucle sur les marcheurs */
      if ( infect[j] == 1 ) { survie[j] -=1 ; } /* le temps qui passe... */
      if ( infect[j] == 1 && survie[j] <= 0 ) { /* ...les malades qui meurent */
         infect[j]=2 ; ninfect -=1 ; nmort   +=1 ;
      }
      if ( infect[j] < 2 ) {                    /* seuls les vivants bougent */
        if ( 1.*rand()/RAND_MAX < 0.5 ) {       /* on bouge en x */
          xnew= x[j] + 2* (floor( 2.*rand()/RAND_MAX ))-1 ; 
          if ( xnew < 0 )   { xnew=0   ; }      /* On ne quitte pas le reseau */
          if ( xnew > N-1 ) { xnew=N-1 ; }
          x[j]=xnew ;                           /* deplacement du marcheur */
        } else {                                /* sinon on bouge en y */
          ynew= y[j] + 2* (floor( 2.*rand()/RAND_MAX ))-1 ; 
          if ( ynew < 0 )   { ynew=0   ; }      /* on ne quitte pas le reseau */
          if ( ynew > N-1 ) { ynew=N-1 ; }
          y[j]=ynew ;                           /* deplacement du marcheur */
        }
      }
      if ( infect[j] == 1 ) {                   /* ce marcheur est malade...  */
        for ( k=0 ; k<M ; k++ ) {               /* ...et peut contaminer les autres */
          if ( k != j && infect[k] == 0 ) {     /* ...qui ne le sont pas deja */
            d2=pow( (x[j]-x[k]),2 )+pow( (y[j]-y[k]),2 ) ;  
            if ( d2 == 0 ) {                    /* sur le meme site: contagion */
              infect[k]=1 ; ninfect +=1 ; survie[k]=L ;
            }
          }
        }
      }
    }
    printf ("iter, malades, morts: %d, %d, %d\n",iter,ninfect,nmort) ; iter+=1 ; 
  }
}