generalrelativity

At OMGPOP, I’m working on a physics-based, real-time multiplayer game. Keeping a user-run simulation synchronized is a difficult (read: impossible) challenge. I’m fairly happy with the result, so I thought I’d share the idea and logic behind the networked simulation synchronization pipeline.

The model on the end user application consists of the current player and a list of all other players in the match. Each player is a subclass of a MovableElement type which defines relevant properties:

MovableElement
{
    //position
    public var x:Number;
    public var y:Number;
 
    //velocity
    public var vx:Number;
    public var vy:Number;
 
    //angular velocity
    public var av:Number;
}

Each user runs a copy of the simulation locally. Because frame rate can vary so drastically from machine to machine, it’s necessary to advance the simulation by a dynamic time-step. This introduced a problem in the collision detection scheme, as a large enough time-step would allow the character to tunnel through the floor without intersection being detected. As per the previous post, collision detection occurs between a circle and a line segment. To guarantee tunneling does not occur (without a true continuous solution, and thus a rewrite), I had to ensure that the time-step was small enough that any change in position was less than the circle’s radius:

max Δt = (r - λ) / ||v||

Where r is the circle’s radius, λ an error term to ensure we’re working with less than the radius, and v is the linear velocity of the circle in question.

Here is the main simulation loop, including solving for the maximum allowed time-step, without any opponent input:

//solve for dt
var dt:Number = synchronizedTime - lastTime;
 
//ActionScript does not throw a divide by zero error- max will be positive infinity if speed = 0
var max:Number = ( currentPlayer.radius - 1e-8 ) / currentPlayer.velocity.magnitude;
 
physicsEngine.step( dt, max );//ensures we never take a sub-step larger than max
 
lastTime += dt;

Now, when a user changes input (through key presses) and thus applies forces to their character, a notification must be sent to everyone with the current physical state of the sending player’s character; alongside a list of keys being pressed and a timestamp. Upon receipt of the character state, the corresponding player is updated on the receiver’s end, and the simulation is rewound to the time the sending player changed inputs, and then re-simulated over the latent time for just that player, according to the new set of inputs (which map to forces applied remotely to all opponents). This looks something like this:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22

private function onUserUpdatedRemotely( event:CharacterUpdateEvent ) : void
{
 
    var s:PlayerState = event.playerState;
    //find the local instance of the remote player who has updated
    var updatedOpponent:MovableElement = playerMap[ event.characterID ];
 
    //copy the state
    updatedOpponent.x = s.x;
    updatedOpponent.y = s.y;
    updatedOpponent.vx = s.vx; 
    updatedOpponent.vy = s.vy;
    updatedOpponent.av = s.av;
 
    //update forces being applied to opponent
    remoteKeyMapper.updateForces( updatedOpponent, event.keyList );
 
    //solve for change in time and advance just this opponent
    var dt:Number = synchronizedTime - event.updateTime;
    physicsEngine.stepSingleElement( updatedOpponent, dt );
 
}

And that’s pretty much it- very simple. All the hard work is writing a clean model and tight logic so that it doesn’t become a heaping pile of updates here and events there, etc.

I’m working on a game at work where a ball interacts physically with a series of capsules. The system was originally written to support circle vs line collision tests. It was then decided that it would be capsules instead of lines.

This post will show how finding the closest point on a line segment to a point can be easily extended to determine if and where a circle and line segment are intersecting; and how that test can be extended again to identify the point of intersection between a circle and capsule.

First, we need a test for the closest point, P, on a line segment, L, to point X.

We can define L as its endpoints A and B:

L = AB

We can derive the vector, v, as that which points from A to B:

v = AB = B - A

Next, we define a vector, w, pointing from A to our point, X.

w = AX = X - A

Now, we’ve expressed the relationship between X and L in terms of vectors. By defining w and v each pointing from A, we’re working in a 2D Cartesian coordinate frame where A is the origin. Because of this, we can take advantage of some geometrical properties of vectors that will allow us to find the closest point on L to X. Specifically, we’re going to project w onto v:

proj(w,v) = dot(w,v) / dot(v,v)

where dot is the dot product:

dot(w,v) = w.x * v.x + w.y * v.y

When envisioning projection, I find it helpful to think of vectors as having a spatial domain. If you compute the dot product of 2 vectors, it represents how much 1 vector exists in the domain of the other. The dot product of a vector and itself is a value which is a maximum of that vector’s domain. Anything greater than this value exceeds the domain. The dot product of a vector and one perpendicular to it is zero, and defines the minimum of that domain (exclusive). So the domain of a vector, u, in practice, ends up being spanned by any scalar value, q, which satisfies the constraint:

0 < q <= dot(u,u)

Now let’s apply this geometric intuition to our problem of finding the closest point on a line segment. What happens if the vector w points in an opposite direction than v? dot(w,v) will be negative. Therefore w cannot exist in the domain of v. If this is the case, it should be evident that A will be the closest point on L to X. Next, let’s consider the easy-to-imagine scenario of w exceeding the maximum of v’s domain, by defining w as pointing in the same direction as v, but having twice its length. Clearly, B in this case will be the closest point to X.

Example 1. Closest point on segment to user’s mouse:

0 <= t <= 1

Where t is the projection. By clamping t to this range, and using it to scale the vector v, we find the point on v closest to X. Because we’re working in the space of A, we need to transform back to world coordinates, by adding A:

P = A + v * clamp( proj(w,v), 0, 1 )

Here’s the implementation in my geometry library (Vector2D is not a native data type), with the variables changed to match the description of the problem:

public static function closestPointLineSegment( X:Vector2D, A:Vector2D, B:Vector2D ) : Vector2D
{
    var v:Vector2D = B.minus( A );
    var w:Vector2D = X.minus( A );
    var wDotv:Number = w.dot( v );
    var t:Number = w.dot( v ) / v.dot( v );
    t = MathUtil.clamp( 0, 1, t );
    return A.plus( v.times( t ) );
}

Circle/Line Collision Test

Now that we’ve found the closest point on a line segment, L, to a point, X, determining intersection of a circle and line segment is as simple as checking whether the distance between X and P is less than the radius of the circle. Of course, X defines the center of the circle.

Example 2. Intersecting a circle and line segment:

Note that in the code below, I’m checking for distance by comparing the squared distance against the squared length of the circle’s radius. This is to avoid an expensive square root call.

public static function testCircleSegment( X:Vector2D, r:Number, A:Vector2D, B:Vector2D ) : Vector2D
{
    var P:Vector2D = closestPointLineSegment( X, A, B );
    P.minusEquals( X );
    if( P.dot( P ) <= r * r ) return P;
    return null;
}

Circle/Capsule Collision Test

Lastly, intersecting a circle and a capsule is as simple as inflating the circle’s radius by the radius of the capsule. You can see in the below example that we’re really just doing a circle/line test on an inflated circle.

Example 3. Intersecting a circle and capsule:

One extremely nice aspect of working with circles is that, not allowing the shapes to penetrate, circle’s can be touching another convex shape at a maximum of 1 point. Due in part to this convenience, the collision normal will always point from the closest point on the convex shape through the circle’s center.

public static function testCircleCapsule( X:Vector2D, r1:Number, A:Vector2D, B:Vector2D, r2:Number ) : Vector2D
{
    var P:Vector2D = closestPointLineSegment( X, A, B );
    P.minusEquals( X );
    if( P.dot( P ) <= ( r1 + r2 ) * ( r1 + r2 ) ) return P;
    return null;
}

Here’s a demo of the auto-gen physics engine (SLRBS- Super Lightweight Rigid Body Simulator) I’ve been working on which allows you to select images to use as rigid bodies. First, the image is analyzed, ignoring transparent pixels, to determine its convex hull. Then mass/center of mass/inertia tensor are solved for. Finally, the rigid body is created and dropped into the world.

I’m sure that in playing with it, you’ll notice all sorts of wonky behavior… but I think it’s fun anyway. Let me know if there’s a class or type of shape/image that tends to be particularly misbehaved!

