24.08.2019 14:47
nickname21 8 4900 Источник

Содержание


Часть 1 — Задача двух тел

Часть 2 — Полу-решение задачи двух тел

Движение в плоскости


Осталось сделать последний штрих. Решить это уравнение:

$ \ddot{\vec{r}} = -\mu \dfrac{\vec{r}}{r^{3}}, $



где $ \mu = G(m_{1} + m_{2}) > 0, $ $ \vec{r} $ — относительное расстояние между телами.

В прошлом выпуске было показано, что при значительном различии масс (например $ m_{1} >> m_{2} $) вектор $ \vec{r} $ можно считать радиусом вектором в новой системе координат связанной с неподвижным массивным телом. Неподвижное оно потому, что центр масс совпадает с ним, а центр масс движется равномерно прямолинейно (тоже в прошлом выпуске доказали).

Но всё это лишь для удобства восприятия нами. Предельный случай, когда массивное тело неподвижно и сосредоточено в начале системы координат, а второе легкое тело вращается вокруг, прямо как спутник на орбите Земли. А $ \vec{r} $ как раз описывает траекторию движения спутника.

Земля конечно не движется равномерно и прямолинейно и тем более не является неподвижной находящейся в центре Вселенной (хотя это наша воля назначить центр). Да и точкой не является, по крайней мере в примере со спутником.

image
Земля и спутник

Но всё таки, если непродолжительное время понаблюдать за Землей со стороны, то нам будет казаться, что она (Земля) пролетела в космосе по прямой, причем равномерно. Вот вы чувствуете центростремительное ускорение вызываемое вращением Земли вокруг Солнца? И я нет. Мы даже не чувствуем (вестибулярным аппаратом), что она (Земля) вращается вокруг своей оси. И сила Кориолиса практически не проявляется. Так может Земля — плоская? Лично я не проверял :D Но давайте будем верить официальной науке, что всё это благодаря масштабам. Вот если взять сделать центрифугу в космосе для создания искусственной гравитации, то по простым расчетам диаметр (или радиус?) должен быть не меньше километра для более менее комфортного жилья. Ибо сила Кориолиса будет создавать непривычные эффекты: бросишь мячик прямо — он улетит куда-то в сторону. Да просто ходить будет неудобно. Станешь вечно пьяным, вечно молодым испытывать непривычные ощущения, как во сне, как в сказке.

А на счёт того, что Земля представлена точкой — тоже довольно близко к правде. Так как объекты в форме шара создают вокруг себя гравитационное поле, как если бы мы сжали всю массу до точки и сконцентрировали её в центре шара. Это мы потом математически докажем, но сейчас и интуитивно, думаю, понятно. Ведь шар абсолютно симметричен. В реальности конечно же если сжать Землю до 1 см в диаметре — она станет черным отверстием. Хотя всё таки черной дырой, извиняюсь — инженерная привычка (у инженеров не бывает дыр, — только отверстия). Все инженера это знают, ибо работает хорошо эволюционный отбор: если начинающий инженер скажет на отверстие «ДЫРКА» в присутствии более опытного инженера, — сразу получит подзатыльник в лучшем случае. В худшем нужно ждать вертуху в щи (шутка).

image
Отверстие и дырка

Перед тем как решать полезно бывает сначала численно по-моделировать, возможно, это позволит заметить какие-то закономерности. Уравнения уже попроще, чем в прошлый раз:

\begin{equation*}
\begin{cases}
\ddot{x} = -\mu \dfrac{x}{(x^{2}+y^{2}+z^{2})^{\frac{3}{2}}},
\\
\ddot{y} = -\mu \dfrac{y}{(x^{2}+y^{2}+z^{2})^{\frac{3}{2}}},
\\
\ddot{z} = -\mu \dfrac{z}{(x^{2}+y^{2}+z^{2})^{\frac{3}{2}}}.
\end{cases}
\end{equation*}
Берем и вбиваем это куда-нибудь, что умеет решать. В мат-пакет, либо прямо в Python :) Я в Python не стану вбивать, хотя могу. Это потому что мне 3d визуализация matplotlib не нравится. Плюс еще с сохранением gif-ок всё не очень просто, было начал искать — нашел лишь сохранение в видео-формате, так еще и устанавливать дополнительно всякого го*на нужно.

Вот и повод поинтересоваться у аудитории Хабра: можно ли делать подобные 3d анимации как в этой статье, но в matplotlib? Чтобы перспектива была такая красивая (там где я делаю — можно перспективу регулировать). Или посоветуйте какие классные библиотеки Python для визуализации, буду благодарен.

Короче вот несколько вариантов с разными начальными условиями:

image
анимация 1
image
анимация 2
image
анимация 3

Если еще 10-100 раз смоделировать, то наверняка мы заметим, что маленький шарик движется всегда в одной плоскости. И эту плоскость можно вычислить, имея начальные условия. Ведь радиус вектор $ \vec{r} $ всегда лежит в этой плоскости, но и скорость $ \vec{v} $ лежит в ней же. А если всегда, то и в начальный момент времени. Значит мы в состоянии вычислить эту плоскость.

В данном случае удобно задавать плоскость с помощью вектора нормали и любой точки лежащей в ней. Точки у нас аж две — начало системы координат и начальное положение спутника, выбирай любую — выберем начало $ (0, 0, 0) $. А вектор нормали легко получить векторно умножив вот это:

$ \vec{h} = [\vec{r}_{0}, \vec{v}_{0}], $


где $ \vec{r}_{0} = \vec{r}(0), \vec{v}_{0} = \dot{\vec{r}}(0) $. $ \vec{h} $ назовём кинетическим моментом системы. Потому что похож на момент импульса:

$ \vec{L} = [\vec{r}, \vec{p}] = [\vec{r}, m\vec{v}] $


Можно даже назвать кинетический момент — удельным моментом импульса:

$ \vec{h} = \dfrac{\vec{L}}{m} $


image
Вектор нормали

Можно ли проверить с помощью уравнений, что $ \vec{h} $ действительно всегда будет параллельным себе же в начальный момент времени?

Ну если это так, то для любого момента времени должно быть выполнено:

$ k(t)\vec{h} = [\vec{r}, \dot{\vec{r}}], $


здесь $ k(t) $ — какая-то функция времени, которая регулирует длину $ \vec{h} $, ведь $ \vec{h} $ в принципе может быть параллелен себе, но по длине изменяться.

Ну давайте продифференцируем это равенство, используя все правила дифференцирования векторных произведений и их свойств:

$ \dot{k}(t)\vec{h} = [\dot{\vec{r}}, \dot{\vec{r}}] + [\vec{r}, \ddot{\vec{r}}] = [\vec{r}, \ddot{\vec{r}}]. $


Как хорошо получается, у нас же:

$ \ddot{\vec{r}} = -\mu \dfrac{\vec{r}}{r^{3}}, $


ну так и подставим туда:

$ \dot{k}(t)\vec{h} = [\vec{r}, -\mu \dfrac{\vec{r}}{r^{3}}] = -\dfrac{\mu}{r^{3}}[\vec{r}, \vec{r}] = 0. $


А это значит, что $ \dot{k}(t) = 0 $. Или же $ k(t) = const $. Короче говоря $ \vec{h} $ постоянен как по направлению так и по длине. Вот и всё, закон сохранения.

Данный факт можно, конечно, было получить и без всех этих рассуждений и физических смыслов. Просто играясь с самими формулами. Брать равенство и векторно умножить его на $ \vec{r} $ (просто потому что справа тогда будет ноль):

$ \ddot{\vec{r}} = -\mu \dfrac{\vec{r}}{r^{3}}, $


$ [\vec{r}, \ddot{\vec{r}}] = -\dfrac{\mu}{r^{3}}[\vec{r}, \vec{r}] = 0, $


учитывая:

$ \dfrac{d}{dt}[\vec{r}, \dot{\vec{r}}] = [\dot{\vec{r}}, \dot{\vec{r}}] + [\vec{r}, \ddot{\vec{r}}] = [\vec{r}, \ddot{\vec{r}}], $


$ \dfrac{d}{dt}[\vec{r}, \dot{\vec{r}}] = 0 $


теперь разок проинтегрировав, получаем уже полученную вектор константу:

$ [\vec{r}, \dot{\vec{r}}] = \vec{h} = \vec{const}. $


Если $ \vec{h} $ — константа, то она может и нулю равняться. А такое может быть только если $ \vec{r} $ и $ \dot{\vec{r}} $ будут параллельны. Это значит, что скорость будет направлена вдоль радиуса вектора. Именно это и происходило, когда на голову Ньютона падало яблоко.

image
Эх, яблочко.

В принципе, в принципе, так можно продолжать играть дальше. Например, попробовать векторно умножить имеющееся у нас уравнение на $ \dot{\vec{r}} $, или использовать скалярное произведение. В общем комбинировать и надеяться на лучшее. Операций у нас не так уж и много, и свойств тоже, — по пальцам пересчитать, так что не грех. На самом деле мы это и будем проделывать позже. И это даст свои плоды, даже откроет еще один способ решения задачи (получения уравнения траектории), но мы решим обеими способами, так статьи длиннее получатся =)).

