行星轨道的数值逼近

问题描述 投票：1回答：1

新的一天新问题......

我按照Eric所说的方式做到了，我创建了一个加速方法！

它看起来像这样：

static Vector2 AccelerationOfTheMoon(Vector2 position)
    {
        double Mm = 7.349 * Math.Pow(10, 22);

        double MagnitudeOfForce()
        {
            double MagF = position.LengthSquared();
            return MagF;
        }

        Vector2 ForceAsVector()
        {
            return position.NegativeUnitVector() * MagnitudeOfForce();
        }

        Vector2 acceleration = ForceAsVector() / Mm;


    }

我看到它的方式，我的方法称为AccelerationOfTheMoon接收带有位置的Vector 2。现在我想使用这个向量，所以我创建了另一个方法，它应该将我的Magnitude作为标量（double）返回。为此我创建了一个方法来计算矢量的大小并对其进行平方。当我在ForceAsVector-Method中调用我的MagnitudeOfForce-Method时，我试图将负单位向量与标量相乘并将其作为Vector2返回。

也许我在这里做得非常糟糕，但我正在努力去理解埃里克所做的一切，但仍然需要你的一些帮助！

谢谢。

c# vector numerical-methods newtons-method orbital-mechanics
1个回答
11
投票

你犯了两个明显的错误和一个糟糕的算法选择。

首先，你已经将月球的位置建模为单个数字而不是3空间中的矢量，所以你在这里建模的是简谐运动，就像弹簧约束到单个维度而不是轨道。首先将位置和速度建模为3空间向量，而不是双倍。

其次，你已经使引力的符号为正，因此它是两个物体之间的排斥力，而不是吸引力。力的标志必须在地球的方向。

第三，你对欧拉算法的实现似乎是正确的，但欧拉在数值上解决轨道问题是一个不好的选择，因为它不是保守的;你可以很容易地进入月亮增加或失去一点动量的情况，然后加起来，并破坏你漂亮的椭圆轨道。

由于月球的轨道是哈密顿量，所以使用辛算法;它旨在模拟保守系统。

https://en.wikipedia.org/wiki/Symplectic_integrator

这种方法和您的欧拉方法基本上是关于找到状态向量：

https://en.wikipedia.org/wiki/Orbital_state_vectors

但是，对于简单的2体系统，有更简单的方法。如果你想要做的是像Kerbal太空计划那样进行模拟，只有你在轨道上运行的物体影响你的轨道，而且有多个卫星的行星是“在轨道上”，你不需要在每个时间单位模拟系统计算出状态向量序列;你可以直接用开普勒元素计算它们：

https://en.wikipedia.org/wiki/Orbital_elements

我们高度精确地了解月球的六个要素;从那些你可以随时计算月球三维空间中的位置，再假设没有任何扰乱它的轨道。 （实际上，月球的轨道是由太阳改变的，因为地球不是球体，是潮汐，等等。）

更新：

首先，因为这是课程作业，您需要引用所有来源，其中包括从互联网获取帮助。你的导师必须知道你的工作是什么，以及你有别人为你做的工作。

你问过如何在两个方面做这项工作;这似乎是错的，但无论如何，做的是课程的工作。

我希望更多初学者接受教育的规则是：制作解决问题的类型，方法或变量。在这种情况下，我们希望表示复杂值的行为，因此它应该是一个类型，并且该类型应该是值类型。值类型是C＃中的struct。所以，让我们这样做。

struct Vector2
{
    double X { get; }
    double Y { get; }
    public Vector2(double x, double y)
    {
        this.X = x;
        this.Y = y;
    }

请注意，向量是不可变的，就像数字一样。你永远不会改变一个向量。当你需要一个新的矢量时，你会做一个新的矢量。

我们需要对矢量执行哪些操作？向量加法，标量乘法和标量除法只是花式乘法。让我们实现这些：

    public static Vector2 operator +(Vector2 a, Vector2 b) => 
      new Vector2(a.X + b.X, a.Y + b.Y);
    public static Vector2 operator -(Vector2 a, Vector2 b) => 
      new Vector2(a.X - b.X, a.Y - b.Y);
    public static Vector2 operator *(Vector2 a, double b) => 
      new Vector2(a.X * b, a.Y * b);
    public static Vector2 operator /(Vector2 a, double b) => 
      a * (1.0 / b);

我们也可以在其他顺序中进行乘法，所以让我们实现：

    public static Vector2 operator *(double b, Vector2 a) => 
      a * b;

使矢量为负数与将其乘以-1相同：

    public static Vector2 operator -(Vector2 a) => a * -1.0;

并帮助我们调试：

    public override string ToString() => $"({this.X},{this.Y})";

}

我们完成了向量。

我们试图模拟轨道状态参数的演变，所以再次制作一个类型。状态参数是什么？位置和速度：

struct State
{
    Vector2 Position { get; }
    Vector2 Velocity { get; }
    public State(Vector2 position, Vector2 velocity)
    {
        this.Position = position;
        this.Velocity = velocity;
    }

现在我们来看核心算法。我们有什么？状态和加速度。我们想要什么？一个新的国家。制作方法：

    public State Euler(Vector2 acceleration, double step)
    {
        Vector2 newVelocity = this.Velocity + acceleration * step;
        Vector2 newPosition = this.Position + this.Velocity * step;
        return new State(newPosition, newVelocity);
    }
}

超。还剩什么？我们需要计算出加速度。制作方法：

static Vector2 AccelerationOfTheMoon(Vector2 position) 
{
  // You do this. Acceleration is force divided by mass,
  // force is a vector, mass is a scalar. What is the force
  // given the position? DO NOT MAKE A SIGN MISTAKE AGAIN.
}

现在你拥有了所需的所有部件。从初始状态开始，您可以重复调用AccelerationOfTheMoon来获取新的加速，然后调用Euler获取新状态，然后重复。

作为初学者，请仔细研究这些技巧。注意我做了什么：我创建了表示概念的类型，以及表示对这些概念的操作的方法。一旦你这样做，该程序变得非常清晰阅读。

考虑一下如何使用这些技术扩展该程序;我们制作了一个硬编码的AccelerationOfTheMoon方法，但如果我们想计算几个物体的加速度呢？再次，制作一个代表Body的类型;一个身体的特点是State和质量。如果我们想解决n体问题怎么办？我们的加速度计算必须将其他物体作为参数。等等。

最新问题
© www.soinside.com 2019 - 2024. All rights reserved.