Как добавление силового поля обратного куба к силовому полю обратного квадрата меняет орбиту?

Начнем с того, что запишем кинетическую энергию для планеты (строго говоря, нужно учитывать и движение солнца, но для ртути$M_{\rm sun} \gg M_{\rm planet}$поэтому я буду игнорировать движение центра масс)

$$T(r,\theta,\phi) = \frac{1}{2}m(\dot{r} + r^2\dot{\theta}^2 + r^2\sin^2\theta\dot{\phi}^2) \tag{1}$$

А потенциальная энергия равна

$$V(r) = -G\frac{Mm}{r^2}+ \alpha\frac{l^2}{2m^2 r^2} \tag{2}$$

Обратите внимание, что${\bf F} = -\nabla V$, это хорошая проверка. Затем лагранжиан можно записать как (ниже вы увидите , что это полностью эквивалентно написанию уравнений Ньютона)

$$L = T - V = \frac{1}{2}m(\dot{r}^2 + r^2\dot{\theta}^2 + r^2\sin^2\theta\dot{\phi}^2) + G\frac{Mm}{r^2} - \alpha\frac{l^2}{2m^2 r^2} \tag{3}$$

Теперь нам просто нужно вывести уравнения движения. В общем, если$q$является обобщенной координатой ($q = \{r,\theta,\phi\}$) затем

$$\frac{{\rm d}}{{\rm d}t}\left(\frac{\partial L}{\partial \dot{q}}\right) - \frac{\partial L}{\partial q} = 0 \tag{4}$$

известное как уравнение Эйлера-Лагранжа$^\color{red}{1}$.

$q = \phi$

$$\frac{{\rm d}}{{\rm d}t}\left(\frac{\partial L}{ \partial \dot{\phi}}\right) - \frac{\partial L}{\partial \phi} = \frac{{\rm d}}{{\rm dt}}(mr^2\sin^2\theta \dot{\phi}) - (0) = 0$$

Из этого вы заключаете, что член в скобках является константой, которая есть не что иное, как угловой момент$l$. Я не буду вдаваться в подробности (это становится довольно длинным), но это гарантирует, что движение будет на плоскости:$\theta = \pi/2$, у нас есть тогда

$$\frac{{\rm d}}{{\rm d}t}(mr^2\dot{\phi}) = 0 ~~~\Rightarrow~~~ l = mr^2\dot{\phi} \tag{5}$$

$q = r$

Для этого случая получаем (не забываем учесть, что$l$является константой)

$$\frac{{\rm d}}{{\rm d}t}(m\dot{r})-\left(-G\frac{M m}{r^2} + \alpha\frac{l^2}{m^2r^3}+ mr \dot{\phi}^2\right) = 0$$

Что сводится к

$$m\frac{{\rm d}}{{\rm d}t}\dot{r} = - G\frac{mM}{r^2} + \frac{l^2}{mr^3}\left(1 - \frac{\alpha}{m} \right) \tag{6}$$

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

$$m\frac{{\rm d}}{{\rm d}t} \dot{r} = -\frac{{\rm d}}{{\rm d}r}V_{\rm eff}(r)$$

для некоторого потенциала$V_{\rm eff}$(пожалуйста, выведите его) и нахождение стабильных точек. Переменная$r$будет больше вокруг такой точки

Я оставлю оставшиеся детали для вас. Удачи!

---

$^\color{red}{1}$Вы можете убедиться, что это сводится к уравнениям Ньютона на простом примере, в 1D$T = m\dot{x}^2/2$так что$L_{1D} = m\dot{x}^2/2-V$и уравнение (4) становится

