轨道力学

Question

有没有人有一个实施轨道力学的例子(最好是在XNA中)？我目前使用的代码如下,但它执行时并没有"感觉正确".物体只是稍微向地球弯曲,无论我调整多少变量,我都无法进入轨道,甚至是部分轨道.

shot.Position += shot.Velocity;  

foreach (Sprite planet in planets)  
{  
  Vector2 directionToPlanet = (planet.Position - shot.Position);  
  directionToPlanet.Normalize();  

  float distance = Vector2.DistanceSquared(shot.Position, planet.Position);  

  float gPull = (float)(planet.gravityStrength * (planet.Mass * shot.Mass) / distance) + planet.gravityField;  
  shot.Position += new Vector2(directionToPlanet.X * gPull, directionToPlanet.Y * gPull);  
} 

编辑 标记Mendelt的答案是正确的,指出我需要更新速度,而不是位置.我还需要将gPull的计算更改为

float gPull = shot.Mass * planet.Mass / distanceSqr * planet.gStr;

Answer 1

在最后一行中,您将更新镜头的位置.你应该更新速度.

您可能需要查看此博客中的代码http://blog.mendeltsiebenga.com/post/Fun-with-planets.aspx没有xna,但是工作轨道机制.(虽然我从来没有摆脱屏幕闪烁)

Answer 2

Newton-Raphson迭代不是一种解决这个问题的稳定方法(也就是说,使用如此简单的微分方程积分器,你无法做到正确).考虑使用第二个(或更高)订单解决方案:Runge-Kutta很好,在这种情况下很容易实现.

从数值分析的角度来看,轨道力学问题减少到求解耦合微分方程组的问题:

x_i'' + G m_i \sum_{i != j} m_j r_ji/(|r_ji|)^3 = 0

其中x三个向量代表物体的位置,m它们是相同物体的质量,r_ji = x_j - x_i是物体j和物体之间的向量位移i.

Answer 3

“跳蛙”方法非常有效和稳定，适用于任何动态粒子/场系统，包括等离子体。对于重力，这很简单。以下是您在单个行星上的单次迭代所做的全部工作（1 体问题，围绕静止太阳的单个行星）：

    public void Push()
    {
        Position += Gravity.dT * Velocity;
        Velocity += Gravity.dT * GravityAccelerationVector(Position);
    }

其中“Gravity.dT”是一个统一的时间步长，以任意的时间度量。我正在使用 System.Windows.Vector，但任何自定义 Vector 类都可以，只要它支持基本的乘法和加法。诀窍是 Position 和 Velocity 不是“同时”的，这对于大多数积分方法来说很常见。相反，它们是交错的。迭代 N 上的位置基于迭代 N - 1/2 的速度更新，但随后迭代 N+1/2 的速度基于迭代 N 的位置更新。

N-body 版本是这样的：

    public static void PushPlanets(Planet[] planets)
    {
        // Position Push at iteration N + 0:
        foreach(var p in planets)
            p.Position += Gravity.dT * p.Velocity; // Velocity from N - 1/2

        // Velocity Push at iteration N + 1/2:
        foreach (var p in planets)
        {
            Vector TotalGravity = new Vector(0,0);
            foreach (var pN in planets)
            {
                if (pN == p) continue;
                TotalGravity += pN.Mass * p.Mass * GravityAccelerationVector(p.Position - pN.Position);
            }
            TotalGravity += Sun.Mass * p.Mass * GravityAccelerationVector(p.Position); // Solar acceleration
            p.Velocity += Gravity.dT * TotalGravity;
        }

在哪里

    public static Vector GravityAccelerationVector(Vector position)
    {
        return Vector.Multiply(-G / position.LengthSquared / position.Length, position);
    }

N 体只是更复杂，因为它不是一个单一的引力源，而是有几个。但代码格式是相同的：每个行星的位置都受到 N-1/2 速度的推动，然后我们根据新位置计算出每个行星上的总重力加速度，然后我们通过该总加速度推动每个行星的速度.

其他方法，即使是高阶方法，通常也不稳定，因为它们同时基于位置和速度线性投影下一步。这总是在为系统增加能量方面出错，轨道会逐渐向外移动。其他方法可以过度补偿这种自然误差并从系统中去除能量。一般来说，人们想要一种能源中性的解决方案。蛙跳法会在轨道相位方面逐渐出错，但不会在整体能量方面出错。