This isn't about understanding angular physics, but more how to actually implement it.
For example, I'm using a simple linear interpolation per frame (with dT) I'm having trouble with the angular units, I want to be able to rotate around arbitrary axes.
(with glm)
Using a vec3 for torque, inertia and angular velocity works excellent for a single axis.
Any more and you get gimbal lock. (i.e. You can rotate around a local x, y or z but superimposing prevents proper results)
Using quaternions I can't get it to function nicely with time, inertia or for an extended period.
Is there any tried-and-true method for representing these features?
The usual solution is to use the matrix representation of rotation. Two rotations in sequence can be expressed by multiplying their respective matrices. And because matrix multiplication is not symmetric, the order of the 2 rotations matters - as it should.