Прикладная математика и механика. T. 87, Номер 5, 2023

Прикладная математика и механика, 2023, T. 87, № 5, стр. 742-756

Применение метода быстрых разложений к построению траектории движения тела переменной массы из его начального положения в заданное конечное в гравитационном поле

А. Д. Чернышов 1, *, М. И. Попов 2, **, В. В. Горяйнов 3, ***, О. Ю. Никифорова 1, ****

1 Воронежский государственный университет инженерных технологий
Воронеж, Россия

2 Воронежский государственный университет
Воронеж, Россия

3 Воронежский государственный технический университет
Воронеж, Россия

* E-mail: chernyshovad@mail.ru
** E-mail: mihail_semilov@mail.ru
*** E-mail: gorvit77@mail.ru
**** E-mail: niki22@mail.ru

Поступила в редакцию 28.09.2022
После доработки 04.08.2023
Принята к публикации 10.08.2023

Полный текст (PDF)

Аннотация

Приводится аналитическое решение задачи о движении космического корабля из начальной точки в конечный пункт за определенное время. Вначале используется метод быстрых синус-разложений. Рассматриваемая здесь космическая задача существенно нелинейная, что порождает необходимость в использовании методов тригонометрической интерполяции, которая по точности и простоте превосходит все известные интерполяции. При этом задача вычисления коэффициентов Фурье интегральными формулами заменяется на решение ортогональной интерполяционной системы. В этой связи рассматривается два случая на отрезке $\left[ {0,a} \right]$: универсальная интерполяция и тригонометрические синус- и косинус-интерполяции. Доказана теорема о быстром убывании коэффициентов разложений, получена компактная формула для вычисления коэффициентов интерполяции. Дается общая теория быстрых разложений. Показано, что в таком случае коэффициенты Фурье с ростом порядкового номера убывают значительно быстрее по сравнению с коэффициентами Фурье в классическом случае. Это свойство позволяет существенно сократить число учитываемых членов в ряде Фурье, существенно увеличить точность расчетов и уменьшить объем вычислений на ЭВМ. Проведен анализ полученных решений задачи движения космического корабля и предложено их сравнение с точным решением тестовой задачи. Приближенное решение по методу быстрых разложений вполне можно принимать за точное, так как используемые из справочников входные данные задачи имеют более высокую погрешность.

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

1. Введение. Расчет траектории движения космического аппарата в атмосфере Земли или других планет имеет огромное значение при подготовке полетов спутников и космических кораблей. Математическая модель движения космического корабля на пассивном и активном участках траектории представляет собой систему нелинейных дифференциальных уравнений. Интегрирование такой системы связано со значительными трудностями в связи с их иррациональной нелинейностью, так как в системе присутствуют дробные степени от неизвестных функций. Задачи, связанные с расчетом траектории движения космических кораблей, обычно решаются численными методами [1–3]. Наиболее успешным из них является метод “раздельного интегрирования” [4]. Полуаналитический метод с использованием вспомогательных эмпирических зависимостей для введенных целевых и базисных функций применяется в работе [5]. По существу такой подход является улучшенным численным методом. Таким образом, актуальным является разработка новых методов для расчета траекторий космических кораблей, способных с высокой точностью при минимальных временных затратах на ЭВМ определить решение в аналитическом виде. Этим условиям удовлетворяет метод быстрых разложений [6], для применимости которого требуется только условие гладкости рассматриваемых функций в заданной области до некоторого заданного порядка. Высокая эффективность метода быстрых разложений апробирована во многих публикациях. Например, решены нелинейные задачи для интегро-дифференциальных уравнений [7], для задач с криволинейными областями [8], для задач с неизвестной границей типа Стефана [9], задачи космического характера [10] и многие другие. Данная статья еще один пример эффективности метода быстрых разложений.

2. Постановка задачи. Космический корабль при работе его реактивных двигателей будем считать телом переменной массы. Запишем уравнения движения корабля в виде системы дифференциальных уравнений для его центра масс:

(2.1)
${\mathbf{\ddot {r}}} + \beta {\mathbf{\dot {r}}} + \frac{\alpha }{{{{r}^{3}}}}{\mathbf{r}} = {\mathbf{P}},$
где ${\mathbf{r}} = \left( {x\left( t \right),y\left( t \right),z\left( t \right)} \right)$ ∈ ${{C}^{6}}\left( {t \in \left[ { - a,a} \right]} \right)$, $r = \sqrt {{{x}^{2}} + {{y}^{2}} + {{z}^{2}}} $, ${\mathbf{P}} = \left( {{{P}_{x}},{{P}_{y}},{{P}_{z}}} \right)$. В систе-ме (2.1) $\alpha = gR_{{\text{З}}}^{2}$ – коэффициент притяжения к Земле единичной массы, ${{R}_{{\text{З}}}}$ – радиус Земли, $g$ – ускорение свободного падения, $\beta $ – коэффициент аэродинамического сопротивления, ${{P}_{x}},{{P}_{y}},{{P}_{z}}$ – составляющие ускорений реактивной и других возможных сил.

В качестве тестового примера, функции ${{P}_{x}},{{P}_{y}},{{P}_{z}}$ в правой части (2.1) подобраны так, чтобы можно было записать точное решение задачи в виде

(2.2)
$x{{{\kern 1pt} }_{*}}\, = a\cos \omega t,\quad y{{{\kern 1pt} }_{*}}\, = a\sin \omega t,\quad z{{{\kern 1pt} }_{*}}\, = {{R}_{{\text{З}}}} + wt$

Траектория (2.2) обеспечивается силами ${{P}_{x}},{{P}_{y}},{{P}_{z}}$:

(2.3)
$\begin{gathered} {{P}_{x}} = - a{{\omega }^{2}}\cos \omega t - \beta a\omega \sin \omega t + \frac{\alpha }{{r_{*}^{3}}}a\cos \omega t \\ {{P}_{y}} = - a{{\omega }^{2}}\sin \omega t + \beta a\omega \cos \omega t + \frac{\alpha }{{r_{*}^{3}}}a\sin \omega t,\quad {{P}_{z}} = \beta w + \frac{\alpha }{{r_{*}^{3}}}\left( {{{R}_{{\text{З}}}} + wt} \right) \\ {{r}_{*}} = \sqrt {x_{*}^{2} + y_{*}^{2} + z_{*}^{2}} \\ \end{gathered} $

Граничные условия для точного решения (2.2) из заданного начального положения $x\left( 0 \right),y\left( 0 \right),z\left( 0 \right)$ в конечное $x\left( {{{t}_{0}}} \right),y\left( {{{t}_{0}}} \right),z\left( {{{t}_{0}}} \right)$ примут вид:

(2.4)
$\begin{gathered} x\left( 0 \right) = a,\quad y\left( 0 \right) = 0,\quad z\left( 0 \right) = {{R}_{{\text{З}}}} \\ x\left( {{{t}_{0}}} \right) = a\cos \omega {{t}_{0}} = {{x}_{{{{t}_{0}}}}},\quad y\left( {{{t}_{0}}} \right) = a\sin \omega {{t}_{0}} = {{y}_{{{{t}_{0}}}}},\quad z\left( {{{t}_{0}}} \right) = {{R}_{{\text{З}}}} + w{{t}_{0}} = {{z}_{{{{t}_{0}}}}} \\ \end{gathered} $

Данная задача является важной и актуальной. Решение задачи (2.1), (2.4) будет получено новым аналитическим методом с высокой точностью при минимальных временных затратах на ЭВМ.

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

