旋转在Java物理引擎中不起作用

Question

我正在制作一个java刚体物理引擎，到目前为止它做得很好，直到我尝试实现旋转。我不知道问题出在哪里。我有使用这些网站上的公式计算凸多边形和圆的惯性矩的方法：

http://lab.polygonal.de/?p=57

http://en.wikipedia.org/wiki/List_of_moments_of_inertia

这是多边形转动惯量的代码：

public float momentOfInertia() {
    Vector C = centerOfMass().subtract(position);   //center of mass
    Line[] sides = sides();                         //sides of the polygon
    float moi = 0;                                  //moment of inertia
    for(int i = 0; i < sides.length; i++) {
        Line l = sides[i];                          //current side of polygon being looped through
        Vector p1 = C;                              //points 1, 2, and 3 are the points of the triangle
        Vector p2 = l.point1;
        Vector p3 = l.point2;
        Vector Cp = p1.add(p2).add(p3).divide(3);   //center of mass of the triangle, or C'
        float d = new Line(C, Cp).length();         //distance between center of mass
        Vector bv = p2.subtract(p1);                //vector for side b of triangle
        float b = bv.magnitude();                   //scalar for length of side b
        Vector u = bv.divide(b);                    //unit vector for side b
        Vector cv = p3.subtract(p1);                //vector for side c of triangle, only used to calculate variables a and h
        float a = cv.dot(u);                        //length of a in triangle
        Vector av = u.multiply(a);                  //vector for a in triangle
        Vector hv = cv.subtract(av);                //vector for height of triangle, or h in diagram
        float h = hv.magnitude();                   //length of height of triangle, or h in diagram
        float I = ((b*b*b*h)-(b*b*h*a)+(b*h*a*a)+(b*h*h*h))/36;     //calculate moment of inertia of individual triangle
        float M = (b*h)/2;                          //mass or area of triangle
        moi += I+M*d*d;                             //equation in sigma series of website
    }
    return moi;
}

这是圆形的：

public float momentOfInertia() {
    return (float) Math.pow(radius, 2)*area()/2;
}

我知道面积函数工作正常，我已经检查过了。我只是不知道如何检查惯性矩方程是否错误。

对于冲突检测，我使用分离轴定理来检查两个多边形和圆的任意组合，在那里它可以找出它们是否在碰撞，碰撞的法向速度，以及碰撞的接触点。这些方法都工作得很漂亮。

我可能还想说位置是如何组织的。每个物体都有一个位置和一个形状，要么是多边形，要么是圆。每个形状都有一个位置，多边形有单独的顶点。所以如果我想找到多边形形状物体顶点的绝对位置，我需要添加物体、多边形和顶点本身的位置。质心方程根据形状处于绝对位置，不考虑物体。质心和转动惯量方法在Shape类中。

对于每个物体，根据物体更新方法中的力和扭矩更新常数，其中dt是增量时间。我还根据旋转的差异旋转多边形，因为顶点是不断变化的。

public void update(float dt) {
    if(mass != 0) {
        momentum = momentum.add(force.multiply(dt));
        velocity = momentum.divide(mass);
        position = position.add(velocity.multiply(dt));
        angularMomentum += torque*dt;
        angularVelocity = angularMomentum/momentOfInertia;
        angle += angularVelocity*dt;
        shape.rotate(angularVelocity*dt);
    }
}

最后，我还有一个CollisionResolver类，它修复了两个碰撞物体的碰撞，包括施加法向力和摩擦力。这是该类唯一一个完成所有这一切的方法：

