Bug with AccelerationBased coupled linear axes motor in 2D

[Schnorg]

There seems to be a problem with the AccelerationBased motor with coupled linear axis. I don't know if it is just a coupled axes problem, but it would be hard to observe errant behaviour on a single axes.

This is using ForceBased
https://github.com/user-attachments/assets/baf9ac51-f7a2-4e66-96ed-23087f0250a9

This is using AccelerationBased (stiffness and damping doubled to account for 2x 1Kg mass on each side)
https://github.com/user-attachments/assets/4d7fa4b6-dc69-40d8-bd02-39802eb5b4ec

I have benchmarked this against a modified version of Godot 2D Physics (where I have added acceleration based and force based controls to the PinJoint2D locked axes). The Godot version acts normally in acceleration mode and similarly to Rapier in force mode.

Further notes on the setup

• The pin joints are all connected to the center of one body (which is offset for the other)
• Each body is the default 1Kg.
• $k (accel mode) = 2.0 * k (Force mode)$ (and similar for $c$), this is to account for resting $M_\text{eff}$ with $\frac{1}{M_\text{eff}} = \frac{1}{m1} + \frac{1}{m2}$.
• I tried different $c$ settings in case damping / stiffness balance was a problem. It didn't fix it.
• I tried it at 240fps and it behaved the same.

See extract of my joint setup code below to indicate settings (from Godot-Rapier-Physics). The joints 2D linear axes are being unlocked, coupled and given a motor that is calculated from Godot softness settings. It works in force mode. See the part in the 'if (softness > 0.0) section'.

            joint.set_motor_model(JointAxis::AngX, MotorModel::AccelerationBased);
            joint.set_motor_max_force(JointAxis::AngX, Real::MAX);
            if angular_limit_enabled {
                joint.set_limits(JointAxis::AngX, [angular_limit_lower, angular_limit_upper]);
            } else {
                joint.limit_axes.remove(JointAxesMask::ANG_X);
            }
            if motor_enabled {
                joint.set_motor(JointAxis::AngX, 0.0, motor_target_velocity, 0.0, 0.0);
                #[cfg(feature = "dim3")]
                joint.set_motor_max_force(JointAxis::AngX, motor_max_force);
            } else if motor_position_enabled {
                joint.set_motor(
                    JointAxis::AngX,
                    motor_target_position,
                    0.0,
                    motor_stiffness,
                    motor_damping,
                );
            } else {
                joint.motor_axes.remove(JointAxesMask::ANG_X);
            }
            if softness > 0.0 {
                let force_based: bool = true;
                // Emulate softness by using ForceBased motor on LinX (coupled with LinY)
                joint.locked_axes = JointAxesMask::FREE_FIXED_AXES;
                joint.coupled_axes = JointAxesMask::LIN_AXES;
                if force_based {
                    joint.set_motor_model(JointAxis::LinX, MotorModel::ForceBased);
                    //joint.set_motor_model(JointAxis::LinY, MotorModel::ForceBased);
                } else {
                    joint.set_motor_model(JointAxis::LinX, MotorModel::AccelerationBased);
                    //joint.set_motor_model(JointAxis::LinY, MotorModel::AccelerationBased);
                }
                joint.set_motor_max_force(JointAxis::LinX, Real::MAX);
                //joint.set_motor_max_force(JointAxis::LinY, Real::MAX);
                // See Catto, Soft Constraints 2011 slides for CFM and ERP definitions
                // softness = CFM / h, where h = time delta
                // bias = ERP
                // Catto 2011 shows CFM and ERP as they relate to k and c
                // plug them in, rearrange, and you get the equations below to
                // figure out k and c. Only works with a ForceBased spring.
                let h: Real = physics_step;
                let mut inv_mass_eff = 1.0 / 1.0 + 1.0 / 1.0;
                if force_based {
                    inv_mass_eff = 1.0;
                }
                // inv_mass_eff *= 0.1;
                let k: Real = inv_mass_eff * bias / (h * h * softness);
                let c: Real = inv_mass_eff * (1.0 - bias) / (h * softness);
                joint.set_motor(JointAxis::LinX, 0.0, 0.0, k, c);
                //joint.set_motor(JointAxis::LinY, 0.0, 0.0, k, c);
            } else {
                joint.locked_axes = JointAxesMask::LIN_AXES;
                joint.coupled_axes = JointAxesMask::FREE_FIXED_AXES;
                joint.motor_axes.remove(JointAxesMask::LIN_X);
            }

Not sure if this is helpful in any way, but In the case where there are 2 degrees of freedom for the spring, in acceleration mode, wouldn't cfm_gain need to be a 2x2 matrix? Unless the impulse equation is rearranged for acceleration mode (to cancel out effective mass terms and eliminate the need for a matrix). Again, not sure if it's helpful, but that's where the math seemed to be leading me when doing a similar thing in Godot (still learning though).