I just tried implementing volume preservation, it doesn't seem to help. Things are less rigid, but still inside-out-able. Since I desired rigid behaviour, I'm forced to use beta = 0 anyway..
This is quite frustrating! Have you tried simulating an almost-planar tetrahedron? It might be difficult to tell (visually) if it has flipped inside out or not.
I feel that the solution doesn't lie in the matrix-fitting steps, since those steps are working on an inverted solution. Adjusting Apq so that the determinant is positive seems to be the correct thing to do, however I don't know how to approach this.
My only thought so far has been, since I have a method which works when the principle axis/axis-of-greatest-spread is axis aligned, I could fit the particles with a line to determine the current principle axis, rotate everything to align the principle axis to one of the world axes, correct Apq (i.e negate one of the columns, depending on which world axis I've aligned everything to), then rotate things back. However, this seems a bit ridiculous.