demo
image pack

This was something I did a while ago and never got around to posting about. It’s a simulation of an incompressible fluid- when you render the density and apply temperature by adding an upward force, it looks like smoke!

It’s fairly computationally intensive. The state of the system has to be solved simultaneously. But because it uses an iterative solver to approximate a solution, it’s easy to tune for performance to realism.

Click and drag to knock the smoke around, release to add density (more smoke).

My implementation is based on this paper by Joe Stam. The leap from the Navier Stokes equations to this compact form is pretty radical. Mr. Stam’s a pretty smart guy.

For the past few months I’ve been working on a general purpose geometry library for representing your interactive elements in a more precise geometrical way. There are some pretty cool features like collision detection management, auto-generation of geometry, and an easy way to bind the model to your DisplayObject hierarchy (or not!). 

It will be released under the MIT license alongside SLRBS- a physics engine that uses the library. The absolute last thing the Flash world needs is another physics engine, but this has some really nice features for beginner developers, prototypers etc. And it does a great job of showing how to use the geometry library.

One of the most frequent questions or requests I get/got with FOAM was something along the lines of “how do I make my MovieClip a rigid body?” It can be hard from my position to know how to tell someone who probably uses the Flash authoring tool a proper way to handle the difference between the model and view. So with SLRBS, I set out to address this issue.

With SLRBS, you have a typical API to create physics simulation, but it also has some streamlined auto-creation methods that make setting it up a breeze. To ensure this, I even downloaded a copy of Flash CS4 to play with it!!! AAAAHH!!! wtf happened to Flash!? 

So, all you do is lay out your MovieClips on the stage, and call a single static method and voila! physics simulation! There’s some pretty cool stuff going on that determine how to handle overlapping DisplayObjects (each demo does it differently), and it’s just a single parameter.

Here’s a screenshot of the stage in a demo where a car is being created:

Car Demo

Because 2 of the tires overlap the car, they get pinned to the car at their center of mass.

 

This is the only ActionScript in the .fla:

1
2
3
4
5
6
7
8
9

import org.generalrelativity.slrbs.SLRBSimulator;
 
var simulator:SLRBSimulator = SLRBSimulator.createSLRBSimulatorFromContainer( this );
addEventListener( Event.ENTER_FRAME, onEnterFrame );
 
function onEnterFrame( event:Event ) : void
{
   simulator.step( 1 / 60 );
}

 

Here’s a second demo’s stage in the Flash authoring tool- this time a ragdoll. The overlap this time gets pinned at the center of overlap instead of the center of mass as in the car demo- this works well for creating ragdoll joints.

Ragdoll Demo

This is a small framework I wrote that addresses some of the issues associated with ActionScript being executed linearly in one thread. If the language offered an API for threading, computationally expensive routines could be run in a separate thread, keeping the application’s UI appropriately responsive.

The idea is to offer a simple means to decompose:

• deep recursion
• large loops
• computationally expensive tasks

I’ve created 2 demos (FP10) to illustrate the issue.

The first demo is an implementation of the A* algorithm. If you select to execute the path-finding normally, you’ll notice the Flex UI becomes unresponsive; if you’re lucky, you might even get the spinning beach ball of death. But when you run the algorithm in a green thread, the application UI remains responsive. Select start and end nodes at diagonally opposite corners to force the algorithm through a lot of work.

The second demo encodes 3 images. This is a little harder to visualize in a demo, so here’s what’s happening: I take an embedded image, encode it as a .png 3 times and then decode and render the ByteArray to show the encoding’s success. As with the path-finding demo, you’ll see the UI become unresponsive when you run the encoding normally. But decomposed to run over a few seconds, everything runs well.

A thread’s chief responsibility, aside from running its processes, is to maintain frame rate. The API offers the user the ability to set execution frequency and processing share (percentage). It polices itself according to these allocation constraints. If a process exceeds its allotment, it is penalized in the next cycle. This works surprisingly well to maintain frame rate and adhere to the user’s request.