А теперь продолжим. Из за того что траектория лежит в плоскости (пускай пока что как-то там в пространстве расположенной), получается задача то у нас — плоская. А это значит, что вместо трех дифференциальных уравнений второго порядка, на самом деле их два. Но в правильно выбранной системе отсчета. Например в такой:

image
Правильная система отсчета

Новый базис связан со старым так (пишу в порядке вычисления):

$ \vec{i^{'}} = \dfrac{\vec{r}_{0}}{|\vec{r}_{0}|} = a_{11}\vec{i} + a_{12}\vec{j} + a_{13}\vec{k} $


$ \vec{k^{'}} = \dfrac{\vec{h}}{|\vec{h}|} = \dfrac{[\vec{r}_{0}, \vec{v}_{0}]}{|[\vec{r}_{0}, \vec{v}_{0}]|} = a_{31}\vec{i} + a_{32}\vec{j} + a_{33}\vec{k} $


$ \vec{j^{'}} = [\vec{k^{'}}, \vec{i^{'}}] = a_{21}\vec{i} + a_{22}\vec{j} + a_{23}\vec{k} $


Или в матричном виде:

$\begin{bmatrix} \vec{i^{'}} \\ \vec{j^{'}} \\ \vec{k^{'}} \end{bmatrix}=\begin{bmatrix} a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & a_{23} \\ a_{31} & a_{32} & a_{33} \end{bmatrix} \begin{bmatrix} \vec{i} \\ \vec{j} \\ \vec{k} \end{bmatrix} $


А можно и так и сяк :)

$ \textbf{e}^{'} = A\textbf{e} $


Соответственно переход наоборот:

$ \textbf{e} = A^{-1}\textbf{e}^{'} $


Любой вектор, например скорость $ \vec{v} $, остается одним и тем же вектором в независимости выбранной нами системы отсчета, то бишь можно расписывать так:

$\vec{v} = v_{x}\vec{i} + v_{y}\vec{j} + v_{z}\vec{k} = v^{'}_{x}\vec{i^{'}} + v^{'}_{y}\vec{j^{'}} + v^{'}_{z}\vec{k^{'}}$


Или. Же.

$ \vec{v} = \begin{bmatrix} v_{x} & v_{y} & v_{z} \end{bmatrix} \begin{bmatrix} \vec{i} \\ \vec{j} \\ \vec{k} \end{bmatrix} =\begin{bmatrix} v^{'}_{x} & v^{'}_{y} & v^{'}_{z} \end{bmatrix} \begin{bmatrix} \vec{i^{'}} \\ \vec{j^{'}} \\ \vec{k^{'}} \end{bmatrix} $


Хорошо статья тянется, как резиновая:

$ \vec{v} = \textbf{v}^{T}\textbf{e} = \textbf{v}^{'T}\textbf{e}^{'} $


Как преобразуются базисы мы в курсе, тогда очень просто найти как преобразуются компоненты векторов при переходе (это я напоминаю аналитическую геометрию, а может кому-то что-то проясниться и мои записи станут полезными):

$ \textbf{v}^{T}\textbf{e} = \textbf{v}^{T}A^{-1}\textbf{e}^{'} = \textbf{v}^{'T}\textbf{e}^{'} $


Ну и очевидно компоненты равны:

$ \textbf{v}^{T}A^{-1} = \textbf{v}^{'T} $


Транспонируем их, а то как неродные:

$ \textbf{v}^{'} = (\textbf{v}^{T}A^{-1})^{T} = (A^{-1})^{T} \textbf{v} $


Это не я. Это всё правила транспонирования.

Карочь, взяли базисы подставили, туда-сюда, сравнили и всё.

Но нас ждет еще один маленький сюрприз. В хорошем смысле слова. У нас матрица поворота. А значит она — ортогональная матрица. Ортогональные матрицы имеют прикольное свойство:

$ A^{-1} = A^{T} $


Тогда переходы получаются — просто прелесть (не бесовская):

$ \textbf{v}^{'} = (A^{-1})^{T} \textbf{v} = (A^{T})^{T} \textbf{v} = A \textbf{v}, $


$ \textbf{v} = A^{T} \textbf{v}^{'}. $


Просто бери матрицу A и переходи в новую систему координат:

$ v^{'}_{x} = a_{11}v_{x} + a_{12}v_{y} + a_{13}v_{z} $


$ v^{'}_{y} = a_{21}v_{x} + a_{22}v_{y} + a_{23}v_{z} $


$ v^{'}_{z} = a_{31}v_{x} + a_{32}v_{y} + a_{33}v_{z} $


Как вызывать вычислять компоненты матрицы A, надеюсь никто не провтыкал (извините за французский). Например:

$ a_{11} = \dfrac{x_{0}}{\sqrt{x_{0}^2 + y_{0}^2 + z_{0}^2}} $


Ну и хватит, пойдем дальше.

Теперь то, из-за того что точка всегда лежит в плоскости $ x^{'}Oy^{'} $, вектор $ \vec{r} $ будет иметь вид:

$ \vec{r} = x^{'}\vec{i^{'}} + y^{'}\vec{j^{'}} + 0\vec{k^{'}} $


Координата $ z^{'} = 0 $ всегда ноль.
Проще стало:

$ \textbf{r}^{'} = \begin{bmatrix} x^{'} \\ y^{'} \\ 0 \end{bmatrix} $


И ускорение:

$ \ddot{\textbf{r}}^{'} = \begin{bmatrix} \ddot{x}^{'} \\ \ddot{y}^{'} \\ 0 \end{bmatrix} $


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

$ \textbf{r} = (x, y) $


И уравнений будет два (уже есть):
\begin{equation*}
\begin{cases}
\ddot{x} = -\mu \dfrac{x}{(x^{2}+y^{2})^{\frac{3}{2}}},
\\
\ddot{y} = -\mu \dfrac{y}{(x^{2}+y^{2})^{\frac{3}{2}}}.
\end{cases}
\end{equation*}
Это ведь мы решаем-решаем, даст Бог решим. А потом главное вернуться поочередно в исходную систему координат. Так что забывать откуда мы пришли — нельзя. После перехода матрицы выбрасывать — нельзя. Лучше заранее найди обратную, как говорится: готовь сани летом, а тележку зимой.

Думаю на этот раз хватит, продолжение будет в следующих статьях. Которые будут длиннее.

Комментарии (8)


  1. nickname21 Автор
    24.08.2019 19:10
    #20546613

    Уже один человек обиделся на опрос и предложил вставить вариант X = 0. Дорогие друзья! Не обижайтесь, я буду продолжать писать даже если 0 человек мне задонатят.


    1. nickname21 Автор
      24.08.2019 19:16
      #20546623

      И где там феня… Пару слов по приколу ))) Хочу пишу так, если не нравиться многим, — я перестану =)


    1. Exosphere
      24.08.2019 21:30
      #20546891
      +2

      А я как модератор, пользуясь случаем, напомню, что вы можете воспользоваться встроенным механизмом Хабра для донейтов. Для этого вам нужно открыть профиль, нажать «Настроить профиль» и прописать реквизиты платёжных систем. Ну а затем напоминать, чтобы пользователи не забывали пользоваться кнопкой в профиле под статьёй.


      1. nickname21 Автор
        27.08.2019 19:09
        #20556039

        К сожалению PayPal.me в Украине не работает. Яндекс и Вэбмани тоже… Печально, не везет


  1. axe_chita
    24.08.2019 21:13
    #20546867
    +1

    Пиши, читать интересно. Вдруг это вырастет и разовьётся до масштабов «Спросите у Итана nickname21»
    П.С.Даже мой вооруженный нейтралитет с математикой хмыкнул и сказал «Ладно читай, эти формулы не картинки из теста Роршаха, а что то интересное и осмысленное»


  1. metaprog
    27.08.2019 00:49
    #20552713

    Есть чудесная игра про освоение космоса — Kerbal space program.


  1. axe_chita
    27.08.2019 05:54
    #20552933

    Ну тогда вспомним и "Orbiter — симулятор орбитальных космических полётов для Windows (выпускается как Freeware с закрытым исходным кодом). Позволяет осуществлять перелеты между планетами и орбитальными станциями в Солнечной системе. Отличается тщательно проработанной реалистичной физикой и богатыми возможностями конфигурирования и написания различных расширений. Использованы реалистичные физические модели движения планет, динамики кораблей и атмосферных явлений. С помощью файлов конфигурации можно модифицировать существующие объекты, добавлять новые космические корабли, орбитальные станции и планетарные базы, а также создать новую солнечную систему.

    Имеется развитое сообщество, разработано множество расширений. Для Орбитера существуют десятки моделей реальных носителей, кораблей и станций. Среди них: Восток/Восход, Mercury, Gemini, Союз, Saturn/Apollo, Skylab, Спейс Шаттл, Энергия/Буран, Салют, Мир, МКС и множество фантастических кораблей и сценариев.

    Может использоваться в образовании, науке и промышленности."


  1. Yurec666
    27.08.2019 16:03
    #20555209

    Столько написано про уравнения движения и ни разу слово тензор не встретилось:(