Будем рассматривать движение космического корабля как движение его центра масс из заданного начального положения в конечное. Поэтому в задаче (2.1), (2.4) неизвестными являются координаты центра масс космического корабля ${\mathbf{r}}\left( t \right)$ во время его полета. Составляющие ускорений сил ${\mathbf{P}}$ в (2.3) подобраны так, что корабль движется по выбранной винтовой линии (2.2). Имея точное решение (2.2) можно путем сравнения с приближенным, полученным ниже, вычислить абсолютную погрешность местоположения корабля, его скорость и ускорение, что позволит сравнить приближенное аналитическое решение с точным.

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

3. Основные характеристики метода быстрых разложений равномерных интерполяций. Основными недостатками всех известных аналитических интерполяций является их невысокая точность и недопустимость почленного дифференцирования, которое приводит к большой погрешности. Сюда относятся вариационные методы, методы наименьших квадратов, сплайн-методы, классические методы тригонометрической интерполяции и некоторые другие. В известных научных публикациях авторы пишут о возможности дифференцирования при использовании их метода. Это голословные утверждения, так как не приводятся строгие математические обоснования для возможности дифференцирования и не приводятся примеры, подтверждающие возможность дифференцирования. На самом деле во всех известных случаях дифференцирование недопустимо, поэтому авторы не приводят тестовые примеры. Лишь только метод дискретного преобразования Фурье позволяет дифференцирование до второго порядка. Формальное использование данного метода к уравнениям четвертого порядка привело к расходящимся рядам и неправильным результатам.

Метод быстрых разложений (тригонометрических интерполяций) на данный момент является единственным методом, эффективно решающим данную проблему. К настоящему времени разработано несколько разновидностей быстрых разложений: синус-, косинус-разложения [6], универсальный метод быстрых разложений [11].

3.1. Полные тригонометрические интерполяции. При рассмотрении многих прикладных задач удобно использовать полные тригонометрические ряды Фурье, например работа [12]. Пусть некоторая функция $f\left( x \right)$ непрерывна и интегрируема на отрезке $\left[ { - a,a} \right]$. Представим ее полной интерполяционной тригонометрической суммой

(3.1)
$f\left( x \right) = {{a}_{0}} + \sum\limits_{n = 1}^{N - 1} {{{a}_{n}}\cos \pi n\frac{x}{a} + {{b}_{n}}\sin \pi n\frac{x}{a}} ;\quad x \in \left[ { - a,a} \right],\quad f\left( x \right) \in L_{p}^{\alpha }$

Отрезок $\left[ { - a,a} \right]$ равномерно разобьем на $2N$ частей точками ${{x}_{j}}$, коэффициенты интерполяции ${{a}_{0}},{{a}_{n}},{{b}_{n}}$ вычислим из интерполяционной системы, полагая в (3.1) $x = {{x}_{j}}$, ${{x}_{j}} = \frac{{a \cdot j}}{N}$; $j = - N \div N - 1$:

(3.2)
$f\left( {{{x}_{j}}} \right) = {{a}_{0}} + \sum\limits_{n = 1}^{N - 1} {{{a}_{n}}\cos \pi n\frac{{{{x}_{j}}}}{a} + {{b}_{n}}\sin \pi n\frac{{{{x}_{j}}}}{a}} ;\quad j = - N \div N - 1$

Из сходимости ряда (3.1) следует, что в системе (3.2) уравнение при $j = 0$ следует опустить, так как оно выполняется тождественно. Поэтому в (3.2) имеем $2N - 1$ уравнений относительно такого же количества неизвестных ${{a}_{0}},{{a}_{n}},{{b}_{n}}$.

Задача о нахождении ${{a}_{0}},{{a}_{n}},{{b}_{n}}$ из системы (3.2) сильно упрощается, если доказать, что дискретная система функций

(3.3)
$\left( {1,\cos \pi n\frac{{{{x}_{j}}}}{a},\sin \pi n\frac{{{{x}_{j}}}}{a}} \right);\quad j = - N \div N - 1$
является ортогональной при суммировании по переменной $j = - N \div N - 1$ в том же смысле, как и непрерывная система $\left( {1,\cos \pi n\frac{x}{a},\sin \pi n\frac{x}{a}} \right)$; $x \in \left[ { - a,a} \right]$ является ортогональной в пространстве Гильберта. Это означает, что надо доказать справедливость следующих равенств

$1)\;\;\sum\limits_{j = - N}^{N - 1} {\cos \pi n\frac{{{{x}_{j}}}}{a}} = 0,\quad \sum\limits_{j = - N}^{N - 1} {\sin \pi n\frac{{{{x}_{j}}}}{a}} = 0;\quad j = - N \div N - 1$
(3.4)
$\begin{gathered} 2)\;\;\sum\limits_{j = - N}^{N - 1} {\cos \pi m\frac{{{{x}_{j}}}}{a}\cos \pi n\frac{{{{x}_{j}}}}{a}} = 0\quad {\text{при}}\quad m \ne n \\ 3)\;\;\sum\limits_{j = - N}^{N - 1} {\sin \pi m\frac{{{{x}_{j}}}}{a}\sin \pi n\frac{{{{x}_{j}}}}{a}} = 0\quad {\text{при}}\quad m \ne n \\ \end{gathered} $
$4)\;\;\sum\limits_{j = - N}^{N - 1} {\sin \pi m\frac{{{{x}_{j}}}}{a}\cos \pi n\frac{{{{x}_{j}}}}{a}} = 0\quad {\text{при}}\quad \forall \left( {m,n} \right) = 1 \div \left( {N - 1} \right)$

В последнем равенстве системы (3.4) не следует требовать выполнение неравенства $m \ne n$, так как в 4) из (3.4) равенство выполняется при любых $m,\,\,n$. Значение индекса $j = N$ в суммах (3.4) исключается в силу предполагаемой периодичности $f\left( x \right)$, представленной зависимостью (3.1), при ее продолжении на область вне отрезка $\left[ { - a,a} \right]$, где она практически не рассматривается.

Для доказательства первой формулы 1) из (3.4) воспользуемся тригонометрической формулой Эйлера для комплексного числа и представим левую часть данной формулы в виде

(3.5)
$\sum\limits_{j = - N}^{N - 1} {\cos \pi n\frac{{{{x}_{j}}}}{a}} = \frac{1}{2}\sum\limits_{j = - N}^{N - 1} {\exp i\pi n\frac{j}{N} + \exp \left( { - i\pi n\frac{j}{N}} \right)} = \frac{1}{2}\sum\limits_{j = - N}^{N - 1} {{{q}_{{nj}}} + {{q}_{{ - nj}}}} ,$
где $i$ – мнимая единица и использовано обозначение

$\exp \left( {i\pi n\frac{j}{N}} \right) = {{q}_{{nj}}},\quad \exp \left( { - i\pi n\frac{j}{N}} \right) = {{q}_{{ - nj}}}$

Тогда равенство (3.5) можно записать в виде

(3.6)
$\sum\limits_{j = - N}^{N - 1} {\cos \pi n\frac{{{{x}_{j}}}}{a}} = \frac{1}{2}\sum\limits_{j = - N}^{N - 1} {{{q}_{{nj}}} + {{q}_{{ - nj}}}} $

При помощи свойств геометрических прогрессий правую часть в (3.6) преобразуем следующим образом

(3.7)
$\begin{gathered} \sum\limits_{j = - N}^{N - 1} {{{q}_{{nj}}} + {{q}_{{ - nj}}}} = \left( {{{q}_{{ - nN}}} + {{q}_{{n\left( {1 - N} \right)}}} + {{q}_{{n\left( {2 - N} \right)}}} + ... + {{q}_{{n\left( {N - 1} \right)}}}} \right) + \\ + \;\left( {{{q}_{{nN}}} + {{q}_{{n\left( {N - 1} \right)}}} + {{q}_{{n\left( {N - 2} \right)}}} + } \right.\left. {... + {{q}_{{n\left( {1 - N} \right)}}}} \right) = {{q}_{{ - nN}}}\frac{{1 - {{q}_{{2Nn}}}}}{{1 - {{q}_{n}}}} + {{q}_{{nN}}}\frac{{1 - {{q}_{{ - 2Nn}}}}}{{1 - {{q}_{{ - n}}}}} \\ \end{gathered} $

