Вы допустили две очевидные ошибки и неудачный выбор алгоритма.
Во-первых, вы смоделировали положение луны как одно число, а не как вектор в 3-пространстве, так что выМоделируем здесь простое гармоническое движение, похожее на пружину, ограниченную одним измерением, а не на орбиту.Начните с моделирования позиций и скоростей в виде 3-пространственных векторов, а не двойных.
Во-вторых, вы сделали знак гравитационной силы положительного , так что это отталкивающий сила между двумя телами, а не привлекательная .Знак силы должен быть в направлении Земли .
В-третьих, ваша реализация алгоритма Эйлера кажется правильной, но Эйлер - плохой выбор для численного решения орбитальных задач.потому что это не консервативно;Вы можете легко попасть в ситуации, когда Луна получает или теряет немного импульса, и это складывает и разрушает вашу красивую эллиптическую орбиту.
Поскольку орбита Луны гамильтонова, используйте вместо этого симплектический алгоритм;он предназначен для моделирования консервативных систем.
https://en.wikipedia.org/wiki/Symplectic_integrator
Этот подход и ваш эйлеровский подход в основном направлены на поиск векторов состояний:
https://en.wikipedia.org/wiki/Orbital_state_vectors
Однако для вашей простой системы из двух тел существуют более простые методы.Если вы хотите сделать симуляцию, подобную космической программе Kerbal, в которой только ваша орбита влияет на вашу орбиту, а планеты с несколькими лунами «на рельсах», вам не нужно моделировать систему на каждую единицу времениразработать последовательность векторов состояний;Вы можете вычислить их непосредственно, учитывая кеплеровские элементы:
https://en.wikipedia.org/wiki/Orbital_elements
Мы знаем шесть элементов для Луны с высокой степенью точности;из них вы можете вычислить положение в 3-пространстве Луны в любое время, опять же, предполагая, что ничто не мешает его орбите.(На самом деле орбита Луны меняется от солнца, от того, что Земля не является сферой, от приливов и т. Д.)
ОБНОВЛЕНИЕ:
Прежде всего, так как это для курсовой работы , вам необходимо процитировать все ваши источники, и это включает в себя получение помощи из Интернета .Ваши инструкторы должны знать, какую работу вы выполняете и какую работу вы выполняете для вас.
Вы спрашивали, как выполнять эту работу в двух измерениях;это кажется неправильным, но что бы там ни было, делайте то, что говорит курсовая работа.
Правило, которое я хотел бы обучить большему числу начинающих: создайте тип, метод или переменную, которая решит вашу проблему ,В этом случае мы хотим представить поведение комплексного значения , поэтому это должен быть тип , а этот тип должен быть типом значения .Типы значений struct
в C #.Итак, давайте сделаем это.
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, чтобы получить новое ускорение, а затем вызывать Эйлера, чтобы получить новое состояние, и повторять.
Как новичок, внимательно изучите эти приемы .Обратите внимание на то, что я сделал: я создал типы для представления концепций и методы для представления операций в этих концепциях.Как только вы это сделаете, программа станет чрезвычайно понятной для чтения.
Подумайте, как бы вы расширили эту программу, используя эти методы;мы сделали жестко AccelerationOfTheMoon
метод, но что, если мы хотим вычислить ускорения нескольких тел?Опять же, делает тип для представления Body
;тело характеризуется State
и массой.Что, если мы хотим решить проблему n-тела?Наше вычисление ускорения должно было бы принять другие тела в качестве параметра.И так далее.