Accurately calculating gravity for large systems
Gravity for small distances

I've seen two simple ways that are used to fix this problem:


Add a small amount to each distance:

double d2 = d2(p);

d2 += EPSILON_DIST;

Use a minimum distance:

double d2 = d2(p);

if (d2 &lt; MIN_DIST)
    d2 = MIN_DIST;



Either way will prevent small distances from causing huge accelerations.  Since d2 is in units of distance squared, your EPSILON_DIST and MIN_DIST should also be squared distances.

More efficient gravity computation

Your gravity calculation could be improved a little.  I answered another almost identical question here, so you can take a look at the solution I came up with there.  The main point is that you can reduce the number of divisions to only one in grav() and zero in move().  Although it doesn't exactly match your Particle class, I think you can extrapolate from it.  Here is the code, pasted from that other question:

public void interact(Body other) {
    double dx = other.getX() - x;
    double dy = other.getY() - y;
    double r  = calculateDistance(dx, dy);
    double inv_r3 = 1.0 / (r * r * r);

    /* Precalculate force component (1/r^2) times direction (dx/r) = dx / r^3 */
    dx *= inv_r3;
    dy *= inv_r3;

    /* calculate accelerations for both bodies */
    vx += other.getMass() * dx;
    vy += other.getMass() * dy;
    other.addToVx(mass * -dx);
    other.addToVy(mass * -dy);
}

public void update() {
    x += vx;
    y += vy;
}