The source is pretty well commented and should be easy to implement! And if nothing else, there’s an object-oriented implementation of the A* algorithm… (based on Ian Millington’s pseudo code in his Artificial Intelligence For Games).

path-finding demo

image encoding demo

source

I started this framework because of (and got the name from) my awesome friend, Roger Braunstein, who shared his elegant and simple solution with me a couple months ago, and referenced the obscure deprecated Java API.

The next release of the Flash Player has support for 2.5D- which boils down to 3D properties, appropriate perspective, but no depth sorting. This method takes a DisplayObjectContainer, transforms its children into world space and sorts them along that z-axis. It’s sloppy because it’s DisplayObject-based and so there will be noticeable “pop” when a DisplayObject is sorted onto the top of the display list in close proximity to another. That said it should be perfect for carousels and other views that aren’t super tight.

Here’s a demo

And here’s the algorithm:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78

public static function simpleZSort3DChildren( doc:DisplayObjectContainer, recurse:Boolean = true ) : void
{
 
    //transforms from local to world oordinate frame
    var transform:Matrix3D = doc.transform.getRelativeMatrix3D( doc.stage );
 
    var numChildren:int = doc.numChildren;
 
    //v = ( n * 3 )- (x,y,z) set for each child
    var vLength:int = numChildren * 3;
 
    var vLocal:Vector.<Number> = new Vector.<Number>( vLength, true );
    var vWorld:Vector.<Number> = new Vector.<Number>( vLength, true );
 
    //insertion point for child’s coordinates into state vector
    var vIndex:int = 0;
 
    for( var i:int = 0; i <numChildren; i++ )
    {
 
        var child:DisplayObject = doc.getChildAt( i );
        if( recurse && child is DisplayObjectContainer ) simpleZSort3DChildren( DisplayObjectContainer( child ), true );
 
        vLocal[ vIndex ] = child.x;
        vLocal[ vIndex + 1 ] = child.y;
        vLocal[ vIndex + 2 ] = child.z;
 
        vIndex += 3;
 
    }
 
    //populates vWorld with children coordinates in world space
    transform.transformVectors( vLocal, vWorld );
 
 
 
    //bubble sorts children along world z-axis
    for( i = numChildren - 1; i > 0; i-- )
    {
 
        var hasSwapped:Boolean = false;
 
        vIndex = 2;
 
        for( var j:int = 0; j < i; j++ )
        {
 
            //z value at that index for each child
            var z1:Number = vWorld[ vIndex ];
 
            vIndex += 3;
 
            var z2:Number = vWorld[ vIndex ];
 
            if( z2> z1 )
            {
 
                //swap
                doc.swapChildrenAt( j, j + 1 );
 
                //swap z values (don’t need to change x and y because they’re not used anymore)
                vWorld[ vIndex - 3 ] = z2;
                vWorld[ vIndex ] = z1;
 
                //mark as swapped
                hasSwapped = true;
 
            }
 
        }
 
        //if there was no swap, we don’t need to iterate again
        if( !hasSwapped ) return;
 
    }
 
 
}

This is the first in a series of quick tutorials on helpful maths, physics and algorithms, and their ActionScript implementations. Today we’re going to learn how to shoot a projectile out of the sky. Sweet.

DEMO

First, let’s define the algorithm’s structure. We’re going to pass the time to collision. The algorithm will return the muzzle velocity that our projectile requires to hit a moving target out of the air after this specified amount of time.

function solveMuzzleVelocity( timeToImpact:Number ) : Vector

Next, we have to define the equations of motion:

x1 is a particle’s position after t seconds. x0 is the original position and the dots represent the first and second positional derivatives (velocity and acceleration, respectively). This equation allows us to solve for position after any specified time (given constant acceleration). Referring to our algorithm’s structure, we know t (timeToImpact) to be defined. This means that we can determine where the target will be when we hit it, which gives us a point of impact:

//point of impact
var poi:Vector = target.position.plus( target.velocity.times( timeToImpact ) );
poi.plusEquals( new Vector( 0, Projectile.G * timeToImpact * timeToImpact * 0.5 ) );

