Где будет тело через время T
Тело движется вокруг тяготеющего центра — планета вокруг Солнца, станция вокруг Земли — и ничего, кроме этого центра, его не тянет. Ускорение направлено к центру и равно μ/r², где r — расстояние, а μ — гравитационный параметр центра:
$$\mathbf{r}'' = -\frac{\mu}{|\mathbf{r}|^3}\,\mathbf{r}.$$По положению и скорости в момент 0 найдите положение в момент T. Единицы любые, но согласованные: километры и секунды, астрономические единицы и сутки.
Второй пример — Земля 17 октября 1926 года (относительно Солнца, по данным JPL Horizons), через год. Среди закрытых тестов — комета Галлея с того же дня, МКС 17 октября 2026 года и астероид Апофис; все — в задаче двух тел, то есть без притяжения остальных планет.
Вход
В первой строке число запросов q (1 ≤ q ≤ 20). В каждой из следующих строк восемь чисел: μ, x, y, z, vx, vy, vz, T. Орбита — эллипс (полная энергия отрицательна), эксцентриситет не больше 0,996, а T не больше двадцати периодов обращения.
Выход
Для каждого запроса — три числа: x, y, z в момент T. Каждое принимается с ошибкой до 10⁻⁶ от модуля правильного ответа (для ответов меньше единицы — с абсолютной ошибкой 10⁻⁶).
Примеры
Пример 1
2 1 1 0 0 0 1 0 3.1415926535897931 1 1 0 0 0 1 0 6.2831853071795862
-1 -2.953085909e-12 0 1 -1.129566406e-11 0
Пример 2
1 0.00029591310798513925 0.91110800447119999 0.40365312087548 7.96530645313e-05 -0.0072492298972944994 0.015666005457084169 2.515955426257e-06 365.25
0.9110958776 0.4036793268 7.965727319e-05
Ограничения
- Время · Python
- 8000 мс
- Время · Go
- 2000 мс
- Память
- 256 МБ
- Проверка
- 9 тестов
Разбор
Открывайте по ступеням: каждая следующая подсказывает больше. Можно сразу к решению — но задача, решённая с одной подсказкой, запоминается лучше.
1Подсказка 1
Перепишите уравнение второго порядка как систему первого: состояние — шесть чисел, положение и скорость, а производная состояния — скорость и ускорение. Дальше подходит любой метод из статьи. Какой шаг выбрать?
2Подсказка 2
Посмотрите на третий тест: комета Галлея десятилетиями ползёт вдали от Солнца, а у Солнца разворачивается за несколько недель. Один шаг на всю орбиту либо слишком мелок вдали, либо слишком груб у Солнца.
3Идея
Нужен шаг, который подстраивается сам: пара Дормана — Принса оценивает ошибку каждого шага по разнице ответов пятого и четвёртого порядка и укорачивает шаг там, где решение меняется быстро. Допуск должен быть с большим запасом: ошибки шагов копятся, а ответ проверяют с точностью 10⁻⁶ после десятков оборотов.
4Разбор
Шаг Дормана — Принса: семь наклонов по таблице, ответ пятого порядка и оценка ошибки — разность весов двух порядков, умноженная на те же наклоны. Ошибку каждой координаты делим на масштаб: начальное расстояние для положений и начальную скорость для скоростей. Если наибольшее отношение не больше допуска 10⁻¹², шаг принимаем, а седьмой наклон, взятый в конце шага, становится первым наклоном следующего. В любом случае новый шаг — старый, умноженный на 0,9·(1/ошибка)1/5, но не больше чем в 5 раз и не меньше чем в 0,2 раза за раз. Последний шаг укорачиваем, чтобы попасть ровно в T.
Почему не хватает РК4 с постоянным шагом. Две тысячи шагов на запрос — щедро для Земли за год, но комете Галлея за 60 лет это шаги по 11 суток, а в перигелии она за 11 суток поворачивает на десятки градусов, и метод проскакивает поворот: координаты ответа уезжают на десятки процентов. МКС за сутки делает 16 оборотов, и две тысячи шагов — это 125 на оборот: ошибка около 2·10⁻⁴, в сотни раз больше допустимой.
Есть и путь вообще без шагов: в задаче двух тел орбита — эллипс, и положение в любой момент находится из уравнения Кеплера, которое решают методом Ньютона или бисекцией (см. статью о бисекции и методе Ньютона). Такое решение точнее и быстрее, но работает только для двух тел; численное — для любых сил.
4+Решение
5Сложность
Число шагов зависит от орбиты: на почти круговой — несколько сотен на оборот при допуске 10⁻¹², на вытянутой большая их часть уходит на окрестность перицентра. Каждый шаг — шесть вычислений ускорения, O(1) каждое.
6Типичные ошибки
- Постоянный шаг. Проходит первые два теста и падает на комете, станции и случайных вытянутых орбитах.
- Эйлер. Даже миллион шагов не дают шести знаков: ошибка у него убывает только как первая степень шага, а энергия орбиты растёт.
- Слабый допуск. Допуск 10⁻⁶ на шаг не означает 10⁻⁶ в ответе: ошибки сотен шагов складываются и раздуваются. Берите на несколько порядков строже, чем требуется в ответе.
- Масштаб ошибки. Абсолютный допуск в километрах и в астрономических единицах — разные вещи. Нормируйте ошибку на размер орбиты.