public static void resolveCollision(Body a, Body b, float dt) {
    //calculate normal vector
    Vector norm = CollisionDetector.normal(a, b);
    Vector normb = norm.multiply(-1);
    //undo overlap between bodies
    float ratio1 = a.mass/(a.mass+b.mass);
    float ratio2 = b.mass/(b.mass+a.mass);
    a.position = a.position.add(norm.multiply(ratio1));
    b.position = b.position.add(normb.multiply(ratio2));
    //calculate contact point of collision and other values needed for rotation
    Vector cp = CollisionDetector.contactPoint(a, b, norm);
    Vector c = a.shape.centerOfMass().add(a.position);
    Vector cb = b.shape.centerOfMass().add(b.position);
    Vector d = cp.subtract(c);
    Vector db = cp.subtract(cb);
    //create the normal force vector from the velocity
    Vector u = norm.unit();
    Vector ub = u.multiply(-1);
    Vector F = new Vector(0, 0);
    boolean doA = a.mass != 0;
    if(doA) {
        F = a.force;
    }else {
        F = b.force;
    }
    Vector n = new Vector(0, 0);
    Vector nb = new Vector(0, 0);
    if(doA) {
        Vector Fyp = u.multiply(F.dot(u));
        n = Fyp.multiply(-1);
        nb = Fyp;
    }else{
        Vector Fypb = ub.multiply(F.dot(ub));
        n = Fypb;
        nb = Fypb.multiply(-1);
    }
    //calculate normal force for body A
    float r = a.restitution;
    Vector v1 = a.velocity;
    Vector vy1p = u.multiply(u.dot(v1));
    Vector vx1p = v1.subtract(vy1p);
    Vector vy2p = vy1p.multiply(-r);
    Vector v2 = vy2p.add(vx1p);
    //calculate normal force for body B
    float rb = b.restitution;
    Vector v1b = b.velocity;
    Vector vy1pb = ub.multiply(ub.dot(v1b));
    Vector vx1pb = v1b.subtract(vy1pb);
    Vector vy2pb = vy1pb.multiply(-rb);
    Vector v2b = vy2pb.add(vx1pb);
    //calculate friction for body A
    float mk = (a.friction+b.friction)/2;
    Vector v = a.velocity;
    Vector vyp = u.multiply(v.dot(u));
    Vector vxp = v.subtract(vyp);
    float fk = -n.multiply(mk).magnitude();
    Vector fkv = vxp.unit().multiply(fk);                               //friction force
    Vector vr = vxp.subtract(d.multiply(a.angularVelocity));
    Vector fkvr = vr.unit().multiply(fk);                               //friction torque - indicated by r for rotation
    //calculate friction for body B
    Vector vb = b.velocity;
    Vector vypb = ub.multiply(vb.dot(ub));
    Vector vxpb = vb.subtract(vypb);
    float fkb = -nb.multiply(mk).magnitude();
    Vector fkvb = vxpb.unit().multiply(fkb);                            //friction force
    Vector vrb = vxpb.subtract(db.multiply(b.angularVelocity));
    Vector fkvrb = vrb.unit().multiply(fkb);                            //friction torque - indicated by r for rotation
    //move bodies based on calculations
    a.momentum = v2.multiply(a.mass).add(fkv.multiply(dt));
    if(a.mass != 0) {
        a.velocity = a.momentum.divide(a.mass);
        a.position = a.position.add(a.velocity.multiply(dt));
    }
    b.momentum = v2b.multiply(b.mass).add(fkvb.multiply(dt));
    if(b.mass != 0) {
        b.velocity = b.momentum.divide(b.mass);
        b.position = b.position.add(b.velocity.multiply(dt));
    }
    //apply torque to bodies
    float t = (d.cross(fkvr)+d.cross(n));
    float tb = (db.cross(fkvrb)+db.cross(nb));
    if(a.mass != 0) {
        a.angularMomentum = t*dt;
        a.angularVelocity = a.angularMomentum/a.momentOfInertia;
        a.angle += a.angularVelocity*dt;
        a.shape.rotate(a.angularVelocity*dt);
    }
    if(b.mass != 0) {
        b.angularMomentum = tb*dt;
        b.angularVelocity = b.angularMomentum/b.momentOfInertia;
        b.angle += b.angularVelocity*dt;
        b.shape.rotate(b.angularVelocity*dt);
    }
}

至于实际的问题，圆圈和多边形都旋转得很慢，而且经常朝着错误的方向旋转。我知道我在那里扔了很多东西，但是这个问题已经困扰了我一段时间，如果我能得到任何帮助，我将不胜感激。

谢谢

Answer 1

这个答案解决了“我只是不知道如何检查转动惯量方程是否错误。”问题的一部分。

有几种可能的方法，其中一些您可能已经尝试过，可以结合使用：

>

单元测试获取惯性矩代码，并将其应用于教程或教科书中已知解决方案的问题。

量纲分析无论如何，我会推荐这一点用于任何科学或工程项目。为了使发布的代码紧凑，您可能删除了注释，但它们很重要。用其单位注释表示物理量的每个变量。根据输入，检查您计算的每个表达式的结果变量的单位是否正确。例如，在以国际单位制表示的经典方程中：F以牛顿为单位，相当于千克。m/（s^2），m以千克为单位，a以m/（s^2）为单位，因此所有这些都是平衡的。注意物理世界坐标和屏幕坐标之间的转换。

程序简化尝试只处理一个非常简单形状的一个实例，您可以手动完成所有计算。由于您的一些问题与旋转无关，圆可能是一个很好的首选，因为它的对称性。调试，将中间结果与纸笔（和计算器）的等效结果进行比较。逐渐添加相同形状的更多实例，然后调试下一个形状的单个实例。。。

考虑到您怀疑惯性计算存在故意错误，请尝试设置与计算略有不同的任意值，并查看它们在显示中的差异。影响是否与您看到的问题相似？如果是这样，请将其作为假设。

一般来说，进行迭代模拟的程序很容易受到累积浮点误差的影响。除非您真的需要节省空间，并且已经对代码的数值稳定性进行了足够的分析，以确保float是可以的，否则我强烈建议使用双精度。这可能不是你当前的问题，但可能会在以后成为一个问题。