% Document MUST be built with XeLaTeX! \documentclass[11pt,a4paper]{article} \usepackage[margin=1.5cm]{geometry} % Поддръжка на български език % Thanks to https://tex.stackexchange.com/a/136926 \usepackage{fontspec} \defaultfontfeatures{Ligatures=TeX} \setmainfont{CMU Serif} \setsansfont{CMU Sans Serif} \usepackage[bulgarian,english]{babel} \usepackage{amsmath} % Математически символи \usepackage{extarrows} \usepackage{amsfonts} % Още математически символи \usepackage{cancel} % Зачертаване \usepackage[dvipsnames]{xcolor} % Цветни кутийки + допълнителни цветове \usepackage{rotating} \usepackage{enumitem} \setlist[enumerate]{label=\arabic*.} \usepackage{listings} \newfontfamily{\lstsansserif}[Scale=1]{DejaVuSansMono} \lstdefinestyle{mystyle}{ basicstyle=\lstsansserif\footnotesize, breakatwhitespace=false, breaklines=true, captionpos=b, keepspaces=true, numbers=left, numbersep=15pt, showspaces=false, showstringspaces=false, showtabs=false, tabsize=2, frame=single, literate = {-}{-}1{"}{\textquotedbl}1{*}{\normalfont{*}}1, } \lstset{style=mystyle} \newcommand{\rad}{\mathrm{rad}} \usepackage[pdfauthor={Камен Димитров Младенов, }, pdftitle={Решения на задачи за изпит по НМ}, pdfsubject={}, pdfkeywords={}, colorlinks=true, urlcolor=blue]{hyperref} \title{\vspace{-1.5cm}Решения на задачи за изпит по НМ} \author{Камен Димитров Младенов, FN} \date{DATE} \begin{document} \maketitle \vfill Първо ще подготвим таблицата на орбиталните елементи. Уебсайтът на NASA \cite{nasa} предоставя следните константи за началото на епохата (1 Януари 2000, 00:00 часа): \begin{table}[!h] \begin{center} \begin{tabular}{c||c|c|c|c|c|c} & a & e & I & L & long.peri. & long.node. \\ & au, au/Cy & rad, rad/Cy & deg, deg/Cy & deg, deg/Cy & deg, deg/Cy & deg, deg/Cy \\ \hline\hline Mercury & 0.38709927 & 0.20563593 & 7.00497902 & 252.25032350 & 77.45779628 & 48.33076593 \\ & 0.00000037 & 0.00001906 & -0.00594749 & 149472.67411175 & 0.16047689 & -0.12534081 \\ \hline Venus & 0.72333566 & 0.00677672 & 3.39467605 & 181.97909950 & 131.60246718 & 76.67984255 \\ & 0.00000390 & -0.00004107 & -0.00078890 & 58517.81538729 & 0.00268329 & -0.27769418 \\ \hline EM Bary & 1.00000261 & 0.01671123 & -0.00001531 & 100.46457166 & 102.93768193 & 0.0 \\ & 0.00000562 & -0.00004392 & -0.01294668 & 35999.37244981 & 0.32327364 & 0.0 \\ \hline Mars & 1.52371034 & 0.09339410 & 1.84969142 & -4.55343205 & -23.94362959 & 49.55953891 \\ & 0.00001847 & 0.00007882 & -0.00813131 & 19140.30268499 & 0.44441088 & -0.29257343 \\ \hline Jupiter & 5.20288700 & 0.04838624 & 1.30439695 & 34.39644051 & 14.72847983 & 100.47390909 \\ & -0.00011607 & -0.00013253 & -0.00183714 & 3034.74612775 & 0.21252668 & 0.20469106 \\ \hline Saturn & 9.53667594 & 0.05386179 & 2.48599187 & 49.95424423 & 92.59887831 & 113.66242448 \\ & -0.00125060 & -0.00050991 & 0.00193609 & 1222.49362201 & -0.41897216 & -0.28867794 \\ \hline Uranus & 19.18916464 & 0.04725744 & 0.77263783 & 313.23810451 & 170.95427630 & 74.01692503 \\ & -0.00196176 & -0.00004397 & -0.00242939 & 428.48202785 & 0.40805281 & 0.04240589 \\ \hline Neptune & 30.06992276 & 0.00859048 & 1.77004347 & -55.12002969 & 44.96476227 & 131.78422574 \\ & 0.00026291 & 0.00005105 & 0.00035372 & 218.45945325 & -0.32241464 & -0.00508664 \\ \end{tabular} \end{center} \end{table} Всички колони със стойности спрямо градуси ще бъдат преобразувани спрямо радиани, използвайки формулата $\rad = \deg \cdot \frac{\pi}{180}$ (с точност до 8 знака след десетичния знак). Ще сменим наименованията на елементите: $\mathrm{a} \rightarrow a_{2000}$, $\mathrm{e} \rightarrow e_{2000}$, $\mathrm{I} \rightarrow i_{2000}$, $\mathrm{L} \rightarrow \lambda_{2000}$, \break $\mathrm{long.peri.} \rightarrow \tilde{\omega}_{2000}$ и $\mathrm{long.node.} \rightarrow \theta_{2000}$. Ще бележим измененията с $\Delta$. Получаваме: \begin{table}[!h] \begin{center} \begin{tabular}{c||c|c|c|c|c|c} & $a_{2000}, \Delta a$ & $e_{2000}, \Delta e$ & $i_{2000}, \Delta i$ & $\lambda_{2000}, \Delta \lambda$ & $\tilde{\omega}_{2000}, \Delta \tilde{\omega}$ & $\theta_{2000}, \Delta \theta$ \\ & au, au/Cy & rad, rad/Cy & rad, rad/Cy & rad, rad/Cy & rad, rad/Cy & rad, rad/Cy \\ \hline\hline Mercury & 0.38709927 & 0.20563593 & 0.12225995 & 4.40259868 & 1.35189358 & 0.84353100 \\ & 0.00000037 & 0.00001906 & -0.00010380 & 2608.79030501 & 0.00280085 & -0.00218761 \\ \hline Venus & 0.72333566 & 0.00677672 & 0.05924827 & 3.17613446 & 2.29689636 & 1.33831572 \\ & 0.00000390 & -0.00004107 & -0.00001377 & 1021.32854958 & 0.00004683 & -0.00484668 \\ \hline EM Bary & 1.00000261 & 0.01671123 & -0.00000027 & 1.75343756 & 1.79660147 & 0.0 \\ & 0.00000562 & -0.00004392 & -0.00022596 & 628.30757790 & 0.00564219 & 0.0 \\ \hline Mars & 1.52371034 & 0.09339410 & 0.03228321 & -0.07947238 & -0.41789517 & 0.86497713 \\ & 0.00001847 & 0.00007882 & -0.00014192 & 334.06130168 & 0.00775643 & -0.00510637 \\ \hline Jupiter & 5.20288700 & 0.04838624 & 0.02276602 & 0.60033114 & 0.25706047 & 1.75360053 \\ & -0.00011607 & -0.00013253 & -0.00003206 & 52.96631189 & 0.00370929 & 0.00357253 \\ \hline Saturn & 9.53667594 & 0.05386179 & 0.04338874 & 0.87186604 & 1.61615531 & 1.98378354 \\ & -0.00125060 & -0.00050991 & 0.00003379 & 21.33653879 & -0.00731244 & -0.00503838 \\ \hline Uranus & 19.18916464 & 0.04725744 & 0.01348507 & 5.46703627 & 2.98371499 & 1.29183904 \\ & -0.00196176 & -0.00004397 & -0.00004240 & 7.47842217 & 0.00712187 & 0.00074012 \\ \hline Neptune & 30.06992276 & 0.00859048 & 0.03089309 & -0.96202600 & 0.78478315 & 2.30006864 \\ & 0.00026291 & 0.00005105 & 0.00000617 & 3.81283674 & -0.00562720 & -0.00008878 \\ \end{tabular} \end{center} \end{table} Без долен индекс и $\Delta$ ще обозначаваме величина след време $t$ от началото на епохата. Тоест, обозначаваме $a = a_{2000} + t \cdot \Delta a$, \, $e = e_{2000} + t \cdot \Delta e$, ... Финално, ще изчислим рожденната ми дата във векове. Аз съм роден на 27 Април 2002 година. Това е 117-ят ден в годината, което е $\frac{117}{365.25} \approx 0.3203285421$ Юлиански години. Общо излиза $2.3203285421$ Юлиански години след началото на епохата, което са $\frac{2.3203285421}{100} = 0.023203285421$ Юлиански века. \paragraph{1. Решение} Ще напишем програма на езика \texttt{R}, която намира координатите и скоростите на планетите в деня на моето раждане. Входните данни на програмата ще бъдат орбиталните елементи в началото на епохата и рожденната ми дата (в Юлиански векове от началото на епохата). Дефинираме функция \texttt{atDay}, която пресмята новата стойност на величина след време $t$: \begin{lstlisting}[language=R] atDay <- function(epochConstant, rate, t) { epochConstant + t * rate } \end{lstlisting} От теорема знаем, че векторът координат $r$ и векторът скорост $v$ се пресмятат чрез: \begin{align*} r = \begin{pmatrix} x \\ y \\ z \\ \end{pmatrix} = Q \begin{pmatrix} a(\cos u - e) \\ a \sqrt{1 - e^2} \sin u \\ 0 \end{pmatrix} \qquad v = \begin{pmatrix} v_x \\ v_y \\ v_z \\ \end{pmatrix} = Q \begin{pmatrix} -\sin u \\ \sqrt{1 - e^2} \cos u \\ 0 \end{pmatrix} \frac{an}{1 - e \cos u} \end{align*} Където $Q$ е матрицата на въртения около нулата: \begin{align*} Q = \begin{pmatrix} \cos\theta & -\sin\theta & 0 \\ \sin\theta & \cos\theta & 0 \\ 0 & 0 & 1 \\ \end{pmatrix} \begin{pmatrix} 1 & 0 & 0 \\ 0 & \cos i & -\sin i \\ 0 & \sin i & \cos i \\ \end{pmatrix} \begin{pmatrix} \cos g & -\sin g & 0 \\ \sin g & \cos g & 0 \\ 0 & 0 & 1 \\ \end{pmatrix} \end{align*} Дефинираме помощни функции \texttt{m3x3} и \texttt{p}. Първата приема редица елементи и ги преобразува в $3 \times 3$ матрица, втората изкарва подаденото на екрана, като всички числа са с 10 цифри прецизност: \begin{lstlisting}[language=R] m3x3 <- function(...) { matrix(c(...), nrow = 3, ncol = 3, byrow = TRUE) } p <- function(x) { formatC(x, format="f", digits = 10) } \end{lstlisting} Дефинираме функция \texttt{Q}, която пресмята матрица на въртене по подадени стойности $\theta$, $i$, $g$ (операторът \texttt{\%*\%} умножава матрици и вектори): \begin{lstlisting}[language=R] Q <- function(theta, i, g) { m3x3( cos(theta), -sin(theta), 0, sin(theta), cos(theta), 0, 0, 0, 1 ) %*% m3x3( 1, 0, 0, 0, cos(i), -sin(i), 0, sin(i), cos(i) ) %*% m3x3( cos(g), -sin(g), 0, sin(g), cos(g), 0, 0, 0, 1 ) } \end{lstlisting} Дефинираме функция \texttt{coordinates}, която изчислява стойността на $r$ вектора по горната формула: \begin{lstlisting}[language=R] coordinates <- function(theta, i, g, a, u, e) { Q(theta, i, g) %*% c( a * (cos(u) - e), a * sqrt(1 - e^2) * sin(u), 0 ) } \end{lstlisting} Дефинираме функция \texttt{speed}, която изчислява стойността на $v$ вектора по горната формула: \begin{lstlisting}[language=R] speed <- function(theta, i, g, a, u, e, n) { Q(theta, i, g) %*% c( -sin(u), sqrt(1 - e^2) * cos(u), 0 ) * a * n / (1 - e * cos(u)) } \end{lstlisting} По дефиниция знаем, че $\tilde{\omega} = g + \theta$ и $\lambda = l + g + \theta = l + \tilde{\omega}$. То, можем да намерим стойностите на $g$ и $l$ чрез формулите: \begin{align*} g = \tilde{\omega} - \theta \qquad l = \lambda - \tilde{\omega} \end{align*} Имплементираме ги във функции \texttt{g} и \texttt{l}: \begin{lstlisting}[language=R] g <- function(omega, theta) { omega - theta } l <- function(lambda, omega) { lambda - omega } \end{lstlisting} Eдин от начините да дефинираме средното движение $n$ е чрез средната аномалия $l$: \begin{align*} l = l_0 + n (t - t_0) \qquad n = \frac{l - l_0}{t - t_0} \end{align*} В нашата задача, $l_0 = l_{2000} = \lambda_{2000} - \tilde{\omega}_{2000}$. Нека да зачетем $t_0$ за началото на епохата и $t$ за време от това начало (във векове), тоест свеждаме $t - t_0$ до $t$. Тогава намираме, че $n$ е изменението на средната аномалия $l$: \begin{align*} n = \frac{l - l_{2000}}{t} = \frac{(l_{2000} + t \cdot \Delta l) - l_{2000}}{t} = \frac{t \cdot \Delta l}{t} = \Delta l \end{align*} Ще покажем, че изменението на средната аномалия е измението на $\lambda$ минус изменението на $\tilde{\omega}$: \begin{align*} \Delta l & = n = \frac{l - l_{2000}}{t} = \frac{(\lambda - \tilde{\omega}) - (\lambda_{2000} - \tilde{\omega}_{2000})}{t} = \frac{(\lambda_{2000} + t \cdot \Delta \lambda) - (\tilde{\omega}_{2000} + t \cdot \Delta \tilde{\omega}) - \lambda_{2000} + \tilde{\omega}_{2000}}{t} \\ & = \frac{t \cdot \Delta \lambda - t \cdot \Delta \tilde{\omega}}{t} = \Delta \lambda - \Delta \tilde{\omega} \end{align*} Под тази формула, $n$ е в единици радиани за век, обаче на нас ни трябва безмерната ѝ стойност. Следователно, ще разделим $n$ на $100$ и $2 \pi$. Имплементираме формулата във функция \texttt{n}: \begin{lstlisting}[language=R] n <- function(dlambda, domega) { (dlambda - domega) / 100 / (2*pi) } \end{lstlisting} Ще дефинираме стойността $u$ използвайки вторият метод за решаване на уравнението на Кеплер: \begin{align*} u = l + e \cdot \sin(l + e \cdot \sin (l + e \cdot \sin( \cdots ))) \end{align*} Нуждаем се от 8 знака прецизност, затова ще използваме следното приближение: \begin{align*} u = l + e \cdot \sin(l + e \cdot \sin (l + e \cdot \sin(l + e \cdot \sin(l + e \cdot \sin(l + e \cdot \sin(l + e \cdot \sin(l + e \cdot \sin(l)))))))) \end{align*} Дефинираме функция \texttt{u}, която пресмята стойността на $u$ чрез рекурсивен алгоритъм: \begin{lstlisting}[language=R] u <- function(l, e, depth=8) { if (depth <= 0) { l } else { l + e * sin(u(l, e, depth - 1)) } } \end{lstlisting} Дължината на векторите $r$ и $v$ са част от решението. Дефинираме функция \texttt{length}, която приема 3-мерен вектор и връща евклидовата му дължина: \begin{lstlisting}[language=R] length <- function(vec) { sqrt(vec[1]^2 + vec[2]^2 + vec[3]^2) } \end{lstlisting} \break \hfill Вече сме готови да имплементираме кода, който обединява всичко това. Дефинираме функция \texttt{planetData}, която приема датата $t$ във векове от началото на епохата, орбиталните елементи в началото на епохата и тяхното изменение. Функцията пресмята координатите $x$, $y$, $z$; радиуса $|r|$; покомпонентните скорости $v_x$, $v_y$, $v_z$ и скоростта на планетата $v$: \begin{lstlisting}[language=R] planetData <- function(t, a2000, da, e2000, de, i2000, di, lambda2000, dlambda, omega2000, domega, theta2000, dtheta) { a <- atDay(a2000, da, t) e <- atDay(e2000, de, t) i <- atDay(i2000, di, t) lambda <- atDay(lambda2000, dlambda, t) omega <- atDay(omega2000, domega, t) theta <- atDay(theta2000, dtheta, t) g <- g(omega, theta) l <- l(lambda, omega) n <- n(dlambda, domega) u <- u(l, e) r <- coordinates(theta, i, g, a, u, e) v <- speed(theta, i, g, a, u, e, n) c(r, length(r), v, length(v)) } \end{lstlisting} Финално дефинираме функцията \texttt{solarSystemAtBirth}, която извиква \texttt{planetData} за всяка планета с подадени константи: \begin{lstlisting}[language=R] solarSystemAtBirth <- function() { birthDate <- 0.023203285421 rbind( planetData(birthDate, # Меркурии 0.38709927, 0.00000037, 0.20563593, 0.00001906, 0.12225995, -0.00010380, 4.40259868, 2608.79030501, 1.35189358, 0.00280085, 0.84353100, -0.00218761), planetData(birthDate, # Венера 0.72333566, 0.00000390, 0.00677672, -0.00004107, 0.05924827, -0.00001377, 3.17613446, 1021.32854958, 2.29689636, 0.00004683, 1.33831572, -0.00484668), planetData(birthDate, # EM Bary 1.00000261, 0.00000562, 0.01671123, -0.00004392, -0.00000027, -0.00022596, 1.75343756, 628.30757790, 1.79660147, 0.00564219, 0.0, 0.0), planetData(birthDate, # Марс 1.52371034, 0.00001847, 0.09339410, 0.00007882, 0.03228321, -0.00014192, -0.07947238, 334.06130168, -0.41789517, 0.00775643, 0.86497713, -0.00510637), planetData(birthDate, # Юпитер 5.20288700, -0.00011607, 0.04838624, -0.00013253, 0.02276602, -0.00003206, 0.60033114, 52.96631189, 0.25706047, 0.00370929, 1.75360053, 0.00357253), planetData(birthDate, # Сатурн 9.53667594, -0.00125060, 0.05386179, -0.00050991, 0.04338874, 0.00003379, 0.87186604, 21.33653879, 1.61615531, -0.00731244, 1.98378354, -0.00503838), planetData(birthDate, # Уран 19.18916464, -0.00196176, 0.04725744, -0.00004397, 0.01348507, -0.00004240, 5.46703627, 7.47842217, 2.98371499, 0.00712187, 1.29183904, 0.00074012), planetData(birthDate, # Нептун 30.06992276, 0.00026291, 0.00859048, 0.00005105, 0.03089309, 0.00000617, -0.96202600, 3.81283674, 0.78478315, -0.00562720, 2.30006864, -0.00008878) ) } \end{lstlisting} Извикваме функцията: \begin{lstlisting}[language=R] p(solarSystemAtBirth()) \end{lstlisting} И получаваме следния резултат: \newpage \hfill \begin{table}[!h] \begin{center}\begin{large} \begin{tabular}{c||c|c|c|c} & $x$ & $y$ & $z$ & $|r|$ \\ \hline\hline Mercury & -0.2569892522 & 0.2164260053 & 0.0412666944 & 0.3385064720 \\ Venus & -0.1166445847 & 0.7094677286 & 0.0164331060 & 0.7191804113 \\ EM Bary & -0.7984835049 & -0.6131132909 & 0.0000033801 & 1.0067193328 \\ Mars & 0.0101799533 & 1.5686596135 & 0.0326113752 & 1.5690315855 \\ Jupiter & -1.8132518650 & 4.8899907784 & 0.0203482204 & 5.2153912785 \\ Saturn & 2.0761255995 & 8.7959937429 & -0.2357956393 & 9.0407634088 \\ Uranus & 16.4964276183 & -11.3052260072 & -0.2558295226 & 20.0001426991 \\ Neptune & 18.9271242195 & -23.3976888193 & 0.0456520201 & 30.0946832084 \\ \end{tabular} \end{large}\end{center} \end{table} \begin{table}[!h] \begin{center}\begin{large} \begin{tabular}{c||c|c|c|c} & $v_x$ & $v_y$ & $v_z$ & $|v|$ \\ \hline\hline Mercury & -1.3864619581 & -1.1839190472 & 0.0305314033 & 1.8234234939 \\ Venus & -1.1642070756 & -0.1972100811 & 0.0645017151 & 1.1825524945 \\ EM Bary & 0.5928042970 & -0.7969898554 & 0.0000043938 & 0.9932823185 \\ Mars & -0.7825575318 & 0.0743119103 & 0.0207870540 & 0.7863527531 \\ Jupiter & -0.4169723661 & -0.1321091661 & 0.0098828616 & 0.4375116648 \\ Saturn & -0.3329338533 & 0.0739849818 & 0.0119514504 & 0.3412646559 \\ Uranus & 0.1274263738 & 0.1777740834 & -0.0009918592 & 0.2187283458 \\ Neptune & 0.1409109730 & 0.1159848250 & -0.0056353427 & 0.1825928231 \\ \end{tabular} \end{large}\end{center} \end{table} \paragraph{2. Решение} Ще напишем програма на езика \texttt{R}, която намира елементите на Делоне и Поанкаре (от първи и втори вид) в деня на моето раждане. За целта се нуждаем и от масите на всички планети. В допълнението към астрономическия алманах, глава 8, страници 19-20 \cite{mass}, намираме: \begin{table}[!h] \begin{center} \begin{tabular}{c|c|c} GMs/GM1 & 6023600.0 & $mass^{−1}$ for Mercury \\ GMs/GM2 & 408523.71 & $mass^{−1}$ for Venus \\ GMs/GMB & 328900.5614 & $mass^{−1}$ for EM-Bary \\ GMs/GM4 & 3098708.0 & $mass^{−1}$ for Mars \\ GMs/GM5 & 1047.3486 & $mass^{−1}$ for Jupiter \\ GMs/GM6 & 3497.898 & $mass^{−1}$ for Saturn \\ GMs/GM7 & 22902.98 & $mass^{−1}$ for Uranus \\ GMs/GM8 & 19412.24 & $mass^{−1}$ for Neptune \\ \end{tabular} \end{center} \end{table} Ще се нуждаем от $\mu = \frac{m_\mathrm{планета}}{m_\mathrm{слънце}}$. Преобразуваме таблицата по формулата $\mu = \frac{1}{\mathrm{GMs/GMX}}$: \begin{table}[!h] \begin{center} \begin{tabular}{c||c} & $\mu$ \\ \hline\hline Mercury & 0.00000016601367952719 \\ Venus & 0.00000244783833966454 \\ EM Bary & 0.0000030404326333266 \\ Mars & 0.00000032271514450539 \\ Jupiter & 0.00095479193842432214 \\ Saturn & 0.0002858859806661029 \\ Uranus & 0.00004366244043351564 \\ Neptune & 0.00005151389020535497 \\ \end{tabular} \end{center} \end{table} Входните данни на програмата ще бъдат орбиталните елементи в началото на епохата, рожденната ми дата (в Юлиански векове от началото на епохата) и съотношенията на масите $\mu$. \newpage \hfill \hfill Първо ще разгледаме елементите на Делоне: $L$, $G$, $\Theta$, $l$, $g$, $\theta$. Искаме да обновим величините от началото на епохата спрямо време $t$ и да изкарваме на екрана стойности с 8 цифри прецизност. Копираме функциите \texttt{atDay} и \texttt{p} от 1 задача, като увеличаваме прецизността на 21 цифри след десетичната запетая: \begin{lstlisting}[language=R] atDay <- function(epochConstant, rate, t) { epochConstant + t * rate } p <- function(x) { formatC(x, format="f", digits = 21) } \end{lstlisting} Величината $\theta$ е орбитален елемент. За $l$ и $g$ копираме функциите \texttt{l} и \texttt{g} от задача 1: \begin{lstlisting}[language=R] l <- function(lambda, omega) { lambda - omega } g <- function(omega, theta) { omega - theta } \end{lstlisting} За величините $L$, $G$ и $\Theta$ знаем следните формули по теорема: \begin{align*} L = m \sqrt{\gamma}\sqrt{a} \qquad G = L \sqrt{1 - e^2} \qquad \Theta = G \cos i \end{align*} Където $a$, $e$ и $i$ са орбитални елементи. За останалите величини ($m$ и $\gamma$) имаме формулите: \begin{align*} \gamma = 1 + \mu \qquad m = \frac{\mu}{\gamma} = \frac{\mu}{1 + \mu} \end{align*} Имплементираме функции по формулите за всяка велична: \begin{lstlisting}[language=R] gamma <- function(mu) { 1 + mu } m <- function(mu) { mu / (1 + mu) } L <- function(m, gamma, a) { m * sqrt(gamma) * sqrt(a) } G <- function(L, e) { L * sqrt(1 - e^2) } Theta <- function(G, i) { G * cos(i) } \end{lstlisting} Обединявайки всичко, получаваме следната функция за намиране на елементите на Делоне по подадени време $t$, съотношение на маси $\mu$, орбитални елементи в началото на епохата и тяхното изменение за век: \begin{lstlisting}[language=R] delaunay <- function(t, mu, a2000, da, e2000, de, i2000, di, lambda2000, dlambda, omega2000, domega, theta2000, dtheta) { a <- atDay(a2000, da, t) e <- atDay(e2000, de, t) i <- atDay(i2000, di, t) lambda <- atDay(lambda2000, dlambda, t) omega <- atDay(omega2000, domega, t) theta <- atDay(theta2000, dtheta, t) m <- m(mu) gamma <- gamma(mu) L <- L(m, gamma, a) G <- G(L, e) Theta <- Theta(G, i) l <- l(lambda, omega) g <- g(omega, theta) c(L, G, Theta, l, g, theta) } \end{lstlisting} \newpage \hfill Дефинираме функцията \texttt{delaunayAtBirth}, която връща всички елементи на Делоне за всички планети с подадени константи: \begin{lstlisting}[language=R] delaunayAtBirth <- function() { birthDate <- 0.023203285421 rbind( delaunay(birthDate, 0.00000016601367952719, # Меркурии 0.38709927, 0.00000037, 0.20563593, 0.00001906, 0.12225995, -0.00010380, 4.40259868, 2608.79030501, 1.35189358, 0.00280085, 0.84353100, -0.00218761), delaunay(birthDate, 0.00000244783833966454, # Венера 0.72333566, 0.00000390, 0.00677672, -0.00004107, 0.05924827, -0.00001377, 3.17613446, 1021.32854958, 2.29689636, 0.00004683, 1.33831572, -0.00484668), delaunay(birthDate, 0.0000030404326333266, # EM Bary 1.00000261, 0.00000562, 0.01671123, -0.00004392, -0.00000027, -0.00022596, 1.75343756, 628.30757790, 1.79660147, 0.00564219, 0.0, 0.0), delaunay(birthDate, 0.00000032271514450539, # Марс 1.52371034, 0.00001847, 0.09339410, 0.00007882, 0.03228321, -0.00014192, -0.07947238, 334.06130168, -0.41789517, 0.00775643, 0.86497713, -0.00510637), delaunay(birthDate, 0.00095479193842432214, # Юпитер 5.20288700, -0.00011607, 0.04838624, -0.00013253, 0.02276602, -0.00003206, 0.60033114, 52.96631189, 0.25706047, 0.00370929, 1.75360053, 0.00357253), delaunay(birthDate, 0.0002858859806661029, # Сатурн 9.53667594, -0.00125060, 0.05386179, -0.00050991, 0.04338874, 0.00003379, 0.87186604, 21.33653879, 1.61615531, -0.00731244, 1.98378354, -0.00503838), delaunay(birthDate, 0.00004366244043351564, # Уран 19.18916464, -0.00196176, 0.04725744, -0.00004397, 0.01348507, -0.00004240, 5.46703627, 7.47842217, 2.98371499, 0.00712187, 1.29183904, 0.00074012), delaunay(birthDate, 0.00005151389020535497, # Нептун 30.06992276, 0.00026291, 0.00859048, 0.00005105, 0.03089309, 0.00000617, -0.96202600, 3.81283674, 0.78478315, -0.00562720, 2.30006864, -0.00008878) ) } \end{lstlisting} \hfill \hfill Елементите на Поанкаре от първи вид са следната редица стойности по дефиниция: \begin{align*} L, \quad L - G, \quad G - \Theta, \quad l + g + \theta, \quad - g - \theta, \quad - \theta \end{align*} Изполвайки вече дефинираната \texttt{delaunayAtBirth}, можем итеративно да конвертираме всеки ред с елементи на Делоне към елементите на Поанкаре: \begin{lstlisting}[language=R] poincare1 <- function() { data <- delaunayAtBirth() result <- c() for (planet in 1:8) { L <- data[planet,1]; G <- data[planet,2]; Theta <- data[planet,3] l <- data[planet,4]; g <- data[planet,5]; theta <- data[planet,6] result <- rbind( result, c(L, L - G, G - Theta, l + g + theta, - g - theta, - theta) ) } result } \end{lstlisting} Елементите на Поанкаре от втори вид са следната редица стойности по дефиниция: \begin{align*} \begin{matrix} L, & \xi = \sqrt{2(L - G)} \cos(g + \theta), & p = \sqrt{2(G - \Theta)} \cos(\theta), \\ l + g + \theta, & \eta = - \sqrt{2(L - G)} \sin(g + \theta), & q = - \sqrt{2(G - \Theta)} \sin(\theta) \end{matrix} \end{align*} \newpage Напълно аналогично на \texttt{poincare1}, дефинираме \texttt{poincare2}: \begin{lstlisting}[language=R] poincare2 <- function() { data <- delaunayAtBirth() result <- c() for (planet in 1:8) { L <- data[planet,1]; G <- data[planet,2]; Theta <- data[planet,3] l <- data[planet,4]; g <- data[planet,5]; theta <- data[planet,6] result <- rbind( result, c(L, sqrt(2*(L-G)) * cos(g + theta), sqrt(2*(G-Theta)) * cos(theta), l + g + theta, -sqrt(2*(L-G)) * sin(g + theta), -sqrt(2*(G-Theta)) * sin(theta)) ) } result } \end{lstlisting} \hfill Извикваме всяка от тези функции: \begin{lstlisting}[language=R] p(delaunayAtBirth()) p(poincare1()) p(poincare2()) \end{lstlisting} И получаваме следните резултати: \hfill \hfill \hfill \hfill \hfill \hfill \begin{table}[!h] \begin{center} \begin{tabular}{c|c|c|c} Делоне & $L$ & $G$ & $\Theta$ \\ \hline\hline Mercury & 0.000000103289225349226 & 0.000000101081776576878 & 0.000000100327287123918 \\ Venus & 0.000002081863094111913 & 0.000002081815303335347 & 0.000002078162453508473 \\ EM Bary & 0.000003040432177217194 & 0.000003040007655895516 & 0.000003040007655849318 \\ Mars & 0.000000398355249808993 & 0.000000396614057982431 & 0.000000396407441380291 \\ Jupiter & 0.002176825400448131106 & 0.002174276008944916582 & 0.002173712615531286918 \\ Saturn & 0.000882731532142888770 & 0.000881450722537746524 & 0.000880621120751847700 \\ Uranus & 0.000191260833739243018 & 0.000191047155475361307 & 0.000191029787586196811 \\ Neptune & 0.000282474574098807182 & 0.000282464148236046821 & 0.000282329368214719945 \\ \end{tabular} \end{center} \end{table} \begin{table}[!h] \begin{center} \begin{tabular}{c|c|c|c} Делоне & $l$ & $g$ & $\theta$ \\ \hline\hline Mercury & 63.583146161762698511666 & 0.508478328661191358862 & 0.843480240260780167461 \\ Venus & 24.577414857910831358367 & 0.958694185509240348608 & 1.338203261100615826251 \\ EM Bary & 14.535505234845919275699 & 1.796732387344969428256 & 0.000000000000000000000 \\ Mars & 8.089562546332686565620 & -1.282573840780286689167 & 0.864858645439424700285 \\ Jupiter & 1.572177054766796810981 & -1.496536886718685899211 & 1.753683424433265081177 \\ Saturn & -0.249041797926948138553 & -0.367680995663244303273 & 1.983666633030800374726 \\ Uranus & 2.656679993326903410633 & 1.692024027566735533412 & 1.291856213215605730937 \\ Neptune & -1.658208241330383803103 & -1.515413999540041523062 & 2.300066580012320560655 \\ \end{tabular} \end{center} \end{table} \begin{table}[!h] \begin{center} \begin{tabular}{c|c|c|c} Поанкаре 1 вид & $L$ & $L - G$ & $G - \Theta$ \\ \hline\hline Mercury & 0.000000103289225349226 & 0.000000002207448772348 & 0.000000000754489452959 \\ Venus & 0.000002081863094111913 & 0.000000000047790776567 & 0.000000003652849826874 \\ EM Bary & 0.000003040432177217194 & 0.000000000424521321678 & 0.000000000000000046198 \\ Mars & 0.000000398355249808993 & 0.000000001741191826562 & 0.000000000206616602140 \\ Jupiter & 0.002176825400448131106 & 0.000002549391503214524 & 0.000000563393413629664 \\ Saturn & 0.000882731532142888770 & 0.000001280809605142246 & 0.000000829601785898824 \\ Uranus & 0.000191260833739243018 & 0.000000213678263881711 & 0.000000017367889164495 \\ Neptune & 0.000282474574098807182 & 0.000000010425862760361 & 0.000000134780021326876 \\ \end{tabular} \end{center} \end{table} \begin{table}[!h] \begin{center} \begin{tabular}{c|c|c|c} Поанкаре 1 вид & $l + g + \theta$ & $-g - \theta$ & $-\theta$ \\ \hline\hline Mercury & 64.935104730684670926166 & -1.351958568921971526322 & -0.843480240260780167461 \\ Venus & 26.874312304520689309584 & -2.296897446609856174859 & -1.338203261100615826251 \\ EM Bary & 16.332237622190888259865 & -1.796732387344969428256 & -0.000000000000000000000 \\ Mars & 7.671847350991824576738 & 0.417715195340861988882 & -0.864858645439424700285 \\ Jupiter & 1.829323592481375992946 & -0.257146537714579181966 & -1.753683424433265081177 \\ Saturn & 1.366943839440607932900 & -1.615985637367556071453 & -1.983666633030800374726 \\ Uranus & 5.640560234109244674983 & -2.983880240782341264349 & -1.291856213215605730937 \\ Neptune & -0.873555660858104765509 & -0.784652580472279037593 & -2.300066580012320560655 \\ \end{tabular} \end{center} \end{table} \begin{table}[!h] \begin{center} \begin{tabular}{c|c|c|c} Поанкаре 2 вид & $L$ & $\xi$ & $p$ \\ \hline\hline Mercury & 0.000000103289225349226 & 0.000014424827544034668 & 0.000025827152487161843 \\ Venus & 0.000002081863094111913 & -0.000006491251753606220 & 0.000019701746607823352 \\ EM Bary & 0.000003040432177217194 & -0.000006527533112570199 & 0.000000009612286302040 \\ Mars & 0.000000398355249808993 & 0.000053937790938604713 & 0.000013187832795805843 \\ Jupiter & 0.002176825400448131106 & 0.002183802985649161126 & -0.000193054628745658468 \\ Saturn & 0.000882731532142888770 & -0.000072301145995106139 & -0.000516837678338016642 \\ Uranus & 0.000191260833739243018 & -0.000645611834774382142 & 0.000051316013457963606 \\ Neptune & 0.000282474574098807182 & 0.000102183215085003621 & -0.000345950770443006518 \\ \end{tabular} \end{center} \end{table} \begin{table}[!h] \begin{center} \begin{tabular}{c|c|c|c} Поанкаре 2 вид & $l + g + \theta$ & $\eta$ & $q$ \\ \hline\hline Mercury & 64.935104730684670926166 & -0.000064860017692109048 & -0.000029016152403852770 \\ Venus & 26.874312304520689309584 & -0.000007310622668714660 & -0.000083171755027473018 \\ EM Bary & 16.332237622190888259865 & -0.000028397780808021289 & -0.000000000000000000000 \\ Mars & 7.671847350991824576738 & 0.000023939472880300536 & -0.000015469785726697444 \\ Jupiter & 1.829323592481375992946 & -0.000574271300256992776 & -0.001043799184507826754 \\ Saturn & 1.366943839440607932900 & -0.001598872025701959104 & -0.001179865410141265995 \\ Uranus & 5.640560234109244674983 & -0.000102673689729541247 & -0.000179171552127487335 \\ Neptune & -0.873555660858104765509 & -0.000102030956457410801 & -0.000387140939560313784 \\ \end{tabular} \end{center} \end{table} \vfill \renewcommand{\refname}{Източници} \begin{thebibliography}{2} \bibitem{nasa} Solar System Dynamics. (2025, Май 1). \textit{Approximate Positions of the Planets.} Jet Propultion Laboratory, National Aeronautics and Space Administration. \url{https://ssd.jpl.nasa.gov/planets/approx_pos.html} \bibitem{mass} Williams, James G., \& Standish, Myles E. (2006). Orbital Ephemerides of the Sun, Moon, and Planets. In Seidelmann, P. K. (Eds.), \textit{Explanatory Supplement to the Astronomical Almanac}. \url{https://www.researchgate.net/publication/232203657_Orbital_Ephemerides_of_the_Sun_Moon_and_Planets} \end{thebibliography} \end{document}