Skip to content

Commit

Permalink
Fix locked axes in gyro torque (#486)
Browse files Browse the repository at this point in the history
# Objective

I think that #485 wasn't quite the right fix for #474 - applying locked axes to an inertia tensor zeros out one or more of its rows, making the matrix singular - taking the inverse then results in a matrix that's all NaNs, making the delta_ang_vel also all NaNs, which causes the angular momentum to be unaffected by anything that timestep (including applied torques!).

The problem with the status quo ante seems to be that a similar thing happened - rather than trying to invert the (singular) locked inverse inertia, the calculated Jacobian was (near-)singular in a lot of cases when axes were locked, and depending on round-off errors you'd either get NaNs or just an outright cancellation of the angular velocity

## Solution

Calculating the delta-ω from gyroscopic effects with an unconstrained inertia matrix and *then* applying the locked axes to that angular velocity increment should be more reliable, and should avoid undesirable NaN values showing up
  • Loading branch information
unpairedbracket authored Aug 19, 2024
1 parent 4d082a7 commit f2dec36
Showing 1 changed file with 3 additions and 7 deletions.
10 changes: 3 additions & 7 deletions src/dynamics/integrator/semi_implicit_euler.rs
Original file line number Diff line number Diff line change
Expand Up @@ -77,13 +77,9 @@ pub fn integrate_velocity(
// However, the basic semi-implicit approach can blow up, as semi-implicit Euler
// extrapolates velocity and the gyroscopic torque is quadratic in the angular velocity.
// Thus, we use implicit Euler, which is much more accurate and stable, although slightly more expensive.
let effective_inertia = locked_axes.apply_to_rotation(inv_inertia.0).inverse();
delta_ang_vel += solve_gyroscopic_torque(
*ang_vel,
rotation.0,
Inertia(effective_inertia),
delta_seconds,
);
let delta_ang_vel_gyro =
solve_gyroscopic_torque(*ang_vel, rotation.0, inv_inertia.inverse(), delta_seconds);
delta_ang_vel += locked_axes.apply_to_angular_velocity(delta_ang_vel_gyro);
}

if delta_ang_vel != AngularVelocity::ZERO.0 && delta_ang_vel.is_finite() {
Expand Down

0 comments on commit f2dec36

Please sign in to comment.