Так как

$\begin{gathered} {{q}_{{ \pm 2Nn}}} = \exp \left( { \pm i2\pi n} \right) = \cos \left( { \pm 2\pi n} \right) + i\sin \left( { \pm 2\pi n} \right) = 1 \\ {{q}_{{ \pm n}}} = \exp \left( { \pm i\frac{{\pi n}}{N}} \right) \ne 1\quad {\text{при}}\quad n = 1 \div N - 1 \\ \end{gathered} $
то из (3.7) получаем доказательство первого равенства из (3.4). Аналогично доказываются и все остальные равенства из (3.4), при помощи тригонометрической формулы Эйлера и получающихся геометрических прогрессий.

Остается вычислить норму для дискретной системы (3.3). Вначале рассмотрим сумму квадратов косинусов

(3.8)
$\begin{gathered} \sum\limits_{j = - N}^{N - 1} {{{{\cos }}^{2}}\pi n\frac{{{{x}_{j}}}}{a}} = \frac{1}{2}\sum\limits_{j = - N}^{N - 1} {\left( {1 + \cos 2\pi n\frac{j}{N}} \right)} = N + \frac{1}{2}\sum\limits_{j = - N}^{N - 1} {\cos 2\pi n\frac{j}{N}} = \\ = N + \frac{1}{4}{{q}_{{ - 2nN}}}\frac{{1 - {{q}_{{4nN}}}}}{{1 - {{q}_{{2n}}}}} + \frac{1}{4}{{q}_{{2nN}}}\frac{{1 - {{q}_{{ - 4nN}}}}}{{1 - {{q}_{{ - 2n}}}}} \\ \end{gathered} $

Для преобразования равенства (3.8) воспользуемся выражениями

(3.9)
${{q}_{{ \pm 2nN}}} = \cos 2\pi n \pm i\sin 2\pi n = 1,\quad {{q}_{{ \pm 4nN}}} = \cos 4\pi n \pm i\sin 4\pi n = 1$

Кроме (3.9) при $n = 1 \div N - 1$ выполняются неравенства

(3.10)
${{q}_{{ \pm 2n}}} = \cos \frac{{2\pi n}}{N} + i\sin \frac{{2\pi n}}{N} \ne 1,\quad {{q}_{{ \pm 4nN}}} = \cos \frac{{4\pi n}}{N} + i\sin \frac{{4\pi n}}{N} \ne 1$

Для доказательства справедливости неравенств (3.10) положим обратное: пусть выполняются равенства

${{q}_{{ \pm 2n}}} = \cos \frac{{2\pi n}}{N} + i\sin \frac{{2\pi n}}{N} = 1,\quad {{q}_{{ \pm 4nN}}} = \cos \frac{{4\pi n}}{N} + i\sin \frac{{4\pi n}}{N} = 1$

Докажем, что данные равенства приведут к противоречиям. Эти равенства возможны, если ${n \mathord{\left/ {\vphantom {n N}} \right. \kern-0em} N} = {{k}_{0}}$ – какое-либо целое число, или для второго случая ${{2n} \mathord{\left/ {\vphantom {{2n} N}} \right. \kern-0em} N} = {{k}_{0}}$ – какое-либо целое число. Но при $n = 1 \div N - 1$ дробь ${n \mathord{\left/ {\vphantom {n N}} \right. \kern-0em} N}$ всегда остается правильной и потому не может равняться целому числу. Во втором случае имеем следующие возможные значение дроби

$\frac{{2n}}{N} = \frac{{2n}}{N} \div \frac{{2\left( {N - 1} \right)}}{N} = \frac{{2n}}{N} \div 2\left( {1 - \frac{1}{N}} \right)$

При $N \geqslant 3$ дробь ${2 \mathord{\left/ {\vphantom {2 N}} \right. \kern-0em} N}$ всегда остается правильной и потому ${{2n} \mathord{\left/ {\vphantom {{2n} N}} \right. \kern-0em} N}$ не может равняться целому числу. Полученное противоречие является обоснованием справедливости неравенств (3.10) при $N \geqslant 3$.

При помощи (3.9) и (3.10) из (3.8) получаем выражение для нормы косинусов

(3.11)
$\sum\limits_{j = - N}^{N - 1} {{{{\cos }}^{2}}\pi n\frac{{{{x}_{j}}}}{a}} = N$

Аналогично доказывается и выражение для нормы синусов

(3.12)
$\sum\limits_{j = - N}^{N - 1} {{{{\sin }}^{2}}\pi n\frac{{{{x}_{j}}}}{a}} = N$

Следует заметить, что нормы (3.11) и (3.12) не зависят от номера $n$. Ортогональные свойства (3.4) и выражения норм (3.11) и (3.12) позволяют получить решение интерполяционной системы (3.2) относительно ${{a}_{0}},{{a}_{n}},{{b}_{n}}$ в явном конечном виде. Для этого вначале левую и правую части (3.2) просуммируем по индексу $j$:

(3.13)
$\sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)} = 2N{{a}_{0}} + \sum\limits_{n = 1}^{N - 1} {{{a}_{n}}\sum\limits_{j = - N}^{N - 1} {\cos \pi n\frac{{{{x}_{j}}}}{a}} + {{b}_{n}}\sum\limits_{j = - N}^{N - 1} {\sin \pi n\frac{{{{x}_{j}}}}{a}} } $

Выше было доказано, что суммы косинусов и синусов по свойству 1) из (3.4) равны нулю, поэтому из (3.13) имеем

(3.14)
${{a}_{0}} = \frac{1}{{2N}}\sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)} {\kern 1pt} ,$

т.е. коэффициент ${{a}_{0}}$ равен среднеарифметическому значению $f\left( {{{x}_{j}}} \right)$ в точках равномерной интерполяции на отрезке $\left[ { - a,a} \right]$.

Для нахождения ${{a}_{n}}$ умножим левую и правую части (3.2) на $\cos \pi m{j \mathord{\left/ {\vphantom {j N}} \right. \kern-0em} N}$ и просуммируем по индексу $j$:

(3.15)
$\begin{gathered} \sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)\cos \pi m\frac{{{{x}_{j}}}}{a}} = \sum\limits_{j = - N}^{N - 1} {{{a}_{0}}\cos \pi m\frac{{{{x}_{j}}}}{a}} + \\ + \;\sum\limits_{n = 1}^{N - 1} {\sum\limits_{j = - N}^{N - 1} {{{a}_{n}}\cos \pi m\frac{{{{x}_{j}}}}{a}\cos \pi n\frac{{{{x}_{j}}}}{a} + {{b}_{n}}\sin \pi n\frac{{{{x}_{j}}}}{a}\cos \pi m\frac{{{{x}_{j}}}}{a}} } \\ \end{gathered} $

Для упрощения равенств (3.15) воспользуемся свойствами (3.4) и (3.11),

(3.16)
$\begin{gathered} \sum\limits_{j = - N}^{N - 1} {\cos \pi m\frac{{{{x}_{j}}}}{a}} = 0,\quad \sum\limits_{j = - N}^{N - 1} {\cos \pi m\frac{{{{x}_{j}}}}{a}\cos \pi n\frac{{{{x}_{j}}}}{a}} = 0\quad {\text{при}}\quad m \ne n,\quad \Rightarrow \\ \Rightarrow \quad \sum\limits_{j = - N}^{N - 1} {\cos \pi m\frac{{{{x}_{j}}}}{a}\cos \pi n\frac{{{{x}_{j}}}}{a}} = N\quad {\text{при}}\quad m = n \\ \end{gathered} $

