Синтез итерационных алгоритмов решения краевых задач и нелинейных уравнений

УДК 629.191

В. И. Миронов, Ю. В. Миронов, Р. М. Юсупов

СИНТЕЗ ИТЕРАЦИОННЫХ АЛГОРИТМОВ РЕШЕНИЯ КРАЕВЫХ ЗАДАЧ И НЕЛИНЕЙНЫХ УРАВНЕНИЙ

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

Ключевые слова: итерационный алгоритм, краевая задача, нелинейные уравнения.

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

К настоящему времени разработано большое число различных универсальных методов решения таких задач, они приведены в литературе, в частности, в работах [1—6]. Вместе с тем в различных областях знаний получены приближенные решения многих упрощенных модельных задач, имеющие ограниченное применение. Так, например, применительно к динамике полета космических аппаратов получено множество аналитических и упрощенных численных алгоритмов решения задач маневрирования в модельных гравитационных полях: однородном, линеаризованном, однородном центральном, квазиньютоновском и ньютоновском. Обычно эти алгоритмы используются для определения начального приближения при решении более сложных задач. Однако они также могут быть использованы в качестве базовых элементов при создании алгоритмов численного решения усложненных задач. Некоторые вопросы синтеза алгоритмов такого рода рассматривались в работах [7].

В настоящей работе рассматривается метод решения краевых задач и нелинейных уравнений — метод приближенного корректирующего оператора (ПКО), который позволяет использовать возможные упрощенные алгоритмы приближенного расчета в схеме численного поиска точного решения полной задачи. Такой подход расширяет конструктивный базис синтеза быстродействующих алгоритмов решения указанных краевых задач и нелинейных уравнений.

Метод приближенного корректирующего оператора. Пусть поведение объекта описывается векторным дифференциальным уравнением

х = ф(хq,¿X х(^ = ^ г е^Т], где х — и-мерный вектор фазового состояния; q — т-мерный вектор управляющих параметров.

Требуется найти вектор q, который переводит управляемый объект из исходного состояния х(г0 ) = х0 на требуемую траекторию движения, обеспечивающую в заданный момент времени Т достижение заданных граничных условий рт = р[х(Т)]. Будем считать, что вектор Рт имеет размерность т и т < п. Предполагается, что функции ф(х, q, г) и р[х(Т)] определены в некоторых заданных областях изменения аргументов, непрерывны и обладают необходимой степенью гладкости, так что обеспечивается единственность решения рассматриваемой задачи.

Предположим, что известен точный нелинейный оператор связи А между вектором начального состояния объекта х0, вектором управляющих параметров q и состоянием объекта в конечный момент времени р[х(Г)]. Оператор А задается интегрированием соответствующей системы дифференциальных уравнений. Тогда требуемое значение вектора q, обеспечивающее достижение заданной терминальной точки рт, будет удовлетворять уравнению

рг = А

Допустим далее, что известен оператор Д, устанавливающий приближенную связь между векторами хо, Рт и q

Рт = А1

Будем также считать, что для приближенного оператора А известен обратный оператор

А-1 по вектору параметров управления q. В практической ситуации оператор А может формироваться путем приближенного учета физических факторов, определяющих движение, либо формального упрощения точного оператора интегрирования А, либо комбинацией этих приемов.

Используя приближенный оператор А1, запишем уравнение (1) в следующей эквивалентной форме:

А = рт - ДА^); М^) = А^) - А

Подвергнем левую и правую части равенства (3) операторному преобразованию А-1:

q = А1-1[рт -ДАя)].

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

q

Значение вектора q в первом приближении находим из уравнения (4) при ДА^(0)) = 0 , так что q(1) = А;-1 [рт ].

Согласно (4), на каждом шаге решения задачи необходимо вычислять значение разностного оператора ДА^( к)). Получим более эффективный алгоритм, в котором исключается необходимость определения этого оператора.

Рассмотрим несколько итераций. Согласно общей схеме расчета, величина q в первом

приближении находится из условия рт = А1^(1)) , так что

q

На втором шаге значение q(2) находится из условия

рт = А^(2)) + А^(1)) - А^(1)). Учитывая (5), получаем уравнение

2рт = А^

отсюда

q(2) = А1-1[2рт - A(q(1))] = q(2) = А-1^ - А(1)]; А(1) = A(q(1)) - рт.

На третьем шаге величина q(3) находится из условия

рт = А^

однако из выражения (6) следует, что

