Орбитальная механика

Есть ли у кого-нибудь пример внедрения Орбитальной Механики (желательно в 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);  
}

Edit Маркировка Mendelt отвечает правильно, указав, что мне нужно обновить скорость, а не позицию. Мне также нужно было изменить вычисление gPull на

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

Ответ 1

В последней строке вы обновляете позицию кадра. Вы должны обновлять скорость.

Возможно, вы захотите взглянуть на код в этом блоге http://blog.mendeltsiebenga.com/post/Fun-with-planets.aspx Нет xna, но рабочая орбитальная механика. (хотя я никогда не избавился от мерцания экрана)

Ответ 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.

Ответ 3

Метод "прыгучей лягушки" очень эффективен и стабилен и хорошо работает для любой динамической системы частиц/полей, включая плазму. Для гравитации это просто. Вот вам все, что вы делаете для одной итерации на одной планете (проблема с одним телом, одна планета вокруг неподвижного солнца):

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

Где "Gravity.dT" является равномерным временным шагом в произвольной мере времени. Я использую System.Windows.Vector, но любой пользовательский класс Vector будет делать, если он поддерживает базовое умножение и добавление. Хитрость заключается в том, что позиция и скорость не находятся в то же время, что очень часто встречается для большинства методов интеграции. Скорее, они пошатнулись. Позиция на итерации 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, тогда мы вычисляем полное ускорение силы тяжести на каждой планете на основе новых позиций, затем мы нажимаем каждую скорость планеты на это полное ускорение.

Другие методы, даже высокого порядка, часто нестабильны, потому что они линейно проецируют следующий шаг, основанный на скорости позиции и в то же время. Это всегда будет ошибочно на стороне добавления энергии в систему, и орбиты будут постепенно двигаться дальше и дальше наружу. Другие методы могут перекомпенсировать эту естественную ошибку и удалить энергию из системы. В общем, требуется энергетически нейтральное решение. Метод скачкообразной лягушки будет постепенно ошибаться по фазе орбит, но не в общей энергии.

Ответ 4

Проходной объект не войдет на орбиту. Одной из характеристик орбиты является то, что вы вернетесь в ту же точку (относительно тела, находящегося на орбите) с той же скоростью. Если вы начали с эффективной бесконечности, вы вернетесь к эффективной бесконечности.

Чтобы войти на орбиту, вам нужно будет изменить скорость в какой-то момент, не связанной с гравитацией, или, возможно, иметь дополнительные большие тела. Точно так же вы не можете запустить объект на орбиту с поверхности: вы должны иметь что-то (например, последний ожог ракеты), когда спутник достигнет желаемой высоты. В противном случае он попытается вернуться к точке запуска, которая находится на поверхности.

Некоторые из моих худших событий отладки были, когда программа была в порядке, и мои тестовые данные или расчеты были отключены. Убедитесь, что вы знаете, что искать.

Ответ 5

A) Мы не знаем, каковы ваши входные значения.

B) Вы можете использовать лучшее приближение, чем Ньютон-Рафсон.

C) Проходящие объекты, как правило, не попадают на орбиту IRL, сила тяжести крайне слаба, для одинаковых слабых скоростей или поистине исключительных масс получается гораздо больше, чем кривизна.