Метод Рунге-Кутты 4 порядка и гравитационное взаимодействие.

[polonij]

Доброго времени суток!
Я пытаюсь реализовать программу для моделирования классического(нерелятивистского) гравитационного взаимодействия N-тел с помощью метода Рунге-Кутты, но немного запутался в его интерпретации, поэтому не совсем уверен в своих рассуждениях. Буду очень благодарен если поможете разобраться.
Моделирование производится на основе Ньютоновских уравнений движения:

$$ \begin{cases} \frac {d \vec{v}} {d t} = \vec{a}, \\ \frac {d \vec{r}} {d t} = \vec{v}. \\ \end{cases} $$

,где $\vec{a}$ и $\vec{v}$ - ускорение и скорость соответственно. Ускорение считается самым простым способом(в планковской системе единиц):

$$\vec{a_n} = \sum^{N}_{k=0} {m_k \frac {(\vec{r_n} - \vec{r_k})} {|r_n - r_k|^3}}$$

, здесь $n$ и $k$ - номера тел.
Проблема заключается в том, что для работы метода в правой части уравнений должна быть функция от
двух переменных, но, на мой взгляд, и $\vec{a}$ и $\vec{v}$ зависят только от времени (не укладывается в голове, чтобы ускорение зависело еще и от скорости, а скорость - от координаты). Если это предположение верно, то все, что остается - это считать функции независимыми от второй переменной, и тогда итерационные формулы будут выглядеть так:

$$\vec{r}_{n+1} = \vec{r_n} + \frac {1} {6} \left(\vec{k_1} + 2\vec{k_2} + 2\vec{k_3} + \vec{k_4}\right), \text{где:}$$

$$\vec{k_1} = \vec{v}\Delta t$$

$$\vec{k_2} = \vec{k_3} = \left(\vec{v} + \frac {1} {2} \vec{a}\Delta t\right)\Delta t$$

$$\vec{k_4} = (\vec{v} + \vec{a}\Delta t)\Delta t$$

И аналогичные слагаемые для скорости:

$$\vec{k_1} = \vec{a}(t)\Delta t$$

$$\vec{k_2} = \vec{k_3}= \vec{a}(t + \frac {\Delta t} {2})\Delta t$$

$$\vec{k_4} = \vec{a}(t + \Delta t)\Delta t$$

в которых $\vec{a}(t)$ рассчитывается исходя из положений тел в указанные моменты времени найденных с помощью метода Эйлера.
Верно ли вышенаписанное? Это все еще метод Рунге-Кутты? Меня смущает, что формула для расчета координат при сворачивании превращается в школьное уравнение равноускоренного движения, а также, что $\vec{k_2}$ и $\vec{k_3}$ всегда равны между собой и из-за чего создается впечатление, что порядок ошибки будет ниже, чем $O(\Delta t^4)$.

[Xaositect]

polonij в сообщении #580096 писал(а):

Проблема заключается в том, что для работы метода в правой части уравнений должна быть функция от
двух переменных, но, на мой взгляд, и $\vec{a}$ и $\vec{v}$ зависят только от времени (не укладывается в голове, чтобы ускорение зависело еще и от скорости, а скорость - от координаты)..

Но Вы же сами перед этим написали, как ускорение зависит от координаты:

Цитата:

$$\vec{a_n} = \sum^{N}_{k=0} {m_k \frac {(\vec{r_n} - \vec{r_k})} {|r_n - r_k|^3}}$$