С учетом (3.16) уравнение (3.15) принимает вид

(3.17)
$\sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)\cos \pi m\frac{{{{x}_{j}}}}{a}} = {{a}_{m}}N$

Отсюда находим в явном виде коэффициент ${{a}_{m}}$:

(3.18)
${{a}_{m}} = \frac{1}{N}\sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)\cos \pi m\frac{{{{x}_{j}}}}{a}} ;\quad m = 1 \div N - 1$

Аналогично после умножения левой и правой частей (3.2) на $\sin \pi m{j \mathord{\left/ {\vphantom {j N}} \right. \kern-0em} N}$ и суммирования по индексу $j$ подобными вычислениями найдем ${{b}_{m}}$:

(3.19)
${{b}_{m}} = \frac{1}{N}\sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)\sin \pi m\frac{{{{x}_{j}}}}{a}} ;\quad m = 1 \div N - 1$

Подставляя ${{a}_{0}},{{a}_{n}},{{b}_{n}}$ из (3.14), (3.18), (3.19) в (3.1), получим формулу для полной равномерной тригонометрической интерполяции $f\left( x \right)$ на отрезке $\left[ { - a,a} \right]$:

(3.20)
$\begin{gathered} f\left( x \right) = \frac{1}{{2N}}\sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)} + \frac{1}{N}\sum\limits_{n = 1}^{N - 1} {\sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)\cos \pi n\frac{{{{x}_{j}}}}{a}} \cos \pi n\frac{x}{a}} + \\ + \;\frac{1}{N}\sum\limits_{n = 1}^{N - 1} {\sum\limits_{j = - N}^{N - 1} {f\left( {{{x}_{j}}} \right)\sin \pi n\frac{{{{x}_{j}}}}{a}} \sin \pi n\frac{x}{a}} ;\quad f\left( x \right) \in L_{p}^{\alpha }\left( {x \in \left[ { - a,a} \right)} \right) \\ \end{gathered} $

3.2. Синус- и косинус-интерполяции. Рассмотрим случай, когда $f\left( x \right)$ антисимметричная функция на отрезке $\left[ { - a,a} \right]$, $f\left( x \right) = - f\left( { - x} \right)$. Такой случай часто используется при рассмотрении инженерных задач, когда материальное тело описывается одной координатой $x \in \left[ {0,a} \right]$. Тогда достаточно $f\left( x \right)$ рассматривать на отрезке $\left[ {0,a} \right]$, формулы (3.14), (3.18), (3.19) упрощаются и выражение (3.20) принимает вид синус разложения

(3.21)
$\begin{gathered} {{a}_{0}} = {{a}_{n}} = 0,\quad {{b}_{n}} = \frac{2}{N}\sum\limits_{j = 1}^{N - 1} {f\left( {{{x}_{j}}} \right)\sin \pi n\frac{{{{x}_{j}}}}{a}} ;\quad n = 1 \div N - 1 \\ f\left( x \right) = \frac{2}{N}\sum\limits_{n = 1}^{N - 1} {\left( {\sum\limits_{j = 1}^{N - 1} {f\left( {{{x}_{j}}} \right)\sin \pi n\frac{{{{x}_{j}}}}{a}} } \right)\sin \pi n\frac{x}{a}} ;\quad x \in \left[ {0,a} \right] \\ \end{gathered} $

Особенно важным следует отметить, что в выражении для $f\left( x \right)$ в (3.21) в сумме по $j$ количество слагаемых равно количеству коэффициентов ${{b}_{n}}$. Данное свойство будет использовано в дальнейшем, что позволит получить замкнутую алгебраическую систему относительно неизвестных.

Если же $f\left( x \right)$ продолжить на отрезок $\left[ { - a,0} \right]$ симметрично, т.е. положить $f\left( x \right) = f\left( { - x} \right)$, то формулы (3.14), (3.18), (3.19) соответственно упрощаются и выражение (3.20) принимает вид косинус разложения

(3.22)
$f\left( x \right) = \frac{1}{N}\sum\limits_{j = 0}^{N - 1} {f\left( {{{x}_{j}}} \right)} + \frac{2}{N}\sum\limits_{n = 1}^{N - 1} {\sum\limits_{j = 0}^{N - 1} {f\left( {{{x}_{j}}} \right)\cos \pi n\frac{{{{x}_{j}}}}{a}} \cos \pi n\frac{x}{a}} ;\quad f\left( x \right) \in L_{p}^{\alpha }\left( {x \in \left[ {0,a} \right]} \right)$

Интерполяционные формулы (3.20)–(3.22) простые, но они имеют существенные недостатки, так как данные суммы при больших $N$ в общем случае медленно сходятся и не допускают почленное дифференцирование. Тогда для достижения высокой точности потребуется большое количество слагаемых, что и является основной причиной редкого их использования в научных работах.

4. Решение задачи о расчете траектории тела переменной массы методом быстрой тригонометрической синус-интерполяции. Для решения некоторой задачи методом быстрых синус-разложений какую-либо неизвестную функцию $f\left( t \right)$ следует предварительно представить суммой специальной граничной функции ${{M}_{{2k}}}$ и ряда Фурье для разности $f\left( t \right) - {{M}_{{2k}}}\left( t \right)$ на заданном отрезке

(4.1)
$f\left( t \right) = {{M}_{{2k}}}\left( t \right) + \sum\limits_{m = 1}^N {{{f}_{m}}\sin \pi m\tau } ;\quad \tau \in \left[ {0,1} \right],$
где $\tau = t{\text{/}}{{t}_{0}}$, ${{t}_{0}}$ ‒ время движения космического корабля.

Граничную функцию ${{M}_{{2k}}}$ зададим рекуррентным соотношением:

(4.2)
${{M}_{{2k}}} = {{M}_{{2(k - 1)}}} + {{x}^{{(2k)}}}(0){{P}_{{2k}}} + {{x}^{{(2k)}}}({{t}_{0}}){{Q}_{{2k}}}$

Здесь ${{P}_{{2k}}}\left( x \right)$, ${{Q}_{{2k}}}\left( x \right)$ – быстрые полиномы, которые в случае граничных условий Дирихле вычисляются по рекуррентным интегральным формулам

${{P}_{0}} = 1 - \tau ,\quad {{Q}_{0}} = \tau $
(4.3)
$\begin{gathered} {{P}_{{2k}}}\left( t \right) = \int\limits_0^t {\int\limits_0^{{{t}_{{\left( 1 \right)}}}} {{{P}_{{2k - 2}}}\left( {{{t}_{{\left( 2 \right)}}}} \right)d{{t}_{{\left( 2 \right)}}}} d{{t}_{{\left( 1 \right)}}}} - \frac{t}{{{{t}_{0}}}}\int\limits_0^{{{t}_{0}}} {\int\limits_0^{{{t}_{{\left( 1 \right)}}}} {{{P}_{{2k - 2}}}\left( {{{t}_{{\left( 2 \right)}}}} \right)d{{t}_{{\left( 2 \right)}}}} d{{t}_{{\left( 1 \right)}}}} \\ {{Q}_{{2k}}}\left( t \right) = \int\limits_0^t {\int\limits_0^{{{t}_{{\left( 1 \right)}}}} {{{Q}_{{2k - 2}}}\left( {{{t}_{{\left( 2 \right)}}}} \right)d{{t}_{{\left( 2 \right)}}}} d{{t}_{{\left( 1 \right)}}}} - \frac{t}{{{{t}_{0}}}}\int\limits_0^{{{t}_{0}}} {\int\limits_0^{{{t}_{{\left( 1 \right)}}}} {{{Q}_{{2k - 2}}}\left( {{{t}_{{\left( 2 \right)}}}} \right)d{{t}_{{\left( 2 \right)}}}} d{{t}_{{\left( 1 \right)}}}} \\ \end{gathered} $
$k = 1 \div p,\quad 0 \leqslant {{t}_{{\left( 2 \right)}}} \leqslant {{t}_{{\left( 1 \right)}}},\quad 0 \leqslant {{t}_{{\left( 1 \right)}}} \leqslant t,\quad 0 \leqslant t \leqslant {{t}_{0}}$