It’s assumed that the only accelerative force is gravity (Projectile.G here), that we know the starting position of the target and our projectile, and the current velocity of the target.
Plugging this point of impact in for our projectile’s equation of motion (as x1), we eliminate all unknowns except the muzzle velocity. Through basic algebraic manipulation, we can isolate and solve for this unknown:

Finally, divide each side of the equation by t:

The algorithm will return this value. Here’s my implementation, in full:

1
2
3
4
5
6
7
8
9
10
11
12
13

private function solveMuzzleVelocity( timeToImpact:Number ) : Vector
{
 
    //point of impact
    var poi:Vector = target.position.plus( target.velocity.times( timeToImpact ) );
    poi.plusEquals( new Vector( 0, Projectile.G * timeToImpact * timeToImpact * 0.5 ) );
 
    var diff:Vector = poi.minus( projectile.position );
    diff.y -= Projectile.G * timeToImpact * timeToImpact * 0.5;
 
    return diff.dividedBy( timeToImpact );
 
}

Super simple! Now you can integrate dynamically and knock some shit out of the sky.

Here’s a verbose implementation that doesn’t use my Vector class and isn’t dependent on my application:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15

private function solveMuzzleVelocity( t:Number, xp:Point, xt:Point, vt:Point, pixelToMeterRatio:Number = 40 ) : Point
{
 
    const g:Number = 9.81;
    var vDotTerm:Number = g * pixelToMeterRatio * t * t * 0.5;
 
    var poiX:Number = xt.x + vt.x * t;
    var poiY:Number = xt.y + vt.y * t + vDotTerm;
 
    var diffX:Number = poiX - xp.x;
    var diffY:Number = poiY - xp.y - vDotTerm;
 
    return new Point( diffX / t, diffY / t );
 
}

t is time to impact, xp holds the x and y coordinates of the projectile’s starting position, xt is the target’s starting position, vt is the target’s starting velocity and pixelToMeterRatio defines how many pixels represent one meter.

I’m sorry that this took so long- especially those of you who’ve written me and patiently waited for my help with rendering… I simply don’t have as much time as I would like. Anyway- this will offer a solution for rendering your FOAM physics simulation with:

• Drawing API (SimpleFoamRenderer)
• BitmapData (BitmapFoamRenderer)
• DisplayObjects (DisplayObjectFoamRenderer)

Here’s a very silly demo- Mouse down anywhere to switch from the DisplayObject renderer to the Bitmap renderer.

demo
source

Note that these new renderers are meant as example only- I do think in most cases the developer will need a more specific setup than I could hope to generically offer. So consider this my implementation and at best a nudge in the right direction.

Hopefully both of these will be fairly self-explanatory, but I’ll briefly go over the DisplayObject renderer and what my thinking was when addressing this issue.

I wanted a means to render a DisplayObject (base class of Sprite, MovieClip etc.) that’s the visual representation of a RigidBody being simulated. I knew that I’d like to be able to pass either the DisplayObject’s Class or the instance itself. I also knew that I’d need an easy way to adjust the position of the DisplayObject- so first, I wrote a datatype to hold these properties- DisplayObjectData.

You’ll see that it has a setter for defining its displayObject property. This is what facilitates passing either a DisplayObject Class or a DisplayObject instance. The other properties:

• offsetX – the amount to offset the DisplayObject from its origin in the x direction
• offsetY – the amount to offset the DisplayObject from its origin in the y direction
• autoCenter – Boolean indicating whether to center the contents about the origin- Because of the way a RigidBody rotates about its center of mass, I figured this would be typical.
• hasBeenDisplayed – a helper Boolean for determining whether to perform the indicated transformations

In Rocket.as, you can see the way I’m embedding my visual assets:

[Embed(source="../../../bin/images/rocket/lunarPod.png")]
private var rocketSpriteClass:Class;
 
[Embed(source="../../../bin/images/rocket/asteroid.png")]
private var asteroidSpriteClass:Class;