Дя(2)) = 2рт - А(я(1)). Тогда после подстановки (8) в (7) находим:

рт = А1(я(3)) + А(1) + А(2); А(2) = Дя(2)) - рг,

откуда следует, что

Я(3) = А1-1[рт - А(1) - А(2)]. Продолжив этот процесс, приходим к следующей вычислительной схеме:

Я( к+1) = АГ1

„(О

; к = 1,2,...

5 • • • 1

Рт-Е А(

_ г=1

А(г) = Л(^г)) - рт = р(я(г), хо, Т) - рт; р(я(г), хо, Т) - Л(^г)). Выражения (10) удобно представить в следующем виде:

,(к +1) _

= М

рт -Е А'

г=1

(о

к = 1,2,3,...

1

где М[•] - Л1 1[^] — приближенный корректирующий оператор, определяющий алгоритм приближенного решения краевой задачи (1).

Условия сходимости данного метода устанавливаются на основе принципа сжимающих отображений [1, 2]. Можно показать, что вычислительный процесс сходится со скоростью

геометрической прогрессии. Очевидно, что чем ближе приближенные операторы Л1 и Л]"1 к

точным операторам Л и Л-1 , тем выше скорость сходимости.

В целом, метод ПКО отличается достаточно высокой экономичностью с вычислительной точки зрения, поскольку значение неизвестного вектора q полностью уточняется на каждом итерационном шаге путем однократного интегрирования дифференциальных уравнений движения объекта.

Рассмотренный выше метод ПКО естественным образом распространяется на задачи решения нелинейных уравнений. Согласно этому методу, для решения уравнения

X = 0

применяется соотношение

qk+1 = М

-Е х , т)

1=0

к = 1,2,3,

где М [•] — ПКО, определяющий алгоритм решения приближенного уравнения X= 0 .

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

Важной особенностью метода ПКО является то обстоятельство, что на каждой итерации значения рт и X вычисляются один раз. Этим обеспечивается высокая экономичность

вычислений. Ниже это будет показано на примере. Значение оператора М [•] может уточняться в ходе вычислительного процесса по результатам каждой итерации или через несколько итераций. В этом случае можно говорить о комбинированных вариантах использования метода ПКО.

С общих методологических позиций очевидно, что в качестве оператора М [•] могут

рассматриваться и операторы различных известных методов численного решения нелинейных уравнений, таких как метод Ньютона и др. При линейном представлении X вида

X = Aq + Ь

из (13) следует

qк+1 = qk - ); м = А"1,

что по форме напоминает модифицированный (упрощенный) метод Ньютона и в частном

"дхП-1

совпадает с ним.

случае при M =

dq

Jq=qo

Оператор M может уточняться в ходе итераций по некоторому правилу

Mk = Mk (q¿ > qjt-b-> qk-s )>

тогда соотношение (13) при s = 1 принимает вид

q k+i = q k- Mk

-1

В частном случае при Mk =

дХ

dq

из (14) следует алгоритм метода Ньютона.

Jq=qk

В целом, метод ПКО выражает достаточно общую идеологию конструирования алгоритмов решения краевых задач и нелинейных уравнений. При использовании более простых корректирующих операторов М может наблюдаться замедление сходимости вычислений в окрестности решения. В этом случае целесообразно предусмотреть специальные меры для сокращения числа итераций. Для этого можно применить параметрическую модификацию метода ПКО либо перейти в ходе итераций к комбинированному использованию метода ПКО с другими методами, например, с одним из методов секущих [1].

Параметрическая модификация метода ПКО может быть представлена в виде

q(k+1) = M

Рт -ak¿Л

г=1

k = 1,2,3,

Рациональным выбором параметров ак можно добиться ускорения сходимости вычислительного процесса. Возможные способы определения этих параметров аналогичны рассмотренным в [6].

Для сокращения числа итераций можно использовать известные методы ускорения сходимости вычислительных процессов, основанные на линейной интерполяции и использовании асимптотических свойств линейно сходящихся последовательностей [1, 2, 5]. Однако эти условия могут нарушаться на одной из стадий вычислительного процесса, что может привести к его расходимости. При комбинировании метода ПКО с другими более предпочтительным является использование методов, основанных на конечно-разностной аппроксимации дХ

матрицы Якоби — в методе Ньютона (14), характерной для методов секущих [1, 4]. Извест-дq