Полиномы ${{P}_{0}}\left( t \right)$, ${{Q}_{0}}\left( t \right)$ являются производящими, через них при помощи двукратных интегралов находятся все остальные с четными индексами. Быстрые полиномы ${{P}_{{2k}}}\left( t \right)$, ${{Q}_{{2k}}}\left( t \right)$ можно вычислить также при помощи производных [6]:

$\begin{gathered} P_{{2k}}^{{''}}\left( t \right) = {{P}_{{2k - 2}}}\left( t \right),\quad {{P}_{{2k}}}\left( 0 \right) = {{P}_{{2k}}}\left( {{{t}_{0}}} \right) = 0;\quad k \ne 0,\quad k = 1 \div p \\ Q_{{2k}}^{{''}}\left( t \right) = {{Q}_{{2k - 2}}}\left( t \right),\quad {{Q}_{{2k}}}\left( 0 \right) = {{Q}_{{2k}}}\left( {{{t}_{0}}} \right) = 0 \\ \end{gathered} $

Значение индекса $2k$ в (4.1) может быть выбрано произвольно, но не ниже порядка старшей производной в дифференциальном уравнении, в которое предполагается подставлять разложение (4.1). Быстрое разложение (4.1) допускает почленное дифференцирование ряда Фурье $2k$ раз, при этом ряды остаются быстро сходящимися. За счет специальной конструкции полиномов ${{P}_{{2k}}}\left( t \right)$, ${{Q}_{{2k}}}\left( t \right)$ в (4.3) граничная функция ${{M}_{{2k}}}\left( t \right)$ значительно увеличивает скорость сходимости ряда Фурье по сравнению с классическим рядом, с увеличением порядка $2k$ скорость сходимости ряда существенно возрастает [13, 14].

Для определенности и простоты дальнейших выкладок выберем граничную функцию в быстром синус-разложении (4.1) невысокого порядка при $k = 1$, т.е. ${{M}_{2}}$. Тогда выражение для ${\mathbf{r}}\left( t \right)$ будет иметь вид

(4.4)
$\begin{gathered} {\mathbf{r}}\left( t \right) = {\mathbf{r}}\left( 0 \right)\left( {1 - \tau } \right) + {\mathbf{r}}\left( {{{t}_{0}}} \right)\tau + {\mathbf{\ddot {r}}}\left( 0 \right)t_{0}^{2}\left( {\frac{1}{2}{{\tau }^{2}} - \frac{1}{6}{{\tau }^{3}} - \frac{1}{3}\tau } \right) + \\ + \;\frac{1}{6}{\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right)t_{0}^{2}\left( {{{\tau }^{3}} - \tau } \right) + \sum\limits_{m = 1}^N {{{{\mathbf{r}}}_{m}}\sin \pi m\tau } , \\ \end{gathered} $
где ${{{\mathbf{r}}}_{m}} = \left( {{{x}_{m}},{{y}_{m}},{{z}_{m}}} \right)$.

Теорема о скорости убывания коэффициентов быстрого синус-разложения.

Пусть $x\left( t \right) \in L_{2}^{4}\left( {t \in \left[ {0,{{t}_{0}}} \right]} \right)$, где $L_{2}^{4}$ – классы функций Соболева–Лиувилля [15]. Тогда коэффициенты ${{x}_{m}}$ быстрого разложения (4.4) для $x\left( t \right)$ с ростом номера $m$ будут убывать по закону ${{x}_{m}}\sim {{\left( {m\pi } \right)}^{{ - 4}}}$.

Покажем это на примере для $x\left( t \right)$. Из (4.4) запишем ряд Фурье для $x\left( t \right)$ разностью между $x\left( t \right)$ и ее граничной функции:

$\begin{gathered} \sum\limits_{m = 1}^N {{{x}_{m}}\sin \pi m\tau } = x\left( t \right) - x\left( 0 \right)\left( {1 - \tau } \right) - x\left( {{{t}_{0}}} \right)\tau - \\ - \;t_{0}^{2}\left[ {\ddot {x}\left( 0 \right)\left( {\frac{1}{2}{{\tau }^{2}} - \frac{1}{6}{{\tau }^{3}} - \frac{1}{3}\tau } \right) + \frac{1}{6}\ddot {x}\left( {{{t}_{0}}} \right)\left( {{{\tau }^{3}} - \tau } \right)} \right] \\ \end{gathered} $

Отсюда для коэффициентов ${{x}_{m}}$ имеем интегральную формулу

(4.5)
$\begin{gathered} {{x}_{m}} = \frac{2}{{{{t}_{0}}}}\int\limits_0^{{{t}_{0}}} {\left( {x{{{\left( t \right)}}^{{^{{}}}}}{\kern 1pt} - x\left( 0 \right)\left( {1 - \tau } \right) - x\left( {{{t}_{0}}} \right)\tau - } \right.} \\ \left. { - \;t_{0}^{2}\left[ {\ddot {x}\left( 0 \right)\left( {\frac{1}{2}{{\tau }^{2}} - \frac{1}{6}{{\tau }^{3}} - \frac{1}{3}\tau } \right) + \frac{1}{6}\ddot {x}\left( {{{t}_{0}}} \right)\left( {{{\tau }^{3}} - \tau } \right)} \right]} \right)\sin \pi m\tau dt \\ \end{gathered} $

К интегралу в (4.5) применим интегрирование по частям четыре раза, получим

${{x}_{m}} = \frac{{2t_{0}^{3}}}{{{{{\left( {\pi m} \right)}}^{4}}}}\int\limits_0^{{{t}_{0}}} {{{x}^{{\left( 4 \right)}}}\left( t \right)\sin m\pi \tau dt} $

Что и требовалось доказать.

Аналогичные оценки получим для коэффициентов Фурье ${{y}_{m}}$ и ${{z}_{m}}$.

Отсюда следует, что решение системы в виде (4.4) гарантирует возможность его почленного двукратного дифференцирования по времени, которые присутствуют в уравнениях движения (2.1). При использовании граничных функций ${{M}_{{2k}}}\left( t \right)$ более высокого порядка $2k \geqslant 4$ степень скорости убывания коэффициентов Фурье возрастает пропорционально порядку граничной функции и равна $2k + 2$. В данной работе ограничимся простейшим вариантом, когда $k = 1$.

Выполняя граничные условия (2.4) и дифференциальные уравнения (2.1), для построения решения в виде (4.4) надо найти следующие постоянные коэффициенты

(4.6)
${\mathbf{r}}\left( 0 \right),{\mathbf{r}}\left( {{{t}_{0}}} \right),{\mathbf{\ddot {r}}}\left( 0 \right),{\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right),{{{\mathbf{r}}}_{n}};\quad n = 1 \div N$

Неизвестные ${\mathbf{r}}\left( 0 \right),{\mathbf{r}}\left( {{{t}_{0}}} \right)$ получаем из граничных условий (2.4). Для нахождения остальных $3N + 6$ коэффициентов (4.6) подставим выражения для координат траектории полета из (4.4) в уравнения движения (2.1):

