Collision detection discrepancies

[rieder]

Related to #88:
There seems to be a discrepancy in which collisions are detected by different codes.

Expected result:
Collision between p0 and p1 is reported if p0.radius + p1.radius > (p1.position - p0.position).length(), and this may only depend on some expected roundoff/accuracy errors.
Possibly, a switch can be used to detect only approaching collisions (negative relative velocity)

Actual result:
Hermite does this perfectly
Rebound as well, but it reports p0 and p1 in reverse
ph4 misses most of the collisions, and not just the ones that are receding from each other
PhiGRAPE misses some of the collisions
BHTree misses some of the collisions

I didn't test other codes so far.

from amuse.lab import *
particles = Particles(5)#,keys=[0,1,2,3,4])
particles.mass = 1|nbody_system.mass
particles.radius = 0.6|nbody_system.length
particles.x = VectorQuantity([0,1,2,3.2,4],nbody_system.length)
particles.y = 0|nbody_system.length
particles.z = 0|nbody_system.length
particles.vx = VectorQuantity([0.1,0.01,0,0,-0.1],nbody_system.speed)
particles.vy = 0|nbody_system.speed
particles.vz = 0|nbody_system.speed

def check_collisions(particles):
    for i in range(len(particles)-1):
        p0 = particles[i]
        for j in range(i+1,len(particles)):
            p1 = particles[j]
            sum_radii = p0.radius+p1.radius
            distance = (p1.position - p0.position).length()
            if (sum_radii > distance):
                print p0.key, p1.key, sum_radii, ">", distance 

print "Manual check"
check_collisions(particles)

print "By code"
codes = [Hermite,Rebound,ph4,PhiGRAPE,BHTree]
codenames = ["Hermite","Rebound","ph4","PhiGRAPE","BHTree"]
for i in range(len(codes)):
    code = codes[i]
    print "code=%s"%(codenames[i])
    gravity = code()
    CD = gravity.stopping_conditions.collision_detection
    CD.enable()
    gravity.particles.add_particles(particles)
    gravity.evolve_model(0.01|nbody_system.time)
    print "time:",gravity.model_time
    print "code collisions found:"
    for pair in zip(CD.particles(0), CD.particles(1)):
        print pair[0].key, pair[1].key
    print "manual collisions found:"
    check_collisions(gravity.particles)
    gravity.stop()

[rieder]

Follow-up:
ph4 returns one collision, never more, probably by design. Once this one is resolved it can return another one.
However, ph4 runs into problems when particles are just touching after the collision is resolved. In that case, the timestep may become 0 and ph4 will be stuck in a loop.
A possible solution to this would be to set a tiny, but non-zero minimum timestep.

[dipto4]

Is there any way to bypass the collision detection issue for either PhiGRAPE or ph4?

[arjenve]

Yes, You can run these codes without collision detection (in fact these do not run with collision detection, it needs to be enabled)

See these lines:

CD = gravity.stopping_conditions.collision_detection
CD.enable()

The last line can be

CD.disable()

In which case no collisions will be detected.

[arjenve]

In general note that the order in which collisions are reported will never be standard, collision detection is closely coupled to the integrator in question.

PH4 can be build to support detecting all collisions, but currently it is closely coupled to the multiples module and will only detect "viable" collisions (i.e. those that are most likely to lead to a new multiples system), see https://github.com/amusecode/amuse/blob/master/src/amuse/community/ph4/src/jdata.cc#L810.

It would be interesting to determine why BHTree and PhiGrape miss some of the collisions (Stepsize and particle order are very important in the collision detection mechanism). Could you debug these? (add some prints around the collision detection code)

[rieder]

@arjenve I'll try that.

[stale]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed in 28 days if no further activity occurs. Thank you for your contributions.