$$\frac{{\rm d}}{{\rm d}t}(m \dot{x}) - (-V'(x)) = 0 ~~~\Rightarrow~~~ m\ddot{x} = -V'(x) = F(x)$$

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

$$\frac{2 \pi}{\sqrt{3+\tfrac{F'(r_0) r_0}{F(r_0)}}}.$$ Когда $F$является обратным квадратом, мы имеем $r F'/F = -2$, так что перигелии точно $2 \pi$друг от друга, иными словами, орбиты смыкаются. Чтобы решить вашу проблему, вы можете посмотреть, как быстро прецессирует орбита Меркурия, и сопоставить это с предложенным членом коррекции обратного куба.

Вывод Предположим, что наша орбита имеет угловой момент$L$. Это означает, что его угловая скорость$L/(m r^2)$, поэтому центростремительная сила$m r \left( \tfrac{L}{m r^2} \right)^2 = \tfrac{L^2}{m r^3}$. Другими словами, мы имеем

$$m r'' = -F(r) + \frac{L^2}{m r^3}.$$ Если бы у нас была идеальная круговая орбита радиусом $r_0$, мы бы хотели иметь $F(r_0) = \tfrac{L^2}{m r_0^3}$.

Теперь предположим, что мы немного искажаем нашу идеальную круговую орбиту, чтобы$r_0+\epsilon(t)$. Если$\epsilon$достаточно мало, чтобы мы могли аппроксимировать формулы их рядом Тейлора до линейного порядка, то

$$m \epsilon'' = -F(r_0) - F'(r_0) \epsilon + \frac{L^2}{m r_0^3} - \frac{3 L^2}{m r_0^4} \epsilon = - \left(F'(r_0)+ \frac{3 L^2}{m r_0^4} \right) \epsilon.$$ Итак, согласно стандартному анализу массы на пружине, $\epsilon$колеблется с периодом $$2 \pi \sqrt{\frac{m}{F'(r_0) + \tfrac{3 L^2}{m r_0^4}}}.$$

Чтобы получить угол, пройденный между колебаниями, умножьте период на угловую скорость$L/(m r_0^2)$получить

$$\theta = 2 \pi \sqrt{\frac{m}{F'(r_0) + \tfrac{3 L^2}{m r_0^4}}} \frac{L}{m r_0^2}=\frac{2 \pi}{\sqrt{3+\tfrac{m F'(r_0) r_0^4}{L^2}}}.$$ Наконец, помните, что $F(r_0) = \tfrac{L^2}{m r_0^3}$переписать это как $$\frac{2 \pi}{\sqrt{3+\tfrac{ F'(r_0) r_0}{F(r_0)}}}.$$

Вопрос 1:

Это ничего качественно не меняет в том смысле, что ваша задача двух тел принимает, формально говоря, такие же решения, как и задача двух тел Кеплера, т. е. коники, но теперь параметризованные$\alpha$.

Задачу двух тел лучше всего решать в системе центра масс системы. В этой системе отсчета это выглядит как задача о центральном силовом поле, примененная к частице с уменьшенной массой.$\mu = \frac{Mm}{m+M}$. Так как в нашем случае$M\gg m$, по предположению движением массы M вокруг центра масс можно смело пренебречь (центр масс совпадает со звездой). В этом приближении приведенная масса$\mu\sim m$и полный угловой момент системы двух тел$\vec{L}=\mu \vec{r}\times \dot{\vec{r}}$с планеты$\vec{L}\sim m \vec{r} \times \dot{\vec{r}}$.

Центральная проблема в системе отсчета с неподвижной звездой имеет силовое поле, определяемое формулой

$$\vec{F} = \frac{-GMm}{r^2} \hat{r} +\alpha\frac{L^2}{m^2r^3} \hat{r} = -\frac{dU}{dr} \hat r$$ происходящие из потенциальной энергии, $$U(r) = \frac{GMm}{r}-\frac{1}{2}{\alpha}\frac{L^2}{m^2r^2}$$

Хорошо известен результат, что в центральном силовом поле угловой момент$L$является первым интегралом движения, т.е. остается постоянным,

$$\dot{\vec{L}}=0.$$ Это означает, что движение в любое время остается плоским, в плоскости, ортогональной угловому моменту. Используя полярные координаты в этой плоскости единичного вектора нормали $\hat{e}_z$, надо, ${L} = L \hat{e}_z = mr^2\dot \phi \hat{e}_z$, $$L=mr^2\dot \phi = \mbox{Cte}.$$

Другой первый интеграл движения есть полная энергия$E$системы, определяемой как сумма кинетической энергии и потенциальной энергии планеты,

$$E=\frac{1}{2} m(\dot r^2 + r^2\dot \phi^2)+U(r)=\frac{1}{2} m {\dot r^2} + \frac{L^2}{2mr^2}+U(r)$$ .

Это уравнение энергии показывает, что движение в центральном поле формально похоже на «линейное или поступательное» движение (т. е. с одной степенью свободы, радиусом$r$), а в поле "действующей" потенциальной" энергии,

$$U_{\rm eff}(r)=\frac{L^2}{2mr^2}+U(r) = (1-\alpha)\frac{L^2}{2mr^2}+\frac{GMm}{r}$$

Теперь настройка$L_\alpha = L\sqrt{1-\alpha}$, получается,

$$U_{\rm eff}(r)=\frac{L^2}{2mr^2}+U(r) = \frac{L_\alpha^2}{2mr^2}+\frac{GMm}{r}$$ .

и полная энергия

$$E=\frac{1}{2} m {\dot r^2} + U_{\rm eff}(r) = \frac{1}{2} m {\dot r^2} + \frac{L_\alpha^2}{2mr^2}+ \frac{GmM}{r}$$

Оба уравнения (эффективный потенциал или полная энергия) описывают движение классической кеплеровской системы с двумя телами, которая имела бы$L_\alpha$вместо$L_0=L$как полный угловой момент. Таким образом, мы приходим к выводу, что замена$L$к$L_\alpha$во всех результатах классической кеплеровской теории двух тел дает правильное решение вашей проблемы. Получается такое же движение, как и для классической кеплеровской системы двух тел (случай$\alpha =0$), но траектории вместо того, чтобы полностью определяться парой$(E,L_0)$теперь управляются$(E, L_\alpha)$. Решения для$\alpha\neq 0$выводятся из кеплеровского случая$\alpha=0$гомотетически по соотношению$\sqrt{1-\alpha}$(возможно, чисто мнимый комплекс, когда$\alpha>1$).

С гравитационным силовым полем (с потенциалом$\frac{GMm}{r}$) --- и вообще говоря, в каждом кулоновском силовом поле (потенциального$-k/r$) можно легко проинтегрировать траекторию как конику уравнения,

$$p/r=1+e\cos\phi$$

где,

$$p=\frac{L^2}{km}=-\frac{L^2}{GMm^2}$$ называется параметром и $e$эксцентриситет коник определяется выражением $$e=\sqrt{1+\frac{2EL^2}{mk^2}}=\sqrt{1+\frac{2EL^2}{G^2M^2m^3}}$$ .

Заменять$L \rightarrow L_\alpha$в этом классическом результате, и вы демонстрируете, что в вашем случае вы получаете (семейство) коник (индексируемых$\alpha$) соответствующего параметра(ов)$p_\alpha$и эксцентриситет(ы)$e_\alpha$.

Когда$\alpha\leq 1$центробежный термин$\frac{L_\alpha^2}{2mr^2}$в эффективном потенциале положителен. (и) если$E<0$, движение ограничено,$e_\alpha$< 1 и вы получите эллипс; (ii) если$E=0$у вас получится круглый круг. (iii) если$E>0$, движение неограниченно,$e_\alpha$> 1 и вы получите гиперболу. Эти результаты остаются в силе для$\alpha=1$, просто это кеплеровский случай нулевого углового момента, дающего вырожденный эллипс, то есть движение, ускоренное по прямой линии к центральной звезде.

Когда$\alpha>1$центробежный термин$\frac{L_\alpha^2}{2mr^2}$в эффективном потенциале становится отрицательным. Я не вижу существенных изменений в предыдущих результатах, за исключением того, что теперь можно получить ограниченное движение по эллипсам, когда$E>0$и неограниченное движение по гиперболам при$E<0$.

Вопрос 2:

Бьюсь об заклад, что$\alpha=0$поскольку включение вашего дополнительного силового термина на самом деле не обязательно. Предположим, вы вводите термин силы, скажем, с некоторым$L$и$\alpha=0.1$. Ваша модель в некотором смысле эквивалентна модели без дополнительной силы, но с импульсом планеты.$L$уменьшается до адекватного значения, т.е. умножается на$\sqrt{1-\alpha}$. Как физик, я предпочел бы это более простое формальное решение.

Это задача с двумя телами, притягивающими друг друга консервативной центральной силой. Вы можете решить эту проблему, рассматривая движение только одного тела уменьшенной массы.$\mu = \frac{m_1 m_2}{m_1 + m_2}$движение в плоскости - в этом случае нетрудно показать, что момент количества движения системы сохраняется, а значит, движение происходит в плоскости - с вектором положения$\vec{r} = x\vec{i} + y\vec{j}$. Вы можете думать об этом движении так, как если бы референтом была одна из частиц, а другая имела бы массу.$\mu$, а вектор положения - их относительное положение. Чтобы найти параметризацию траектории, в полярных координатах плоскости, в которой происходит движение, решите для$u = 1/r$следующее уравнение:

$$\frac{d^2 u}{d \theta} + u = - \frac{\mu}{l^2}\frac{1}{u^2}F(1/u)$$ Здесь, $l$обозначает величину углового момента (который я предполагаю постоянным в этом случае).

В качестве источников по этому вопросу у вас есть превосходная «Классическая механика» Гольдштейна и «Классическая динамика частиц и систем» Мариона и Торнтона (главы 3 и 8 соответственно). Теперь, обозначая вашу силу через$F(r) = -\frac{k_1}{r^2} - \frac{k_2}{r^3}$(векторная сила$\vec{F}(r) = F(r) \frac{\vec{r}}{r}$), дифференциальное уравнение принимает вид:

$$\frac{d^2 u}{d \theta} + \alpha u = \beta$$ Где $\alpha$и $\beta$константы в функции $\mu$, $l$, $k_1$и $k_2$.

Если$\alpha > 0$, решение имеет вид:$\frac{1}{r} = Acosh(\alpha\theta + B) + \frac{\beta}{\alpha}$, где$A$и$B$определяются начальными условиями.

Если$\alpha = 0$, решение имеет вид$\frac{1}{r} = \beta \theta^2 + A\theta + B$

Если$\alpha < 0$,$\frac{1}{r} = Acos(\alpha\theta + B) + \frac{\beta}{\alpha}$.