🚀 go-pugleaf

Here we are, as promised, I present you with mssrs Henon-Heiles.
Here's the code for iterating the Henon-Heiles system, from 
CMJ fall 1992.

double	x = 0, y = 0.1
		xdot, ydot = 0.2,		/* first derivatives */
		E = 1.0/6.0,			/* max allowable energy */
		xdotdot, ydotdot,		/* second derivatives */
		xdot_squared,			/* for solution of xdot */
		delta_t = 0.01;			/* for integration */

xdot_squared = (2.0*E)-(y*y)+((4.0/3.0)*y*y*y)-(ydot*ydot);
if (xdot_squared < 0.0)	{
	printf("ERROR: xdot_squared < 0.0 : Bad initial conditions.");
	exit(1);
} 
else
	xdot = sqrt(xdot_squared);
	
for(counter = 0; counter < ITERATIONS; counter++)	{
	xdotdot = -x - (2.0 * x * y);
	ydotdot = -y - (x * x) + (y * y);
	
	xdot = xdot + (xdotdot * delta_t);
	ydot = ydot + (ydotdot * delta_t);
	x = x + (xdot * delta_t);
	y = y + (ydot * delta_t);
	
	/* output a poincare section */
	if (xdot < 0.001 && xdot > -0.001)	
		output((fabs(xdot) * 10000.0),(fabs(ydot) * 10.0));
	/* the multiplications here were done to */
	/* put the output into a more usable range */
}