но, что такие методы обеспечивают сверхлинейную сходимость и имеют порядок, по крайней мере, 1 +2^ ~ 1,618. Линейная интерполяция Xв Ят может привести к необходимости использования ряда других методов секущих [1, 3].

Проведенный анализ и накопленный опыт решения прикладных задач позволяют рекомендовать для комбинированного использования метода ПКО следующий способ конечно-

Гахт1

разностной аппроксимации оператора M = — на каждом шаге, начиная с некоторой итерации, который применительно к решению уравнения (12) определяется выражением

M =

dXn 1

0Л-1,

dq

Q = [AqbЛq2,..., Aqra]; л = [Xb^Xm]; = q,-+1-q,•

Значения Лqг и X,, необходимые для проведения вычислений, определяются по результатам первых т итераций, согласно основной процедуре метода ПКО (13). При переходе к новой итерации проводится циклический сдвиг элементов матриц Q и Л, при котором элементы Лql и Xi исключаются, все остальные элементы сдвигаются влево на одну позицию, а на место сдвинутых элементов Лqm и Xm ставятся элементы Лqm+l и Xm+i. Для ускорения

вычислений можно применить пошаговую аппроксимацию обратных матриц Л-1 на основе метода Шульца быстрых обращений [1, 4].

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

Пример. Рассмотрим особенности применения метода ПКО на примере решения задачи расчета импульсной программы полета космического аппарата из некоторого исходного состояния, определяемого значениями его фазовых координат х0, y0, zo, VXo, Vy0, VZq на момент

to, в требуемое конечное состояние xT, yT, zT за заданное время T.

Будем считать, что полет происходит в нормальном гравитационном поле Земли. C учетом [8] уравнения движения в абсолютной геоцентрической экваториальной системе отсчета представим в виде

x = Vx; у = Vy; z = Vz;

Vx = -ax; Vy = -ay; Vz = (2bc - a)z;

—3 2 -2 2 -2

a = b[a00 + c(d -1)]; b = R0r ; c = 1,5a20R0 r ; d = 5z r ;

J20 = -0,001 082 627; a00 = 62 564 951 м2/с2;

a20 =-67 889,273 м2/с2; R0 = 6371 км; r = (x2 + у2 + z2)1/2.

В качестве ПКО воспользуемся приближенным аналитическим решением этой задачи, соответствующим движению в однородном центральном гравитационном поле. В этом случае динамика объекта описывается упрощенными уравнениями:

r = V; V = -ш2г , где ш — угловая скорость орбитального движения спутника в момент t0;

Ш = Ш(Х0) = д/п)r0-3 ; r0 = д/xo + Уо + zo ; П0 = 398 600,44км3/с2.

Для такой модели задача импульсного полета имеет следующее решение:

rT - cosшТг0-ш-1б1пшТУ0

о V = ш-;

sin шТ

5У = (5УТ 5У )1/2; а = 5 V _15У,

где 5У — вектор требуемого импульса скорости; 5 У — его модуль; а — вектор направляющих косинусов импульса 5У.

Совокупность приведенных соотношений определяет значение ПКО М [•]. Решение исходной задачи производится в соответствии с методом ПКО:

к

Я(к+1} = М[гт -ЕАг°)(Т)]; Я = [5У,5У,а]т ; к = 1,2,...

г=1

На каждой итерации производится вычисление невязок Аг(г)(Т) путем однократного интегрирования приведенной выше полной системы дифференциальных уравнений движения объекта в нормальном гравитационном поле. При этом каждый раз изменяются начальные

условия интегрирования по значениям элементов вектора скорости У0к+1 = У0 + 5Ук.

Приведем некоторые результаты численных расчетов. Пусть в исходном состоянии космический аппарат находится на экваториальной круговой орбите с высотой к = 900 км и имеет следующие начальные значения фазовых координат:

Хо = 7243 км; у0 = 634 км; 2о = 0;

УХ =-0,645 км/с; Уу = 7,376 км/с; У, = 0.

Х0 ' ' у0 ' ' 20

Требуется определить импульсное управление, обеспечивающее попадание космического аппарата за время Т = 500 с в заданную точку Хт = 5762 км; ут = 4596 км; 2т = 0.

Результаты расчета импульсной программы управления по итерациям приведены в таблице. Здесь 5р(Т) — значение модуля координатного промаха относительно терминальной точки.

Номер итерации 5р(Т),км а Х а у 5УХ, м/с 5Уу, м/с 5У, м/с

0 165,446 -0,4317 0,9020 -515 1,075 1,192

