#include <stdio.h>
#include <math.h>

double a = 0.0002;
double T = 12;

double f_s(double t, double s, double i)
{
   return - a*s*i;
}

double f_i(double t, double s, double i)
{
   return  a*s*i - i/T;
}

int main()
{
   double h=0.05, t0=0.0, tk=200, t;
   double n = 1000;    /* populacja */
   double i0 = 1;      /* jeden zainfekowany */ 
   double s0 = n - i0; /* liczba podatnych */ 
   double i = i0, s = s0;
   double k_s1, k_s2, k_i1, k_i2;
   
   printf("# t      S         I       R\n");
   for (t=t0; t<= tk; t=t+h)
   {
      printf("%6.2f  %6.1f  %6.1f  %6.1f\n", t, s, i, n-s-i);
      
      k_s1 = f_s(t, s, i);
      k_i1 = f_i(t, s, i);
      
      k_s2 = f_s(t+h, s + k_s1 * h, i + k_i1 * h);
      k_i2 = f_i(t+h, s + k_s1 * h, i + k_i1 * h);
      
      s = s + 0.5 * (k_s1 + k_s2 ) * h;
      i = i + 0.5 * (k_i1 + k_i2 ) * h;      
   }

   return 0;
}