Java implementation of an N-Body simulation; main loop code.

Simulation.java

// Lets go
for (int i = 0; i < bodies.size(); i++) {
	Body b = bodies.get(i);
	for (int w = i + 1; w < bodies.size(); w++) {
		Body with = bodies.get(w);
		double dist2 = // The distance between the two bodies
			with.getPosition().dist(
					b.getPosition());
		double mag =	// The force between the two bodies
				Constants.GRAVITATIONAL_CONSTANT
						* ((b.getMass() * with
								.getMass()) / dist2);
		Vector gravity = // A unit vector pointing from body "b" to body "with"
				with.getPosition()
						.clone()
						.subtract(
								b.getPosition())
						.normalize();
		// Multiply the unit vector by the force divided by the mass the force
			// is acting on to get an acceleration vector.  Add these multiplied
			// by the time change (deltaT) to the velocities.
		b.getDirection().add(
				gravity.multiply(mag * deltaT
						/ b.getMass()));
		with.getDirection().add(
				gravity.multiply(-mag * deltaT
						/ with.getMass()));
	}
}

Vector.java

package com.pi.mathz;

public class Vector {
	public double x, y, z;

	public Vector(double x, double y, double z) {
		this.x = x;
		this.y = y;
		this.z = z;
	}

	public double getAngleOn(Plane p) {
		switch (p) {
		case XY:
			return Math.atan2(x, y);
		case XZ:
			return Math.atan2(x, z);
		case YZ:
			return Math.atan2(y, z);
		}
		throw new IllegalArgumentException("Not a plane!");
	}

	public Vector(double x, double y) {
		this(x, y, 0);
	}

	public double magnitude() {
		return dist(new Vector(0, 0, 0));
	}

	public Vector translate(double x, double y, double z) {
		this.x += x;
		this.y += y;
		this.z += z;
		return this;
	}

	public Vector translate(double x, double y) {
		return translate(x, y, 0);
	}

	public Vector subtract(Vector p) {
		this.x -= p.x;
		this.y -= p.y;
		this.z -= p.z;
		return this;
	}

	public Vector add(Vector p) {
		this.x += p.x;
		this.y += p.y;
		this.z += p.z;
		return this;
	}

	public Vector multiply(double scalar) {
		this.x *= scalar;
		this.y *= scalar;
		this.z *= scalar;
		return this;
	}

	public double dist(Vector p) {
		return Math.sqrt(distSquare(p));
	}

	public double distSquare(Vector p) {
		double dX = p.x - x, dY = p.y - y, dZ = p.z - z;
		return (dX * dX) + (dY * dY) + (dZ * dZ);
	}

	@Override
	public Vector clone() {
		return new Vector(x, y, z);
	}

	@Override
	public int hashCode() {
		return (int) x << 24 ^ (int) y << 12 ^ (int) z << 6;
	}

	@Override
	public boolean equals(Object o) {
		if (o instanceof Vector) {
			Vector p = (Vector) o;
			double xD = p.x - x, yD = p.y - y, zD = p.z - z;
			return xD == 0 && yD == 0 && zD == 0;
		}
		return false;
	}

	@Override
	public String toString() {
		return "(" + x + "," + y + "," + z + ")";
	}

	public Vector translate(Vector trans) {
		return translate(trans.x, trans.y, trans.z);
	}

	public Vector normalize() {
		double dist = dist(new Vector(0, 0, 0));
		if (dist != 0) {
			x /= dist;
			y /= dist;
			z /= dist;
		}
		return this;
	}

	public static Vector normalize(Vector p) {
		return p.clone().normalize();
	}

	public static double dotProduct(Vector u, Vector v) {
		return (u.x * v.x) + (u.y * v.y) + (u.z * v.z);
	}

	public static Vector crossProduct(Vector u, Vector v) {
		return new Vector((u.y * v.z) - (u.z * v.y), (u.z * v.x)
				- (u.x * v.z), (u.x * v.y) - (u.y * v.x));
	}

	public static Vector negative(Vector p) {
		return new Vector(-p.x, -p.y, -p.z);
	}

	public Vector reverse() {
		x = -x;
		y = -y;
		z = -z;
		return this;
	}

	public Vector abs() {
		x = Math.abs(x);
		y = Math.abs(y);
		z = Math.abs(z);
		return this;
	}
}