If this is unfamiliar with you, simply note that it’s analogous to the classes you’re instantiating as MovieClips in Flash CS3. As you’ll see in a moment, it’s not necessary to pass the Class, we can pass the DisplayObject instance as well.

Here’s how I’m passing the class to FOAM to tie it to to a RigidBody. Note that the instance of DisplayObjectData simply defines Foam.addElement’s renderData argument.

foam.addElement( rocket, true, true, new DisplayObjectData( rocketSpriteClass, 0, -5 ) );

Compare this with how I’m creating the asteroid (I wanted to scale the asteroid, or else I’d have created it in the same way as the rocket- it worked out well for the purpose of demonstration):

var asteroidSprite:DisplayObject = new asteroidSpriteClass();
asteroidSprite.scaleX = asteroidSprite.scaleY = size / 64;
var asteroid:Circle = new Circle( sx, sy, size, size * 4, svx, svy, 0.5, 0.25, 0, -0.1 + Math.random() * 0.2 );
foam.addElement( asteroid, true, true, new DisplayObjectData( asteroidSprite ) );

There’s one last important thing to note in the Rocket example class. Because I don’t want to use the default FOAM renderer, I have to tell it which renderer to use:

foam.setRenderer( new DisplayObjectFoamRenderer() );

Now let’s move on to DisplayObjectFoamRenderer. Note that I’ve chosen to extend SimpleFoamRenderer so that I don’t have to rewrite drawing routines for objects I don’t have graphics for (such as the mouse dragger bungee). You’ll see my use of the DisplayObjectData properties in the addRenderable method.

So- ultimately all of this boils down to 3 lines of code which line up our DisplayObject with our RigidBody:

displayObject.x = ISimulatable( renderable.element ).x;
displayObject.y = ISimulatable( renderable.element ).y;
if( renderable.element is IBody ) displayObject.rotation = IBody( renderable.element ).q * 180 / Math.PI;

It places the DisplayObject at the x and y coordinates of the RigidBody’s center of mass. Because rotation is given in radians, we have to convert to an angle which is what DisplayObject.rotation expects.

And that’s pretty much it; for the most part a FOAM renderer is bookkeeping. It might suit you far better to just maintain reference to your DisplayObject and update it every frame just as I am in the above 3 lines of code… simply set render to false in Foam.addElement.

This article walks through my implementation of modeling orbit with FOAM. The content lays a groundwork which could be used to simulate (a ridiculously simplified version of) our solar system for instance. It has 2 main goals:

1. Give the reader a solid and quick understanding of how one might use FOAM specifically.
2. Exemplify the difference between differential equation solvers.

Orbit Demo
FOAM source

When wanting to create physically realistic animation using FOAM, the first step will often be deciding how to develop the forces required to create the desired motion. Forces in FOAM are primarily handled via modular implementations of IForceGenerator- which simply dictates that the method generate generates a force for its passed element. Gravity in the context of our atmosphere and friction are 2 force generators included in the source.

But how would we go about modeling the force an extremely massive object exerts on another?

Firstly, we need to define an equation which ultimately yields a force that we can apply. Newtonian gravity offers the simple force equation:

• F is the force we’re solving for
• G is a gravitational constant
• m1/m2 are the masses of the bodies in question
• r is the distance between the bodies

It should be clear from looking at this equation that the larger the mass, the shorter the distance, the greater the force- which is about as complex as we need to go to get realistic looking results. One simplification we’ll make is assuming that the body exerting the gravitational force is so much larger than the body being acted on, that we can neglect the force acting back on it.

This is simple enough, but we need to translate Newton’s equation into a form our physics engine can execute- We need to create a new force generator. We’ll look at what’s going on in its constructor and generate method; the class in its entirety can be found here.

First, the constructor takes the element we want to exert the force and gravitationalConstant- this is a magic number that should be tweaked per simulation:

1
2
3
4
5

public function GravitationalForceGenerator( source:ISimulatable, gravitationalConstant:Number = 1.2 )
{
   this.source = source;
   g = gravitationalConstant;
}

I’ve left out some error checking for the sake of brevity here- reference the complete class. Next we implement generate:

