stable gravity integration?

[realMakc]

I played with foam gravity demo to see how RK4 actually behaves and, while it does extremely well, it also has it's limits. Up to source mass 1e5 it works without notable error, but higher mass starts to suck planet into tighter orbits, until ~6.1651e5 where it has bifurcation point (planet gets thrown out).

This, one would think, can be ignored, but when you have 9000 particles interacting with two masses, chances for some of them to come close enough to masses are almost 100%

EDIT: I just re-wrote foam Euler integrator with solution loop that looks like this:

override public function step( dt:Number ) : void
{
			
	//get the original state
	var state:Array = _ode.state;
	//get the derivative based on that state
	_ode.getDerivative( state, derivative ); 

	// here is where things get different - get max d state[i]
	var ds:Number, max_ds:Number = -1; var i:int = -1;
	for (i = 0; i < state.length; i++)
	{
		ds = derivative[ i ] * dt;
		max_ds = Math.max (max_ds, Math.abs (ds));
	}

	// force state variables to change in amounts comparable to dt
	var N:int = Math.ceil (max_ds / Math.abs (dt));
	for (var k:int = 0; k < N; k++)
	{
		// derivative
		_ode.getDerivative( state, derivative ); 

		// advance the state by the derivative
		for (i = 0; i < state.length; i++)
		{
			state[ i ] += derivative[ i ] * dt / N;
		}
	}
			
}

Its accuracy easily outperforms RK4 at mass 5e5, but the number of sub-iterations that it takes makes player freeze when orbiter get very close to mass, so this is not feasible solution. So I keep thinking... where is newblack when one needs him?