摩擦冲量已经落在圆盘内，为何仍需切向质量块求解

min F(J) = 0.5 J^T K J + v_free^T J
subject to |J| <= R

(K + lambda I) J = -v_free
|J| = R

public readonly record struct Vec2(double X, double Y);
public readonly record struct TangentMass(double A, double B, double C);

public static class FrictionBlock
{
    public static Vec2 Solve(TangentMass k, Vec2 freeVelocity, double radius)
    {
        if (!double.IsFinite(k.A) || !double.IsFinite(k.B)
            || !double.IsFinite(k.C) || k.A <= 0 || k.C <= 0
            || !double.IsFinite(k.A*k.C-k.B*k.B) || k.A*k.C-k.B*k.B <= 0)
            throw new ArgumentException("有效逆质量矩阵必须正定");
        if (!double.IsFinite(radius) || radius < 0
            || !double.IsFinite(freeVelocity.X) || !double.IsFinite(freeVelocity.Y))
            throw new ArgumentException("速度与半径必须为有效有限值");
        if (radius == 0) return new(0, 0);
        Vec2 unrestricted = At(k, freeVelocity, 0);
        if (Norm(unrestricted) <= radius) return unrestricted;
        double low = 0, high = 1;
        while (Norm(At(k, freeVelocity, high)) > radius)
        {
            high *= 2;
            if (!double.IsFinite(high)) throw new ArithmeticException("区间扩展溢出");
        }
        // high 始终位于可行侧，返回值不因停止条件越出圆盘。
        for (int i = 0; i < 80; i++)
        {
            double middle = low + (high-low)*.5;
            if (Norm(At(k, freeVelocity, middle)) > radius) low = middle;
            else high = middle;
        }
        return At(k, freeVelocity, high);
    }

    private static Vec2 At(TangentMass k, Vec2 v, double lambda)
    {
        double a = k.A+lambda, c = k.C+lambda;
        double determinant = a*c-k.B*k.B;
        if (!double.IsFinite(determinant) || determinant <= 0)
            throw new ArithmeticException("矩阵求解超出数值范围");
        Vec2 result = new((-c*v.X+k.B*v.Y)/determinant,
            (k.B*v.X-a*v.Y)/determinant);
        if (!double.IsFinite(result.X) || !double.IsFinite(result.Y))
            throw new ArithmeticException("冲量超出数值范围");
        return result;
    }

    public static double Norm(Vec2 v) => Math.Sqrt(v.X*v.X+v.Y*v.Y);
}

Vec2 freeVelocity = new(
    currentVelocity.X - k.A*oldImpulse.X - k.B*oldImpulse.Y,
    currentVelocity.Y - k.B*oldImpulse.X - k.C*oldImpulse.Y);
Vec2 next = FrictionBlock.Solve(k, freeVelocity, radius);
Vec2 delta = new(next.X-oldImpulse.X, next.Y-oldImpulse.Y);
// 调用方按接触雅可比将 delta 分别施加到两个刚体，并保存 next。

dotnet build $WORK_DIR/FrictionLab.csproj -c Release
dotnet run --project $WORK_DIR/FrictionLab.csproj -c Release --no-build

build_exit=0 warnings=0 errors=0 build_elapsed_s=5.68
run_exit=0 passed=21 failed=0 skipped=0
samples=10000 boundary_angles=256 elapsed_ms=301.2997
max_budget_excess=0
max_kkt_residual=3.4684476073050936e-15
max_rotation_error=1.5710734600032811e-15
max_dense_objective_gap=0