(4.7)
$\begin{gathered} {\mathbf{\ddot {r}}}\left( 0 \right)\left( {1 - \tau } \right) + {\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right)\tau - \sum\limits_{m = 1}^N {{{{\mathbf{r}}}_{m}}{{{\left( {\frac{{\pi m}}{{{{t}_{0}}}}} \right)}}^{2}}\sin \pi m\tau } + \alpha \frac{{\mathbf{r}}}{{{{r}^{3}}}} + \beta \left\{ {\frac{{{\mathbf{r}}\left( {{{t}_{0}}} \right) - {\mathbf{r}}\left( 0 \right)}}{{{{t}_{0}}}} + } \right. \\ \left. { + \;{{t}_{0}}\left[ {{\mathbf{\ddot {r}}}\left( 0 \right)\left( {\tau - \frac{1}{2}{{\tau }^{2}} - \frac{1}{3}} \right) + \frac{1}{2}{\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right)\left( {{{\tau }^{2}} - \frac{1}{3}} \right)} \right] + \sum\limits_{m = 1}^N {{{{\mathbf{r}}}_{m}}\frac{{\pi m}}{{{{t}_{0}}}}\cos \pi m\tau } } \right\} = {\mathbf{P}}, \\ \end{gathered} $
где ${\mathbf{P}} = \left( {{{P}_{x}},{{P}_{y}},{{P}_{z}}} \right)$.

В дальнейшем вследствие нелинейности системы (4.7) с дробными степенями относительно неизвестных (4.6), удобно использовать поточечный метод, который в литературе также называют методом коллокаций [16], либо методом тригонометрической интерполяции [17]. При рассмотрении нелинейных краевых задач наиболее просто применять тригонометрическую синус-интерполяцию, так как ряд Фурье быстро сходится и потому в ряде достаточно удерживать небольшое количество неизвестных слагаемых.

Вначале в соответствии с методом быстрых синус-разложений [6] в системе (4.7) следует положить $t = 0$

(4.8)
$\begin{gathered} {\mathbf{\ddot {r}}}\left( 0 \right) + \alpha {{{\mathbf{H}}}_{*}} + \beta \left( {\frac{{{\mathbf{r}}\left( {{{t}_{0}}} \right) - {\mathbf{r}}\left( 0 \right)}}{{{{t}_{0}}}} - \frac{1}{3}{{t}_{0}}\left( {{\mathbf{\ddot {r}}}\left( 0 \right) + \frac{1}{2}{\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right)} \right) + \sum\limits_{m = 1}^N {{{{\mathbf{r}}}_{m}}\frac{{\pi m}}{{{{t}_{0}}}}} } \right) = {{\left. {\mathbf{P}} \right|}_{{t = 0}}} \\ {{{\mathbf{H}}}_{*}} = \left( {\frac{a}{{{{{\left( {{{a}^{2}} + R_{{\text{З}}}^{2}} \right)}}^{{{3 \mathord{\left/ {\vphantom {3 2}} \right. \kern-0em} 2}}}}}},0,\frac{{{{R}_{{\text{З}}}}}}{{{{{\left( {{{a}^{2}} + R_{{\text{З}}}^{2}} \right)}}^{{{3 \mathord{\left/ {\vphantom {3 2}} \right. \kern-0em} 2}}}}}}} \right) \\ \end{gathered} $
и затем в (4.7) взять $t = {{t}_{0}}$:

(4.9)
${\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right) + \alpha \frac{{{\mathbf{r}}\left( {{{t}_{0}}} \right)}}{{{{r}^{3}}\left( {{{t}_{0}}} \right)}} + \beta \left( {\frac{{{\mathbf{r}}\left( {{{t}_{0}}} \right) - r\left( 0 \right)}}{{{{t}_{0}}}} + \frac{{{{t}_{0}}}}{3}\left[ {\frac{1}{2}{\mathbf{\ddot {r}}}\left( 0 \right) - {\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right)} \right] + \sum\limits_{n = 1}^N {{{{\left( { - 1} \right)}}^{n}}{{{\mathbf{r}}}_{n}}\frac{{\pi n}}{{{{t}_{0}}}}} } \right) = {{\left. {\mathbf{P}} \right|}_{{t = {{t}_{0}}}}}$

В (4.8) и (4.9) имеем 6 уравнений. Для нахождения оставшихся $3N$ неизвестных ${{x}_{j}},{{y}_{j}},{{z}_{j}}$, $j = 1 \div N$ отрезок $\left[ {0,{{t}_{0}}} \right]$ разделим на $N + 1$ равных частей интерполяционными точками ${{t}_{j}} = {{{{t}_{0}}} \mathord{\left/ {\vphantom {{{{t}_{0}}} {\left( {N + 1} \right)}}} \right. \kern-0em} {\left( {N + 1} \right)}}$ и запишем уравнения (4.7) в каждой расчетной точке при $t = {{t}_{j}}$, $j = 1 \div N$:

(4.10)
$\begin{gathered} {\mathbf{\ddot {r}}}\left( 0 \right)\left( {1 - {{\tau }_{j}}} \right) + {\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right){{\tau }_{j}} - \sum\limits_{n = 1}^N {{{{\mathbf{r}}}_{n}}{{{\left( {\frac{{\pi n}}{{{{t}_{0}}}}} \right)}}^{2}}\sin \pi n{{\tau }_{j}}} + {{\left. {\frac{{\alpha {\mathbf{r}}}}{{{{r}^{3}}}}} \right|}_{{t = {{t}_{j}}}}} + \beta {{t}_{0}}\left[ {\frac{{{\mathbf{r}}\left( {{{t}_{0}}} \right) - {\mathbf{r}}\left( 0 \right)}}{{t_{0}^{2}}}\mathop + \limits_{_{{_{{}}}}}^{^{{}}} } \right. \\ + \;{\mathbf{\ddot {r}}}\left( 0 \right)\left( {{{\tau }_{j}} - \frac{1}{2}\tau _{j}^{2} - \frac{1}{3}} \right) + \frac{1}{2}{\mathbf{\ddot {r}}}\left( {{{t}_{0}}} \right)\left( {\tau _{j}^{2} - \frac{1}{3}} \right) + \left. {\sum\limits_{n = 1}^N {{{{\mathbf{r}}}_{n}}\frac{{\pi n}}{{t_{0}^{2}}}\cos \pi n{{\tau }_{j}}} } \right] = {{\left. {\mathbf{P}} \right|}_{{t = {{t}_{j}}}}} \\ {{t}_{j}} = {{j{{t}_{0}}} \mathord{\left/ {\vphantom {{j{{t}_{0}}} {\left( {N + 1} \right)}}} \right. \kern-0em} {\left( {N + 1} \right)}},\quad {{\tau }_{j}} = {{t}_{j}}{\text{/}}{{t}_{0}},\quad j = 1 \div N \\ \end{gathered} $

Полученную систему (4.8)–(4.10) решаем в Maple. Подставляя найденные коэффициенты (4.6) в (4.4), имеем приближенное аналитическое решение. Отметим, что в полученном решении коэффициенты быстрого разложения находятся численно, но само решение задачи записано в аналитическом виде. Точность вычислений можно повысить увеличением порядка граничной функции и/или увеличением количества членов ряда Фурье. Следует отметить, что возрастание порядка граничной функции приводит к увеличению количества неизвестных постоянных коэффициентов, для определения которых необходимо записать дополнительные алгебраические уравнения. Например, если в быстром разложении (4.1) использовать граничную функцию четвертого порядка ${{M}_{4}}$ $\left( {k = 2} \right)$, то перечень (4.6) увеличится на шесть новых неизвестных

(4.11)
$\ddddot x\left( 0 \right),\;\;\ddddot x\left( {{{t}_{0}}} \right),\;\;\ddddot y\left( 0 \right),\;\;\ddddot y\left( {{{t}_{0}}} \right),\;\;\ddddot z\left( 0 \right),\;\;\ddddot z\left( {{{t}_{0}}} \right)$
(по сравнению с вышеописанным решением для случая граничной функции второго порядка ${{M}_{2}}$). Алгебраические уравнения для нахождения неизвестных (4.11) получаются следующими действиями:

1) дифференцируем уравнения (4.7) два раза,

2) в полученные уравнения последовательно подставим $t = 0$ и $t = {{t}_{0}}$.

Записанные шесть уравнений добавляем к системе (4.8)–(4.10).

В общем случае, при использовании граничной функции ${{M}_{{2k}}}$, для определения ее неизвестных коэффициентов будем иметь $2k - 2$ алгебраических уравнений, которые получаются при дифференцировании уравнений (4.7) $2k - 2$ раз и последовательном подставлением в уравнения, полученные от четных производных значений $t = 0$ и $t = {{t}_{0}}$.

Движение космического корабля без потери массы возможно лишь при выключенном двигателе. Обобщим постановку задачи. Предположим, что топливо расходуется по линейному закону, поэтому масса корабля выражается функцией $m(t) = {{m}_{0}}(1 - \lambda t)$, где ${{m}_{0}}$ – стартовая масса, $\lambda $ – коэффициент пропорциональности. В этом случае уравнения движения примут вид:

(4.12)
${{m}_{0}}(1 - \lambda t){\mathbf{\ddot {r}}} + \beta {\mathbf{\dot {r}}} + \alpha {{m}_{0}}(1 - \lambda t)\frac{{\mathbf{r}}}{{{{r}^{3}}}} = {\mathbf{P}},$
где в силу ${\mathbf{P}} = \left( {{{P}_{x}},{{P}_{y}},{{P}_{z}}} \right)$ могут входить реактивная сила $\dot {m}\,{{{\mathbf{v}}}_{{\mathbf{r}}}}$, сила давления солнечного ветра, сила сопротивления космической пыли и т.д.

Разделим каждое уравнение на ${{m}_{0}}$ и обозначим ${{\beta }_{0}} = {\beta \mathord{\left/ {\vphantom {\beta {{{m}_{0}}}}} \right. \kern-0em} {{{m}_{0}}}}$ и ${{{\mathbf{P}}}_{0}} = {{\mathbf{P}} \mathord{\left/ {\vphantom {{\mathbf{P}} {{{m}_{0}}}}} \right. \kern-0em} {{{m}_{0}}}}$, окончательно получим:

(4.13)
$(1 - \lambda t){\mathbf{\ddot {r}}} + {{\beta }_{0}}{\mathbf{\dot {r}}} + \alpha (1 - \lambda t)\frac{{\mathbf{r}}}{{{{r}^{3}}}} = {{{\mathbf{P}}}_{0}}$

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

${{{\mathbf{P}}}_{0}} = \left( {1 - \lambda t} \right){{{\mathbf{\ddot {r}}}}_{*}} + {{\beta }_{0}}{{{\mathbf{\dot {r}}}}_{*}} + \alpha (1 - \lambda t)\frac{{{{{\mathbf{r}}}_{*}}}}{{r_{*}^{3}}};\quad {{{\mathbf{r}}}_{*}} = \left( {a\cos \omega t,a\sin \omega t,{{R}_{З}} + wt} \right)$

Ход решения задачи (4.13), (2.4) аналогичен описанному выше для задачи (2.1), (2.4).

5. Полученные результаты и их анализ. Для вычислительных экспериментов определим значения параметров, входящих в систему (2.1): ${{t}_{0}} = 30$ с, $\omega = {\pi \mathord{\left/ {\vphantom {\pi {3600}}} \right. \kern-0em} {3600}}$ с–1, $w = 2000$ м/с, $a = 100$ м, $g = 9.8$ м/с2, ${{R}_{{\text{З}}}} = 6372$ км, $\beta = 0.1$. Все расчеты выполнены в системе компьютерной алгебры Maple.

Абсолютные погрешности траектории корабля, его скорости и ускорения вычислим по формулам

(5.1)
$\begin{gathered} \Delta s = \sqrt {{{{\left( {x{{{\kern 1pt} }_{*}} - \;x} \right)}}^{2}} + {{{\left( {y{{{\kern 1pt} }_{*}} - \;y} \right)}}^{2}} + {{{\left( {z{{{\kern 1pt} }_{*}} - \;z} \right)}}^{2}}} \\ \Delta {v} = \sqrt {{{{\left( {\dot {x}{{{\kern 1pt} }_{*}} - \;\dot {x}} \right)}}^{2}} + {{{\left( {\dot {y}{{{\kern 1pt} }_{*}} - \;\dot {y}} \right)}}^{2}} + {{{\left( {\dot {z}{{{\kern 1pt} }_{*}} - \;\dot {z}} \right)}}^{2}}} \\ \Delta a = \sqrt {{{{\left( {{{{\ddot {x}}}_{*}} - \;\ddot {x}} \right)}}^{2}} + {{{\left( {\ddot {y}{{{\kern 1pt} }_{*}} - \;\ddot {y}} \right)}}^{2}} + {{{\left( {{{{\ddot {z}}}_{*}} - \;\ddot {z}} \right)}}^{2}}} \\ \end{gathered} $

В табл. 1 приведены результаты расчетов абсолютных погрешностей (5.1) приближенного аналитического решения задачи (2.1), (2.4) методом быстрой тригонометрической синус-интерполяции при различных граничных функциях и количестве членов ряда Фурье. Мы видим, что при увеличении количества членов ряда Фурье в 2 раза, погрешность уменьшается как минимум на один порядок, а при увеличении порядка граничной функции на 2, точность увеличивается сразу на 6–8 порядков. Абсолютная погрешность скорости ракеты убывает медленнее, а ускорения еще медленнее по сравнению с погрешностью траектории при увеличении количества членов ряда быстрого синус-разложения. Однако с увеличением порядка граничной функции все три погрешности ведут себя одинаково.

Таблица 1.

Абсолютная погрешность вычислений (5.1)

Исследуемый объект N Граничная функция
${{M}_{2}}$ ${{M}_{4}}$ ${{M}_{6}}$ ${{M}_{8}}$
$\Delta s$ 5 2.3 × 10–10 9.2 × 10–17 3.4 × 10–23 1.1 × 10–29
10 2.2 × 10–11 5.3 × 10–18 1.3 × 10–24 2.7 × 10–31
20 1.6 × 10–12 1.6 × 10–19 1.7 × 10–26 1.6 × 10–33
40 1.1 × 10–13 3.4 × 10–21 1.4 × 10–28 5.6 × 10–36
80 6.8 × 10–15 6.6 × 10–23 8.1 × 10–31 1.2 × 10–38
$\Delta {v}$ 5 1.4 × 10–10 7.2 × 10–17 3.3 × 10–23 1.2 × 10–29
10 2.3 × 10–11 6.5 × 10–18 1.8 × 10–24 4.3 × 10–31
20 3.2 × 10–12 3.5 × 10–19 3.9 × 10–26 4.3 × 10–33
40 4.2 × 10–13 1.4 × 10–20 5.8 × 10–28 2.5 × 10–35
80 5.5 × 10–14 5.3 × 10–22 6.7 × 10–30 8.9 × 10–38
$\Delta a$ 5 5.8 × 10–11 3.9 × 10–17 2.2 × 10–23 9.6 × 10–30
10 1.7 × 10–11 5.6 × 10–18 1.8 × 10–24 4.8 × 10–31
20 4.6 × 10–12 5.2 × 10–19 6.6 × 10–26 7.7 × 10–33
40 1.2 × 10–12 4.1 × 10–20 1.8 × 10–27 7.7 × 10–35
80 3.1 × 10–13 2.9 × 10–21 3.7 × 10–29 5.2 × 10–37

Таблица 2 содержит сведения о таких параметрах вычислительного процесса как время расчета и количество значащих цифр. Видно, что время расчета резко возрастает при увеличении количества членов ряда Фурье, а при увеличении порядка граничной функции увеличение времени расчета незначительно.

