Funny, would have expected the infinite loop to occur around the curves, not the flat features. Seems to be the other way around.
Oh, happy Christmas every one. Spending it on my own with geometry this year.
[EDIT] Finally, seem to have cracked it .
For some reason in my implementation, there are certain faces on the Minkowski difference that I can't seem to get the simplex to resolve correctly when the shapes are disjoint to find the simplex closest to the origin, which is used when checking for the minimum separation vector between two shapes overlapping by only their margin.
This provides a far more accurate MSV than using sampling vectors around a unit cube, which is the fall back when the shapes overlap by greater than the margin - this is because GJK fails to find the closest simplex to the origin when the shapes intersect.
However, it seems that if I fail to compute a MSV with the GJK, I can then just fall back on the unit-sphere sampling method and this seems to return a consistent result along the face, so you still get a smooth sliding between the shapes.
The idea is to just step around a unit sphere creating normals and then project the furthest point on the Minkowski sum in that direction onto the vector. This provides the closest point on the Minkowski sum from the origin to the surface of the sum. Obviously you can only use vectors around a sphere to a fixed level of granularity although the operation is pretty cheap per vertex - assuming the vectors are precalculated, it is basically a dot product and a vector multiplication.
This will be even further improved by implementing an idea from Bullet - in addition to the unit-sphere vector checks, provide some "preferred" vectors based on the shape - these would be the negated face normals of polyhedrons for example. If the unit-sphere vectors do not directly coincide with these vectors, these vectors should "win" the test.
Huge amount of scope for improvement but all in all, I'm very happy to be sitting here on Christmas day 2012 having finally achieved something I've wanted to accomplish for about ten years - correctly detecting and resolving a collision between two convex 3D polyhedrons.
All gets a bit easier now. Trivially implemented a support function for a cylinder and a slightly rotated cube (implemented as a polyhedron). The red pyramid is controlled by the keyboard and correctly collides with and slides around the three green volumes.
The next step is to provide each volume with its own transformation but this is trivial to implement without considering optimisations, which I'm not for now. Just take the normal into the support function, inverse transform it by the shape transform, get the result then transform the result by the transform.
I'm already doing something similar with the position, adding it to the result in the support so the point is calculated in local space but transformed to world space before returning. Obviously in the case of position, this has no effect on a normal so there is no need to inverse transform it first. Once we add in rotations for the volumes, it becomes necessary but easy enough to do.
I'm going to wait till this is fully working before I start to worry about making it more efficient. At the moment it has a lot of scope for improvements in this regard but that will be fun for later.
It will be interesting to see if I can build a reasonable character controller based on a capsule using this system - something I never managed using Bullet. I'm thinking somehow that we want a capsule sitting on a line segment perhaps but that is a problem for another day.