1 22,415 -0,6534 0,7570 -848 0,982 1,297

2 2,806 -0,6304 0,7763 -801 0,986 1,270

3 0,348 -0,6333 0,7739 -807 0,986 1,274

4 0,043 -0,6330 0,7742 -806 0,986 1,273

Данные таблицы свидетельствуют о хорошей сходимости рассмотренной реализации метода ПКО.

Работа выполнена при поддержке РФФИ (проект № 09-08-00259).

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

1. Вержбицкий В. М. Основы численных методов. М.: Высш. школа, 2002. 840 с.

2. КрасносельскийМ. А. и др. Приближенное решение операторных уравнений. М.: Наука, 1969. 455 с.

3. Ортега Дж., Рейнболдт В. Итерационные методы решения нелинейных систем уравнений со многими неизвестными. М.: Мир, 1975. 560 с.

4. Островский А. М. Решение уравнений и систем уравнений. М.: ИЛ, 1963. 383 с.

5. ТраубД. Ф., ВожьняковскийX. Общая теория оптимальных алгоритмов. М.: Мир, 1983. 382 с.

6. Федоренко Р. П. Приближенное решение задач оптимального управления. М.: Наука, 1978. 488 с.

7. Миронов В. И. Конструктивный метод решения краевых задач управляемого движения // Алгоритмы и программы исследования систем управления. Вып. 6. Л.: ВИКИ им. А. Ф. Можайского, 1980. С. 70—74.

8. Космические траекторные измерения / Под ред. П. А. Агаджанова, В. Е. Дулевича, А. А. Коростелева. М.: Сов. радио, 1969. 504 с.

Вячеслав Иванович Миронов —

Юрий Вячеславович Миронов —

Рафаэль Мидхатович Юсупов —

Сведения об авторах д-р техн. наук, профессор; Санкт-Петербургский институт информатики и автоматизации РАН; E-mail: mironuv@yandex.ru д-р техн. наук; Санкт-Петербургский институт информатики и автоматизации РАН; ст. научный сотрудник; E-mail: mironuv@yandex.ru

д-р техн. наук, профессор; Санкт-Петербургский институт информатики и автоматизации РАН; член-корреспондент РАН; E-mail: spiiran@iias.spb.su

Рекомендована институтом

Поступила в редакцию 28.05.09 г.

УДК 621.629.78

В. Е. Косенко, Д. И. Марарескул, В. И. Ермоленко, В. И. Лавров, А. Н. Арапочкин, А. И. Косынкин, А. Б. Симонов

НАВИГАЦИОННОЕ ОБЕСПЕЧЕНИЕ ПОСАДКИ ВОЗДУШНЫХ СУДОВ С ПРИМЕНЕНИЕМ ГЛОНАСС-ТЕХНОЛОГИЙ

Рассматриваются подходы к реализации системы обеспечения захода на посадку и посадки воздушных судов с использованием ГЛОНАСС-технологий. Исследована возможность повышения точности навигационного обеспечения за счет совмещения дифференциальной станции с дополнительным наземным источником ГЛОНАСС-подобного сигнала.

Ключевые слова: спутниковая радионавигационная система, навигационная аппаратура потребителей, дифференциальный режим, контрольно-корректирующая станция, псевдоспутник.

Введение. Анализ международных требований к навигационному обеспечению полетов гражданской авиации показал, что эксплуатируемые спутниковые радионавигационные системы ГЛОНАСС и GPS удовлетворяют требованиям по точности навигационного определения воздушных судов при полетах по маршруту и в аэродромных зонах [1]. Вместе с тем проведенные исследования показали, что в классическом варианте дифференциальной навигации [2], когда воздушное судно (ВС) использует радионавигационные сигналы от спутников и дифференциальные поправки от наземной контрольно-корректирующей станции (ККС), возможно обеспечить требуемый уровень точности навигационных определений лишь для захода на посадку по первой (грубой) категории, и только по плановым координатам. При этом остается нерешенным вопрос о выполнении требований ICAO (International Civil Aviation Organization) по целостности системы, т.е. вероятности оповещения воздушного судна о возможных сбоях в работе той или иной радионавигационной точки за время, не превышающее заданного.

Если ККС дифференциальной подсистемы ГЛОНАСС совместить с наземным излучателем радионавигационного сигнала, аналогичного сигналам навигационных спутников, то за счет существенного улучшения геометрического фактора точность дифференциальных навигационных определений может быть доведена до требований второй и даже третьей