1
2
3
4
5
6
7

public function generate( element:ISimulatable ) : void
{
   //find the difference vector
   var diff:Vector = source.position.minus( element.position );
   //add our solved force along the difference vector to the supplied element
   element.addForce( diff.getUnit().times( g * source.mass * element.mass / diff.dot( diff ) ) );
}

Line 6 contains our equation for universal gravitation as discussed above. The change worth noting here is that the force is applied in the direction of the exerting body (GravitationalForceGenerator.source). We do this by normalizing the bodies’ difference vector, and scaling it by the force found in Newton’s equation. Note also that the dot product of a vector and itself is its squared magnitude- this is how we’re finding the square of the bodies’ distance- the equation’s right-hand side denominator.

Now let’s look at how I’ve used this force generator in the context of a FOAM application. Remember that we’re also going to examine the difference between IODESolvers.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21

//create a source of the gravitational pull
var source:Circle = new Circle( 400, 300, 60, 10000 );
 
//create an element influenced by source's gravitational pull to be integrated with the Euler IODESolver
var eulerOrbital:Circle = new Circle( 700, 300, 30, 100, 0, -7 );
//create an element influenced by source's gravitational pull to be integrated with the RK4 IODESolver
var rungaKuttaOrbital:Circle = new Circle( 700, 300, 30, 100, 0, -7 );
 
//create the gravitational force generator
var gravity:GravitationalForceGenerator = new GravitationalForceGenerator( source );
 
//add this to each orbital
eulerOrbital.addForceGenerator( gravity );
rungaKuttaOrbital.addForceGenerator( gravity );
 
//add each element to FOAM (no collisions)
foam.addElement( source, false, true, { color:SimpleOrbit.SOURCE_COLOR, alpha:1 } );
//specify the use of our Euler solver
foam.addElement( eulerOrbital, false, true, { color:SimpleOrbit.EULER_COLOR }, new Euler( eulerOrbital ) );
//RK4 is default, no need to explicitly set
foam.addElement( rungaKuttaOrbital, false, true, { color:SimpleOrbit.RK4_COLOR } );

The important aspects here are:

• we create a new GravitationalForceGenerator and add it as a force generator to each of our orbiting bodies, setting our source Circle as its source.
• the orbitals are created with identical masses, positions and initial velocities (tangent to desired orbit at location)
• one orbital gets added with FOAM’s default solver, RK4- the other with an explicitly defined Euler solver.

The whole SimpleOrbit class can be found here.

If you take a look at the demo, you’ll notice that slowly but surely the orbital being integrated with our Euler solver drifts out of desirable orbit. It will continue to spiral outward. If we set the number of iterations we allow the engine to use per frame to 1, you would see that this error accumulates faster. If instead, we set the iterations to 100 for instance, its orbit would match the more precise RK4-integrated orbital much longer.

So why is this happening? To better explain, I’ll steal some images and description from an earlier post of mine- reference the Euler and RK4 classes also.

At every step in time in which we evaluate state and derivative, the forces acting on our orbiting bodies is different. If we were able to do this evaluation an infinite number of times every frame, Euler would provide the correct answer. This is obviously impossible and so even though the forces acting on the body continually change, we’re limited to apply the force found at a single moment in time across the entire interval over which we’re integrating.

The line segment piercing the starting point represents the derivative. Because Euler simply advances the equation via this derivative, it should be easy to see why it is so error prone. Compare this with how RK4 works:

Because it uses the slope at multiple points over the interval to evaluate the ultimate derivative, our approximation using RK4 is entirely more accurate. This accuracy comes at a cost- RK4 is very much more expensive to calculate- although if you tried the different iteration suggestion I made, you’d see it is more accurate at 3 iterations than Euler is at 100. This is one of FOAM’s strong points- it lets you choose the solver best suited for a specific task within a simulation.

Hopefully this has been helpful in understanding both how to solve a problem with FOAM and what the different integrators are useful for. Future releases will include many more solvers.

Sorry for the super lame demo graphics- another early article will detail the implementation of a decent FOAM renderer.