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

void sdot(double *s1, const double *s)
{
  double x, y, vx, vy;
  double ax, ay;
  const double g = 9.8;

  x = s[0]; y = s[1]; vx = s[2]; vy = s[3];
  ay = - g * (x*x)/(x*x + y * y) - (vx * vx + vy * vy) * y / (x * x + y * y);
  ax =   g * (x*y)/(x*x + y * y) - (vx * vx + vy * vy) * x / (x * x + y * y); 
  s1[0] = vx;
  s1[1] = vy;
  s1[2] = ax;
  s1[3] = ay; 
}

void aggiorna(double *s, double dt)
{
  double s1[4];
  unsigned int i;

  sdot(s1, s);

  for (i = 2; i < 4; i++) {
    s[i] = s[i] + s1[i] * dt;
  }
  for (i = 0; i < 2; i++) {
    s[i] = s[i] + (s1[i] + s[i + 2]) / 2 * dt;
  }
}

int main()
{
  double dt = 0.01, t;
  double s[4] = {5, 0, 0, 0};

  t = 0;
  while (t < 10) {
    printf("%f %f %f\n", t, s[0], s[1]);

    aggiorna(s, dt);
    t = t + dt;
  }

  return 0;
}