Таблица 2.

Параметры вычислительного процесса

Параметр N Граничная функция
${{M}_{2}}$ ${{M}_{4}}$ ${{M}_{6}}$ ${{M}_{8}}$
Время расчета, с 5 2 3 4 5
10 3 4 6 8
20 9 11 14 20
40 50 55 72 81
80 330 408 482 1512
Минимальное количество значащих цифр в расчете 5 15 22 29 35
10 20 24 30 37
20 25 25 35 40
40 40 35 40 45
80 60 70 75 80

Проанализировав данные табл. 1 и 2, можно сделать вывод, что для достижения высокой точности решения необходимо увеличить порядок граничной функции, а количество членов ряда Фурье оставить минимальным.

Для численного расчета траектории полета космического корабля с включенным двигателем зададим $\lambda = {1 \mathord{\left/ {\vphantom {1 {\left( {2{{t}_{0}}} \right)}}} \right. \kern-0em} {\left( {2{{t}_{0}}} \right)}}$, ${{\beta }_{0}} = 0.1$. Остальные параметры, входящие в задачу (4.13), (2.4) возьмем такими же как при решении задачи (2.1), (2.4). В табл. 3 приведены абсолютные погрешности (5.1) расчета методом быстрой тригонометрической синус-интерполяции с граничной функцией ${{M}_{4}}$.

Таблица 3.

Погрешности активного участка полета

N Абсолютные погрешности Время расчета, с
$\Delta s$ $\Delta {v}$ $\Delta a$
5 9.1 × 10–17 7.1 × 10–17 3.2 × 10–17 3
10 5.3 × 10–18 6.2 × 10–18 5.5 × 10–18 6
20 1.5 × 10–19 3.6 × 10–19 5.2 × 10–19 14
40 3.5 × 10–21 1.4 × 10–20 4.1 × 10–20 60
80 6.1 × 10–23 5.1 × 10–22 2.9 × 10–21 402

Анализ табличных данных показывает ту же динамику изменения погрешности, что и в случае с постоянной массой (см. табл. 1), причем с тем же порядком точности для граничной функции ${{M}_{4}}$. При увеличении количества членов ряда Фурье в 2 раза, погрешность уменьшается как минимум на один порядок. Погрешности скорости и ускорения имеют тот же порядок, что и погрешность траектории для 5, 10 и 20 членов ряда Фурье, и отличаются на порядок для 80 членов.

Таким образом, метод быстрой тригонометрической синус-интерполяции показывает отличные результаты для расчета траектории полета космического корабля как с включенным двигателем, так и с выключенным.

Заключение. С помощью метода быстрой тригонометрической синус-интерполяции найдено приближенное аналитическое решение задачи о движении космического корабля из начальной точки в конечный пункт за определенное время. Исследование полученного решения выявило, что для достижения высокой точности решения в методе быстрой тригонометрической синус-интерполяции необходимо увеличить порядок граничной функции (в четыре раза), а количество членов ряда Фурье оставить минимальным (пять членов). Преимуществом метода быстрой тригонометрической синус-интерполяции является и то, что решение получено в аналитическом виде, что позволяет проводить различные исследования свойств траектории. В системе (2.1) можно было бы учесть силы притяжения к Солнцу и к Луне. Усложнения были бы не принципиальные, временные затраты на ЭВМ при этом существенно не изменятся.

Метод быстрых разложений является очень эффективным при решении краевых задач не только для дифференциальных уравнений в частных производных [18, 19], но и для систем обыкновенных дифференциальных уравнений.

Список литературы

  1. Карагодин В.В. Приближенные методы расчета внеатмосферного активного участка траектории // Тр. МАИ. 2013. Вып. 66. http://trudymai.ru/published.php?ID=40267

  2. Аппазов Р.Ф., Сытин О.Г. Методы проектирования траекторий носителей и спутников Земли. М.: Наука, 1987. 440 с.

  3. Беневольский С.В. Математические модели движения для синтеза методов наведения перспективных баллистических ракет // Оборон. техн. 2007. № 3–4. С. 12–16.

  4. Хайрер Э., Нёрсетт С., Ваннер Г. Решение обыкновенных дифференциальных уравнений. Нежесткие задачи. М.: Мир, 1990. 512 с.

  5. Беневольский С.В., Козлов П.Г. Полуаналитический метод восстановления траекторий ЛА по обобщенным проектным параметрам и параметрам системы управления и перспективы его использования // Электр. науч.-технич. Изд. “Наука и образование”. 2011. № 10. http://technomag.edu.ru/doc/216895.html

  6. Чернышов А.Д. Метод быстрых разложений для решения нелинейных дифференциальных уравнений // ЖВММФ. 2014. Т. 54. № 1. С. 13–24.

  7. Чернышов А.Д., Горяйнов В.В. Решение одного нелинейного интегро-дифференциального уравнения методом быстрых разложений // Вестн. ЧГПУ им. И.Я. Яковлева. Сер.: Механ. пред. сост. 2012. № 4(12). С. 105–112.

  8. Чернышов А.Д. Решение нелинейного уравнения теплопроводности для криволинейной области с условиями Дирихле методом быстрых разложений // ИФЖ. 2018. Т. 91. № 2. С. 456−468.

  9. Чернышов А.Д. Решение двухфазной задачи Стефана с внутренним источником и задач теплопроводности методом быстрых разложений // ИФЖ. 2021. Т. 94. № 1. С. 101–120.

  10. Чернышов А.Д., Горяйнов В.В., Чернышов О.А. Применение метода быстрых разложений для расчета траекторий космических кораблей // Изв. вузов. Авиац. техн. 2015. № 2. С. 41–47.

  11. Chernyshov A.D., Saiko D.S., Kovaleva E.N. Universal fast expansion for solving nonlinear problems // J. Physics: Conf. Ser. 2020. V. 1479. Art. no. 012147.

  12. Горячева И.Г., Горячев А.П. Контактные задачи о скольжении штампа с периодическим рельефом по вязкоупругой полуплоскости // ПММ. 2016. Т. 80. № 1. С. 103–116.

  13. Чернышов А.Д., Горяйнов В.В., Лешонков О.В., Соболева Е.А., Никифорова О.Ю. Сравнение скорости сходимости быстрых разложений с разложениями в классический ряд Фурье // Вестн. ВГУ. Сер.: Сист. анализ и информ. технол. 2019. № 1. С. 27–34.

  14. Чернышов А.Д., Горяйнов В.В. О выборе оптимального порядка граничной функции в быстром разложении // Вестн. ВГУ. Сер.: Сист. анализ и информ. технол. 2011. № 1. С. 60–65.

  15. Ильин В.А. Спектральная теория дифференциальных операторов. М. Наука, 1991. 368 с.

  16. Исаев В.И., Шапеев В.П., Идимешев С.В. Варианты метода коллокаций и наименьших квадратов повышенной точности для численного решения уравнения Пуассона // Вычисл. технол. 2011. Т.16. № 1. С. 85–93.

  17. Бахвалов Н.С., Жидков Н.П., Кобельков Г.М. Численные методы. М. Наука, 1987. 800 с.

  18. Горяйнов В.В., Попов М.И., Чернышов А.Д. Решение задачи о напряжениях в остром клиновидном режущем инструменте методом быстрых разложений и проблема согласования граничных условий // Изв. РАН. МТТ. 2019. № 5. С. 113–130.

  19. Чернышов А.Д., Попов В.М., Горяйнов В.В., Лешонков О.В. Исследование контактного термического сопротивления в конечном цилиндре с внутренним источником методом быстрых разложений и проблема согласования граничных условий // ИФЖ. 2017. Т. 90. № 5. С. 1288–1297.

Дополнительные материалы отсутствуют.