/*
cc seir.c -o bseir -lm
*/

#include <stdlib.h>
#include <stdio.h>

#define	HEALTHY		0
#define	EXPOSED		1
#define	INFECTED	2
#define	IMMUNE		3

#define	EXP_TIME	3
#define	INF_TIME	3
#define	IMM_TIME	3

#define	INF_RATE	(0.35)
#define	INIT_EXP	(0.1)

#define	MAXTIME		100
#define	CELLS		1000

int	state[CELLS], timer[CELLS];
double	inf_rate, scnt, ecnt, icnt, rcnt;

void InitStuff();
void SpreadDisease(void);
void StoE(int, int, int);
void ChangeStates(void);
void Measurements(void);

int main(void)
{
	int i;

	InitStuff();

	Measurements();

	for(i=0;i<MAXTIME;i++)
	{
		SpreadDisease();
		ChangeStates();

		Measurements();
		printf("%d %f %f %f %f\n",
			i, scnt, ecnt, icnt, rcnt);
	}
	return (0);
}

void InitStuff(void)
{
	int i;
	unsigned long seed=145739853;
	srand48(seed);

	inf_rate = INF_RATE;
	for(i=0;i<CELLS;i++)
	{
		state[i] = HEALTHY;
		timer[i] = 0;
		if(drand48()<INIT_EXP)
		{
			state[i] = EXPOSED;
			timer[i] = EXP_TIME * drand48() + 1;
		}
	}
}


void SpreadDisease(void)
{
	int i;

	if(state[0]==HEALTHY) StoE(0,CELLS-1,1);
	for(i=1;i<CELLS-1;i++)
		if(state[i]==HEALTHY) StoE(i,i-1,i+1);
	if(state[CELLS-1]==HEALTHY) StoE(CELLS-1,CELLS-2,0);
}

void StoE(int cell, int cellL, int cellR)
{
	int cnt=0;

	if(state[cellL]==INFECTED) cnt = 1;
	if(state[cellR]==INFECTED) cnt++;
	if(cnt)
	{
		if(drand48()<(double)cnt*inf_rate)
		{
			state[cell] = EXPOSED;
			timer[cell] = EXP_TIME;
		}
	}
}

void ChangeStates(void)
{
	int i;

	for(i=0;i<CELLS;i++)
	{
		switch(state[i])
		{
			case IMMUNE:
				timer[i]--;
				if(timer[i]==0)
					state[i] = HEALTHY;
				break;

			case INFECTED:
				timer[i]--;
				if(timer[i]==0)
				{
					state[i] = IMMUNE;
					timer[i] = IMM_TIME;
				}
				break;

			case EXPOSED:
				timer[i]--;
				if(timer[i]==0)
				{
					state[i] = INFECTED;
					timer[i] = INF_TIME;
				}
				break;

			default:
				break;
		}
	}
}

void Measurements(void)
{
	int i;

	scnt = ecnt = icnt = rcnt = 0.0;

	for(i=0;i<CELLS;i++)
	{
		if(state[i]==HEALTHY) 		scnt += 1.0;
		else if(state[i]==EXPOSED) 	ecnt += 1.0;
		else if(state[i]==INFECTED) 	icnt += 1.0;
		else if(state[i]==IMMUNE) 	rcnt += 1.0;
	}
	scnt /= (double)(CELLS);
	ecnt /= (double)(CELLS);
	icnt /= (double)(CELLS);
	rcnt /= (double)(CELLS);
}