/* Second-order predictor-corrector method x(t). */

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

/* Specify the ODE: d2x/dt2 = acc(x, v, t). */

float acc(float x, float v, float t)
{
    return -x;
}

void output(float x, float v, float t)
{
    printf("%f %f %f\n", t, x, v);
}

main()
{
    float x, v, t, dt;

    /* Initial conditions: */

    t = 0;
    x = 1;
    v = 0;

    /* Time step: */

    dt = 0.1;

    while (t < 10) {

	float xpred, vpred, a, apred;

	output(x, v, t);

	/* Predict x and v. */

	a = acc(x, v, t);

	xpred = x + v*dt + 0.5*a*dt*dt;
	vpred = v + a*dt;

	/* Re-evaluate the acceleration. */

	apred = acc(xpred, vpred, t+dt);

	/* Correct v and complete the step. */

	v = vpred + 0.5*(apred-a)*dt;	/* = v + 0.5*(a+apred)*dt;  */
	x = xpred;			/* NO correction to x, note */
	t += dt;
    }

    output(x, v, t);
}