Сер. 10. 2009. Вып. 2
ВЕСТНИК САНКТ-ПЕТЕРБУРГСКОГО УНИВЕРСИТЕТА
УДК 656.022+004.021 М. В. Сотникова
АЛГОРИТМЫ ФОРМИРОВАНИЯ МАРШРУТОВ ДВИЖЕНИЯ СУДОВ С УЧЕТОМ ПРОГНОЗА ПОГОДНЫХ УСЛОВИЙ
1. Введение. Задача выбора маршрута судна имеет очень важное практическое значение, так как от ее решения зависит как время движения, которое на длительных трансокеанских переходах достигает нескольких суток, так и расход топлива, определяющий экономическую эффективность перевозок. На практике маршрут судна обычно прокладывается судоводителем на основе имеющейся карты района плавания, данных о погодных условиях, характеристиках судна и в значительной мере личного опыта. Естественно, что при таком эмпирическом подходе выбранный маршрут далеко не всегда является наилучшим по отношению к длительности движения и его экономичности. В связи с этим особую актуальность приобретает разработка алгоритмов автоматизированного решения задачи о выборе оптимальных маршрутов или маршрутов, которые близки к ним в определенном смысле.
В дальнейшем под маршрутом судна будем понимать траекторию движения его центра масс на земной поверхности и распределение значений линейных скоростей вдоль этой траектории.
В настоящей работе предложены формализованная постановка задачи выбора маршрута в виде задачи конечномерной оптимизации и алгоритм ее решения. Данная постановка учитывает, во-первых, долгосрочный прогноз погодных условий на район плавания, в том числе наличие опасных зон, запретных для движения судна, во-вторых, географические ограничения, определяемые береговыми линиями и зонами мелководья, в-третьих, динамические параметры судна, позволяющие оценить влияние ветра и волнения на скорость его движения.
Необходимо отметить, что поиск оптимального маршрута на длительных переходах является исключительно сложной задачей оптимизации. Это определяется высокой размерностью задачи и наличием большого количества ограничений. Принципиальную роль играют вычислительные затраты на поиск оптимального маршрута. В статье предлагается оригинальный подход к разработке алгоритмов формирования маршрутов, близких к оптимальным решениям, отражающий существенную ограниченность вычислительных затрат. Реализация подхода демонстрируется на конкретном примере.
Среди опубликованных работ имеется ряд статей, связанных с рассматриваемой тематикой. В частности, одним из альтернативных известных подходов к решению задачи выбора маршрута является метод изохрон, впервые предложенный в работе [1] и в дальнейшем развитый в [2]. Следует отметить, что формирование маршрутов
Сотникова Маргарита Викторовна — аспирантка кафедры компьютерных технологий и систем факультета прикладной математики-процессов управления Санкт-Петербургского государственного университета. Научный руководитель — докт. физ.-мат. наук, проф. Е. И. Веремей. Количество опубликованных работ: 5. Научное направление: оптимизация динамических систем с прогнозирующими моделями. E-mail: s_margosha@mail.ru.
(§М. В. Сотникова, 2009
на базе данного подхода требует введения дополнительных упрощающих предположений, поскольку в противном случае вычислительные алгоритмы становятся слишком громоздкими для практического применения.
Другой подход связан с непосредственным использованием идей классического вариационного исчисления [3]. Однако его применение для решения практических задач возможно только в случае предельно упрощенных постановок с примитивной математической моделью движения судна и без учета статических и динамических ограничений.
Для формирования маршрутов могут также применяться методы, базирующиеся на теории динамического программирования [4]. Но известная проблема «проклятия размерности» не позволяет рассматривать практические задачи построения маршрута с разумными вычислительными затратами при необходимости постоянного пересчета решения с поступлением уточненной информации о погодных данных.
2. Вариационная задача о формировании маршрута. Пусть заданы географические координаты начальной А (фо, Ао) и конечной В (ф\, Ах) точек траектории, где ф и А - широта и долгота соответственно. Простейшим примером задачи оптимизации маршрута является поиск траектории минимальной длины, соединяющей указанные точки на поверхности Земли.
Отметим, что без учета ограничений данная проблема представляет собой классический вариант простейшей основной задачи вариационного исчисления о геодезических линиях. Известно, что ее решением является ортодромия - дуга большого круга, проходящая через две заданные точки. Однако при учете запретных для траектории областей, определяемых островами, зонами мелководья и международными соглашениями, необходимо рассматривать соответствующую более сложную вариационную задачу при наличии ограничений.
Оптимизация маршрута предполагает не только нахождение траектории, но и скорости движения по ней, что приводит к необходимости учитывать как динамические особенности судна, так и воздействия на него внешней среды - в первую очередь ветра, волнения и океанских течений. В этом случае в рассмотрение вводятся соответствующие дифференциальные уравнения движения, что еще больше усложняет задачу оптимизации, превращая ее в вариационную задачу Лагранжа.
Для детальной постановки рассматриваемой задачи будем задавать маршрут судна с помощью двух функций р (Ъ) и V (Ъ), отражающих величины текущего курсового угла и заданной линейной скорости судна вдоль определенного курсом направления соответственно в момент времени Ъ. Курсовой угол отсчитывается от направления на север по часовой стрелке.
Далее будем различать заданную и фактическую скорости судна. Заданная скорость отражает пожелания судоводителя, реализуемые через однозначную связь с оборотами винта, и соответствует скорости хода судна на тихой воде. Естественно, что фактическая скорость зависит от заданной, однако не равна ей, поскольку на движение влияют внешние факторы, обусловленные погодными условиями. Фактическая скорость удовлетворяет следующей системе дифференциальных уравнений, которую можно трактовать как математическую модель движения судна:
'IV fl Ы, п (^ , Рситт, Ргю ауе: Р'шгис0 ,
ф = У*2 (ы,р) , (1)
А = Ь(Ы: р) •
Здесь V - фактическая скорость, п (V) - зависимость частоты вращения гребного винта от заданной скорости на тихой воде, ф - широта, А - долгота, Рсигг = Рсигг (Ъ, Ф, А),
Pwave = Pwave (Ъ, ф, А), Р^ив, = Pwind (Ъ, Ф, А) - внешние силы, определяющие воздействие течения, волнения и ветра соответственно на фактическую скорость движения, которые считаются известными исходя из долгосрочного прогноза погоды.
Обозначим в качестве По район плавания судна, включающий начальную и конечную точки траектории, трактуя его как сферический прямоугольник:
По = {(ф, А) | Ао < А < Ах, фо < ф < Ф1} •
Прогноз погоды задается для интервала времени [То,Т\], включающего моменты отправления судна из начальной точки и прибытия в конечную точку, т. е. на весь период плавания.
В каждый момент времени Ъ € [То,Т] на интервале задания прогноза и в любой точке (ф, А) € По, принадлежащей району плавания, в качестве исходных данных доступны значения следующих параметров:
а) ветра (направление и скорость);
б) волнения (высота волны, период, направление);
в) течения (направление и скорость).
Отметим, что первое уравнение в системе (1) описывает динамику движения судна с учетом погодных условий. Для его составления необходимо задать ряд технических параметров судна, включая водоизмещение, длину, ширину, осадку и др. Оставшиеся два уравнения системы (1) определяют динамику изменения географических координат с учетом известных геодезических соотношений [5].
Пусть £о - начальный момент времени, т. е. время отправления судна из точки А (фо,Ао). Тогда начальными условиями для системы (1) будут
V (Ъо) = vо, ф (Ьо) = фо, А (Ьо) = Ао- (2)
Можно считать, что в начальный момент заданная и фактическая скорости совпадают. Из системы (1), (2) следует, что функции р (Ъ) и V (Ъ) однозначно задают траекторию и скорость движения по ней, т. е. маршрут движения судна.
Пусть Ъх - нефиксированный конечный момент времени, зависящий от выбора функций р (Ъ) и V (Ъ). Тогда терминальное ограничение, соответствующее попаданию в момент времени Ъх в окрестность конечной точки В (фх,Ах), определяется неравенством
(ф1 _ ф (Ъ1)) +(А1 _ А (Ъ1)) < £, Ъ1 ^ ^, (3)
в котором £ > 0 - заданная малая величина. Причем время ^ ограничено конечным моментом задания прогноза Т1.
Помимо ограничения (3) на выбор функций р (Ъ) и V (Ъ) налагаются дополнительно:
1) ограничения на величину заданной скорости
vmin ^ v ^ vmax,
Vг € [Ъо,^, (4)
где vmax и vmin - максимальная и минимальная скорости судна на тихой воде соответственно;
2) статические ограничения, связанные с зонами мелководья, международными соглашениями, островами и береговыми линиями;
3) динамические ограничения, обусловленные запретными для движения судна зонами с опасными погодными условиями.
Каждое из статических ограничений представляет собой ограниченное замкнутое связное множество точек на земной поверхности. Пересечение данных множеств с районом плавания По будем обозначать П® , к = 1, Ж. Тогда допустимая по отношению к статическим ограничениям область плавания определяется множеством точек
N
Оя = По\и П9к ■
к = 1
Следовательно, функции р (г) и V (г) должны обеспечивать выполнение условия
(Ф (г) ,л (г)) е о,я, шг € [го,г1 (5)
в каждый момент времени г.
Статические ограничения могут быть заданы набором замкнутых контуров, представленных наборами точек с соответствующими географическими координатами. Каждый такой контур ограничивает зону, запретную для движения судна. Например, замкнутый контур Г°к, являющийся границей множества 0°к, может быть задан в следующем виде:
Гдк = {(Ф1,Х1), (ф2,Л2) ,■■■, (фп,Л„)} ,
причем точки задаются последовательно, в соответствии с направлением обхода либо по часовой стрелке, либо против нее. Для любой отдельно рассматриваемой задачи статические ограничения задаются не на весь земной шар, а на зону плавания Оо.
Динамические ограничения требуют, чтобы точка (ф (г), Л (г)) не попадала в запретные по погоде области в любой момент г е [го, £1]. Необходимо отметить, что динамические ограничения изменяются с течением времени в соответствии с прогнозом погоды. Пусть П“ (£),;?’ € 1, М, - опасная зона (в пересечении с районом плавания По)- В каждый фиксированный момент времени она представляет собой замкнутое ограниченное связное множество. Аналогично допустимому множеству Оя, введем в рассмотрение допустимое по погодным условиям множество, зависящее от времени:
м
Ош (г) = Оо\ у оа (г).
3=1
Следовательно, в каждый момент времени должно выполняться условие
(Ф (г), Л (г)) е Ош (г), шг е [го,г1 ■ (6)
В формализованном виде опасные зоны задаются аналогично статическим ограничениям.
Для характеристики качества конкретного выбранного маршрута введем в рассмотрение следующие функционалы:
1) время движения из начальной точки в конечную точку
Ъ (р (■) ^(')) = ! А = г1 (р (г), V (г)) - го; (7)
tо
2) расход топлива
(Р (') ^(')) = У Р (Р (г)) ^, (8)
tо
где Р (V (г)) - расход топлива в единицу времени.
Таким образом, задача о поиске оптимального маршрута может быть представлена в виде следующей вариационной задачи Лагранжа на условный экстремум: необходимо найти такие функции р (г) и V (г), которые доставляют минимум функционалу (7) или
(8) при наличии связи (1) с учетом граничных условий (2), (3) и ограничений (4)-(6).
Поставленная задача представляется исключительно сложной в вычислительном плане прежде всего в связи с тем, что отсутствуют явные аналитические представления почти для всех функций, входящих в состав указанных функционалов, связей, условий и ограничений. Их задание в основном осуществляется алгоритмически, что крайне затрудняет поиск оптимального решения методами, базирующимися на необходимых условиях экстремума.
Альтернативным вариантом могут служить прямые методы вариационного исчисления, но при этом необходимо интегрировать систему дифференциальных уравнений (1) на всем времени перехода, что порождает значительные вычислительные трудности.
Дополнительные проблемы связаны с тем, что при формировании маршрутов приходится постоянно учитывать поступление новых прогнозов погодных условий и повторять решение поставленной вариационной задачи заново через фиксированный промежуток времени. Это налагает существенные ограничения на допустимое время проведения соответствующих расчетов.
Отмеченные обстоятельства определяют целесообразность в отказе от непосредственного поиска точного решения с переходом к различным вариантам постановок задач приближенной оптимизации маршрута. Одним из наиболее естественных путей подобного перехода служит схема существенной дискретизации задачи, которая по существу соответствует практике судовождения. Естественно, что для подтверждения состоятельности получаемых при этом приближенных решений необходимо сформировать систему нижних оценок для функционалов времени и расхода топлива.
3. Задача оптимизации маршрута в конечномерном варианте. Одним из существенных моментов, определяющих целесообразность и возможность перехода к дискретному варианту задачи формирования маршрута, служит тот факт, что используемый прогноз погодных условий имеет естественное дискретное представление во времени и в пространстве.
Прежде всего прогноз задается для узлов сетки, покрывающей зону плавания Оо на земной поверхности и имеющей фиксированный шаг по широте Дф и долготе ДЛ. Погода в узлах сетки задается с фиксированным временным шагом ДТ, т. е. в моменты То,То + ДТ, ■■■,Т1. В каждом узле сетки доступны перечисленные выше параметры погодных условий.
При необходимости информация о погоде в остальных точках района плавания в любой фиксированный момент времени может быть получена путем линейной интерполяции значений в узлах сетки и между двумя ближайшими моментами времени задания прогноза.
Так же как и прогноз, опасные зоны по погодным условиям задаются на район плавания Оо с фиксированным временным шагом ДТ и с возможностью линейной интерполяции.
Будем полагать в дальнейшем, что траектория движения состоит из конечного числа локсодромических участков, т. е. частей траектории с постоянным значением курсового угла (рис. 1). На рисунке приняты следующие обозначения: А (ф0, Ао) и В (фі, Аі) -начальная и конечная точки траектории, М(ф,А) - текущая позиция судна с географическими координатами ф и Л, М^, і Е 1,р + 1, - точки поворота траектории, причем М\ = А и Мр+1 = Б, ер і , Т^, г Е 1,р, - значения курсового угла, длины и заданной
скорости на г-м локсодромическом участке.
л В
(фр ■> Зр)>
/// (ф/Л), V;-м ///
Рис. 1. Траектория движения и текущее положение судна
Пусть задано общее число локсодромических участков траектории р. Тогда в качестве варьируемых переменных для задачи конечномерной оптимизации принимаются следующие величины:
1. Группа переменных, определяющая траекторию движения (<^, *%), г Е 1,р. Отметим, что число независимых переменных здесь равно 2 • (р — 1), поскольку длина и курсовой угол последнего участка однозначно определяются положением предпоследней точки Мр. Введем в рассмотрение вектор переменных
г = { (Р1,51), (Р2,52), (рр,БР)}еЕ2р,
задающих траекторию движения, которую будем обозначать 7 (г), определяя ее непосредственную зависимость от выбора вектора г. Траектория 7 (г) содержит множество точек (ф, Л), принадлежащих кривой, состоящей из участков, зависимых от вектора г. При этом будем предполагать, что точки А (фо, Ло) и В (ф1, Л1) принадлежат траектории 7 (г) для любого г.
2. Переменные, определяющие величину заданной скорости на каждом участке траектории У^ г Е 1,р. Введем вектор переменных, характеризующих скорость движения по траектории:
V = (У1 V,..., Ур) Е Ер.
Пара совместно рассматриваемых векторов г и V однозначно устанавливает маршрут движения судна.
N А
Очевидно, что первое уравнение системы (1) имеет положение равновесия, если заданная скорость и внешние силы постоянны на определенном промежутке времени. Данное положение равновесия ш = находится из уравнения
fl (£, ш, п ('У ) , Рсигг , Pwave, Р'шъпс1) 0 (9)
и приближенно соответствует фактической скорости движения, если внешние силы мало изменяются на рассматриваемом промежутке времени. Следовательно, уравнение (9) можно использовать для расчета фактической скорости вместо интегрирования системы дифференциальных уравнений (1). При этом фактическая скорость на траектории, состоящей из нескольких участков, должна вычисляться заново в соответствии с (9) в точках поворота траектории и на локсодромических участках с фиксированным временным шагом смены прогноза. Отметим, что уравнение (9) имеет смысл только в точках, принадлежащих допустимой области .
Пусть известны вектор г, определяющий траекторию 7 (г), и вектор V заданных скоростей движения по траектории. Причем будем считать, что 7 (г) С , т. е. траектория допустима по отношению к статическим ограничениям.
Сформируем алгоритм вычисления времени движения 1т = 1т (г, V) от начальной точки до момента попадания в заданную окрестность конечной точки и для расчета суммарного расхода ,1р = 1р (г, V) топлива для выбранного маршрута. При реализации алгоритма будем считать, что расход топлива в единицу времени определяется известной функцией ц = ц (У) для конкретной силовой установки судна, где У - заданная скорость хода.
Пусть в начальный момент времени £ = £о судно находится в точке А (фо,Ло). Обозначим через £ время, оставшееся до очередного момента смены прогноза.
Алгоритм вычисления состоит из следующих этапов:
1. Рассчитывается фактическая скорость ш на основании уравнения (9) при известных заданной скорости У1, курсовом угле рч и внешних силах, соответствующих моменту £ = £о, на данном участке. Положим, что при этом 1т = 0 и ,1р = 0.
2. Определяется расстояние, которое сможет пройти судно с полученной фактической скоростью ш до момента смены прогноза 5 = шЬ. Вычисляется расстояние Ь от текущей точки до следующей точки поворота траектории и сравниваются величины 5 и Ь.
3. Если 5 < Ь, то определяется следующее положение точки М(ф, Л) на траектории в соответствии с формулами прямой геодезической задачи. При этом точка М(Ф, Л) остается на текущем участке г траектории. Время £ попадания в новую точку соответствует моменту смены прогноза и вычисляется по формуле £ = £ + £. Пересчитываются значения суммарного времени в пути 1т = 1т + Ь и расхода топлива ,1р = ,1р + Ъц (V,). Для следующего этапа алгоритма полагается £ = АТ, где АТ - шаг задания прогноза, и определяется заново фактическая скорость ш в новой точке для момента времени £, заданной скорости У, и курсового угла р, из уравнения (9). Затем осуществляется переход к п. 2.
4. Если 5 ^ Ь, то в качестве следующего положения точки М(ф,Л) принимается точка поворота траектории. В данном случае точка М(ф,Л) переходит с текущего участка траектории г на следующий участок г + 1. Обозначим через V' = ^ время, необходимое для попадания в новую точку. Тогда время £ попадания в новую точку находится по формуле £ = £ + г. Далее пересчитываются значения суммарного времени в пути 1т = 1т + Ь и расхода топлива ,1р = ,1р + Ь'ц (У,). Если данная точка
поворота траектории совпадает с точкой В (ф1, Л1) либо находится в ее заданной малой окрестности, то вычисления заканчиваются. При этом время попадания в окрестность конечной точки равно текущему моменту £1 = £. Если нет, то для следующего этапа алгоритма берется время £ = £ — £ ', оставшееся до смены прогноза, и рассчитывается заново фактическая скорость ш в новой точке для момента времени £, заданной скорости У,+ 1 и курсового угла р,+ 1 по уравнению (9). Далее осуществляется переход к п. 2 и вычисления повторяются.
Таким образом, приведенный алгоритм позволяет однозначно определить величины функционалов 1т = 1т (г, V) и ,1р = ,1р (г, V) для заданного маршрута движения.
Для постановки оптимизационных задач, связанных с указанными функционалами, преобразуем описанные выше ограничения (3)-(6), относящиеся к вариационному варианту, к соответствующей дискретной форме.
Прежде всего заметим, что терминальное ограничение (3) с учетом принятых варьируемых переменных может быть представлено в виде
(ф1 — ф (£1, г, V))2 + (Л1 — Л (£1, г, V))2 <£, £1 < Т1. (10)
Суть данного ограничения состоит в том, что текущая точка М(ф, Л) должна попасть в окрестность конечного положения В (ф1 ,Л1) не позднее конечного момента задания прогноза Т1 .
Ограничение (4) на величину заданной скорости соответствует следующим ограничениям на переменные, определяющие заданную скорость движения на каждом участке траектории:
^тт ^ ^ ^ %ахз ^ = 1 :Р- (И)
Теперь рассмотрим ограничение (5), которое можно записать так:
7 (г) С , (12)
поскольку траектория движения 7 (г) полностью определяется заданием вектора г. Заметим, что для проверки допустимости траектории, т. е. выполнения включения (12), необходимо установить пересекается ли она с контурами, ограничивающими запретные статические области.
Наконец, перейдем к рассмотрению ограничений, обусловленных погодными условиями. Положение точки М (ф (£, г, V) ,Л(£, г, V)) на траектории в момент времени £ зависит от выбора траектории 7 (г) и вектора заданных скоростей V. Следовательно, определяемое опасными погодными условиями ограничение (6) имеет вид
М (ф (£, г, V), Л (£, г, V)) С Пш (£), У£ Е [£о,£] ■ (13)
Для проверки ограничения (13) следует рассмотреть функцию Та1 (г, V), значение которой равно времени нахождения в опасных зонах для выбранного маршрута (г, V). Если для фиксированного маршрута значение данной функции равно нулю, то он удовлетворяет ограничению (13), которое в этом случае можно представить в форме
Та1 (г, v)=0■ (14)
Будем считать, что момент отправления принадлежит первому шагу задания прогноза £о Е [То, То + АТ]. На каждом шаге задания прогноза опасные зоны определяются набором ограничивающих их контуров. Следовательно, для каждого шага задания
прогноза можно установить набор сегментов, которые получаются в пересечении траектории движения 7 (г) с опасными зонами.
Вычисление суммарного времени нахождения судна в опасных зонах для заданного маршрута (г, V) можно осуществить с помощью следующего алгоритма:
1. Полагается, что начальное значение суммарного времени нахождения в опасных зонах нулевое: Та1 (г, V) = 0.
2. Для г-го шага задания прогноза, соответствующего интервалу времени [То + (г — 1) АТ, То + гАТ], определяется сегмент траектории, который проходит судно на данном интервале времени. Этот сегмент может содержать точки траектории, принадлежащие нескольким соседним участкам.
3. Устанавливается набор сегментов, соответствующих пересечению траектории 7 (г) с опасными зонами на указанном интервале времени.
4. Определяется, пересекается ли построенный на шаге 2 сегмент с полученным на-
бором сегментов. Если есть пересечение, то рассчитывается время нахождения в опасных зонах по формуле Та1 (г, V) = Та1 (г, V) + £а1, где £а1 ^ АТ - время нахождения
в опасной зоне на данном шаге.
5. Если пересечения нет и сегмент траектории, проходимый судном на данном шаге прогноза, включает конечную точку или конечная точка сегмента находится в окрестности точки В (ф1, Л1), то вычисления заканчиваются. В противном случае увеличивается номер шага г = г + 1 и выполняется переход к п. 2.
На базе ограничений (10)—(12) и (14) можно сформировать допустимое множество П Е Е3р в пространстве переменных (г, V) Е Е3р для конечномерной задачи оптимизации в виде
о_ / (г, V) I Р (£1, г, V) <£,£1 < Т1; 7 (г) С Пй
; ^ш1п ^ ^ ^ ^шахз ^ 1 \ /1 к \
П^ Та1 (г, V)=0■ I (15)
В (15) функция р (£1, г, V) определяется в соответствии с ограничением (10).
Если множество П пусто, то не существует маршрута, обеспечивающего выполнение всех ограничений. В противном случае на множестве допустимых маршрутов необходимо выбрать тот, который доставляет минимум функционалу (7) или (8).
Таким образом, задачу формирования маршрута можно представить как задачу конечномерной оптимизации
1т (г, V) —> шт (16)
(г^)ЕПСЕ3р
в случае минимизации времени перехода или
1р (г, V) —> шт (17)
(г^)ЕПСЕ3р
в случае оптимизации расхода топлива.
Задачи (16), (17) являются очень сложными для непосредственного привлечения численных методов поиска решения. Это в первую очередь связано с высокой размерностью варьируемых векторов, включающих две группы разнородных переменных. Кроме того, сложность вызвана и алгоритмическим заданием существенно нелинейных функций, определяющих критерии и допустимые множества для оптимизации.
Для упрощения решения задач формирования маршрута разделим переменные, задающие траекторию и скорости движения по ней, и изменим постановки задач оптимизации.
Пусть множество допустимых заданных скоростей зависит от выбора вектора г Е Е2р, определяющего траекторию, и имеет вид
V (г) = {V € Ер\р(г1,г,у) < е, г-! < Т\; г;т;п < У < утах, г е 1 ,р; Таг (г, у) = 0} . (18)
Введем допустимое множество
Ид = {г Е Е2р| 7 (г) С Пд, V (г) = 0} . (19)
Заметим, что множество Ид состоит из векторов г, определяющих траектории, по которым возможно движение без нарушения ограничений, формирующих множество V (г). На множестве векторов г Е Ид зададим функцию
Т (г) = шт 1т (г, V), (20)
у£У(г)
где 1т (г, V) вычисляется по приведенному выше алгоритму. Значением функции Т (г) для вектора г является минимально возможное время перехода от начальной до конечной точки с учетом допустимого множества заданных скоростей V (г).
Следовательно, вопрос о минимизации времени перехода в рамках принятых соображений может быть поставлен в следующей форме:
Т (г) —> шт ,
или
1т (г*, V*) = шт шт 1т (г, V), (21)
ге^ у£У(г)
здесь (г*, V*) - обозначение для оптимального решения. Очевидно, что задача (21) может быть трактована как некоторое приближение к решению задачи (16) о минимальном времени движения.
Аналогично, по отношению к задаче (17) введем в рассмотрение функцию
Е (г) = шт 1р (г, V).
у£У(г)
Значением функции Е (г) для вектора г является минимально возможный расход топлива за время перехода от начальной до конечной точки с учетом допустимого множества заданных скоростей V (г). Соответственно задачу о минимальном расходе топлива можно в данном случае трактовать как
Е (г) —> шт ,
ге^
или
1р (г**, V**) = шт шт 1р (г, V), (22)
геБ^ у£У(г)
где (г**, V**) - оптимальное решение.
Дополнительное обоснование целесообразности перехода к задачам (21), (22) состоит в том, что целевые функции и ограничения для задач (16), (17) определены только на допустимых траекториях, т. е. только для тех векторов г Е Е2р, для которых выполнено условие 7 (г) Е Пд. Отсюда следует, что для этих задач можно применять только такие численные методы, которые генерируют минимизирующие последовательности векторов {г^}, к = 0,1, 2,..., для которых выполняются условия
7 (г^) Е Пд, к = 0,1, 2,... . При этом обязательным этапом является решение отдельной оптимизационной задачи о поиске начального вектора го Е Е2р такого, что 7 (го) Е Пд.
В отличие от указанной ситуации, оптимизационные задачи в вариантах (21), (22) свободны от подобных трудностей, поскольку на внутренних шагах численного метода можно использовать и векторы г Е Е2р, для которых 7 (г) Е Пд. Это достигается разделением переменных, позволяющим для любой недопустимой траектории 7 (г) строить ближайшую к ней траекторию 7 (г') Е Пд, что существенно повышает эффективность поиска.
Таким образом, задачи выбора маршрутов с уменьшением времени движения или с экономией расхода топлива далее будем рассматривать в вариантах (21) и (22) на допустимых множествах (18), (19).
4. Алгоритмы формирования маршрутов движения судна. Для решения задач (21), (22) могут применяться комбинации различных методов конечномерной оптимизации, например метод деформируемого многогранника для оптимизации по переменным, определяющим траекторию, и метод последовательного квадратичного программирования для оптимального выбора распределения скоростей вдоль траектории.
Однако вычислительные затраты, необходимые для решения задачи формирования маршрутов трансокеанских переходов при наличии сложных статических и динамических ограничений, в этом случае оказываются очень большими. Поэтому данный подход практически применим только в случае простых статических и динамических ограничений для переходов с относительно малыми длительностями. В частности, если отбросить динамические ограничения, его можно использовать для построения нижних оценок функционалов времени и расхода топлива.
Приводимые ниже алгоритмы базируются на представлении допустимых множеств в задачах (21), (22) конечными совокупностями допустимых траекторий. Ясно, что при этом не гарантируется достижение экстремума, однако можно получить определенное приближение к нему с вполне реальными вычислительными затратами на формирование маршрута.
Таким образом, для реализации предлагаемых алгоритмов вначале необходимо выполнить предварительный этап вычислений, состоящий в построении набора траекторий, покрывающих область плавания и допустимых по отношению к статическим ограничениям. Далее на каждой траектории из набора с помощью соответствующего алгоритма решается задача оптимизации распределения скоростей в варианте (21) или (22).
Последовательно рассмотрим алгоритмы построения наборов траекторий и оптимизации распределения скоростей.
4.1. Алгоритм построения набора траекторий. Очевидно, что задача о формировании допустимых траекторий движения либо вообще не имеет решения, либо имеет бесчисленное множество таких решений, причем последний вариант наиболее естественный, поскольку районы плавания задаются достаточно опытными судоводителями.
Отсюда следует и большое разнообразие возможных методов построения наборов траекторий. Самый простой способ, рассматриваемый ниже, заключается в изначальном равномерном геометрическом покрытии района плавания с последующей корректировкой траекторий на допустимость с учетом статических ограничений. Другие более
сложные варианты построения могут, например, дополнительно учитывать известный прогноз погоды для района плавания.
Исходно задаваемыми параметрами алгоритма являются: число траекторий в наборе, максимальная длина одного участка траектории и ширина набора траекторий, определяющая охватываемую данным набором область. Алгоритм состоит из следующих этапов:
1) формируется траектория, соединяющая начальную и конечную точки по линии прямого курса;
2) строится набор траекторий, расположенных выше и ниже траектории прямого курса. Построение осуществляется таким образом, чтобы заполнить всю область между крайними траекториями в определенном смысле равномерно;
3) проверяется допустимость каждой траектории из полученного в п. 2 набора. При этом учитывается, что каждой из них соответствует вектор г* Е Е2р, г е 1,ЛГ, где N - число траекторий в наборе, р - число участков каждой траектории. Иными словами, проверяется выполнение условий 7 (г*) € Пд, * = 1, Ы: если траектория 7 (г*) допустима, то вектор г* оставляется без изменений. В противном случае строится ближайшая допустимая траектория и вектор г* заменяется соответствующим вектором г*, определяющим новую траекторию 7 (г*).
В результате на выходе алгоритма получается набор допустимых по отношению к статическим ограничениям траекторий {7 (г*)}^=1 и соответствующих им векторов
Гг, г = 1,М.
4.2. Алгоритм формирования маршрута, оптимального по времени перехода. Исходными данными для алгоритма является набор траекторий {7 (г* )}^=1, допустимых по отношению к статическим ограничениям, т. е. удовлетворяющих условиям 7 (г*) € Пд, * = 1, N.
Далее на каждой траектории 7 (г*) Е Пд, определяемой вектором г*, решаем вспомогательную задачу (20) о минимизации времени перехода на допустимом множестве
V (г*). Это задача нелинейного программирования по отношению как к целевой функции, так и к ограничениям. Как показывает опыт, наиболее эффективным для ее решения является стандартный метод последовательного квадратичного программирования (SQP), реализованный, например, в пакете МА^АВ.
Однако заметим, что решение вспомогательной задачи существует не для любой траектории 7 (г), поскольку допустимое множество заданных скоростей V (г) может быть пусто. Например, в сложных погодных условиях может оказаться, что по данной траектории нельзя пройти, не попав в опасную погодную зону.
Если допустимое множество V (г*) не пустое, то в результате решения получаем оптимальный вектор V* Е V (г*) и соответствующее значение времени перехода ,1т (г*, V*). При этом маршрут (г*, V*) допустим.
Теперь можно сформировать приближение к решению задачи (21), т. е. выбрать в качестве предлагаемого оптимального маршрута (г*, V*) допустимый маршрут с наименьшим временем перехода Лт (г*, V*) = шт Лт (г*, V*).
*е1 С{1,...,М }
Если допустимых маршрутов нет, то алгоритм представляет в качестве решения вариант с минимальным временем нахождения в опасной зоне. Для поиска допустимых маршрутов в подобной ситуации можно либо каким-либо образом изменить исходный набор траекторий, либо ослабить ограничения, определяющие опасные области по отношению к прогнозу погодных условий.
4.3. Алгоритм формирования маршрута, оптимального по расходу топлива. Для приближения к решению задачи (22) можно, в принципе, применить алгоритм решения задачи (21), представленный выше, однако с введением дополнительного ограничения на величину расхода топлива. При этом используются те же исходные данные, что и для предшествующего алгоритма, т. е. набор допустимых траекторий {7 (г*)}^1 и соответствующих им векторов г*, г = 1, N. Кроме того, задается также величина J_max, определяющая ограничение на допустимый расход топлива.
Для реализации алгоритма введем новое допустимое множество заданных скоростей на траектории 7 (г):
V* (г) = {v| v е V (г), JF (г, v) < Jmx} ,
учитывающее дополнительное ограничение. Любой вектор из этого множества определяет распределение скоростей на траектории 7 (г), при котором выполняются все ограничения, формирующие множество V (г), и, кроме того, расход топлива не превосходит заданного предела J^ax.
Теперь для каждой траектории 7 (г*) из исходного набора можно решить задачу о минимизации времени движения, аналогичную задаче (20), но на допустимом множестве V* (г*).
Естественно, что множество V* (г) может оказаться пустым, например, если по траектории y (г) невозможно осуществить переход, потратив меньшее количество топлива, чем заданное в ограничении, или если, в силу сложных погодных условий, время в опасной зоне не удается сделать нулевым.
Если множество V* (г*) не пустое, то в результате решения задачи минимизации находится оптимальный вектор v**. Соответствующий маршрут является допустимым (г*, v**) е П, и для него можно вычислить соответствующее время перехода Jt (г*, v**).
Анализируя весь набор траекторий, можно сформировать приближение к решению задачи (22). В качестве оптимального маршрута (г**, v**) выбирается допустимый маршрут с наименьшим временем перехода JT (г**, v**) = min JT (г*, v**).
Таким образом, результатом работы алгоритма является наилучший среди допустимых маршрут, для которого время перехода наименьшее при расходе топлива, не превышающем заданного предела.
Если допустимых маршрутов нет, то алгоритм представляет в качестве решения вариант с минимальным временем нахождения в опасной зоне или с наименьшим расходом топлива. Для поиска допустимых маршрутов в подобной ситуации можно либо каким-либо образом изменить исходный набор траекторий, либо ослабить ограничения, определяющие опасные погодные области, либо увеличить значение допустимого расхода топлива.
5. Пример решения задачи формирования маршрута движения. Пусть начальная и конечная точки имеют координаты A (52.5, —12.1) и B (8.4, —57.9) соответственно. На рис. 2 представлены траектории, соединяющие начальную и конечную точки вдоль дуги большого круга (Great Circle) и по линии прямого курса (Rhumb Line). Длина данных траекторий составляет 6395.5 и 6446.4 км соответственно.
Для района плавания, показанного на рис. 2, заданы прогноз погоды, статические (Static) и динамические (Dynamic) ограничения. На рисунке представлены динамические ограничения, соответствующие первому шагу задания прогноза. Существенные для данного примера статические ограничения определяются, в частности, Азорскими островами, расположенными вблизи траектории прямого курса.
+30с
+15
Рис. 2. Траектории вдоль дуги большого круга и линии прямого курса
Прогнозируемые погодные данные задаются на 10 суток с момента отправления судна с шагом АТ = 6 ч. Опасными зонами плавания считаются те, в которых скорость ветра превышает 11 м/с (21.4 узла).
По доступным параметрам рассматриваемого судна (контейнеровоз водоизмещением 40 000 т) может быть составлено уравнение
(9) для определения его фактической скорости движения. Максимальная величина заданной скорости составляет vmax = 24 узла, минимальная - vmin = 5 узлов.
На рис. 3, а показан набор допустимых траекторий, построенный в соответствии с алгоритмом, приведенным в п. 4.1. Набор состоит из 21-й траектории, каждая из которых имеет 24 участка.
Далее с целью определения оптимального по времени маршрута для каждой траектории решается задача минимизации времени перехода на допустимом множестве заданных скоростей. Ввиду сложных динамических ограничений, решение этой задачи существует только для 11 траекторий из набора, показанных на рис. 3, б. Ниже приведен вектор JT, компонентами которого служат значения времени перехода в сутках для соответствующих 11 допустимых маршрутов:
Jt = [ 7.62 6.39 6.87 8.52 7.90 6.24 6.94 6.29 8.15 6.33 7.34].
Видно, что оптимальному маршруту соответствует минимальное время перехода JT* = 6.24 суток. Данный маршрут (Time optimal route) особо выделен на рис. 3, б. Приведем вектор v* распределения заданных скоростей для оптимального маршрута
* Г 20.00 20.00 20.08 24.00 24.00 24.00 24.00 24.00 24.00 24.00 24.00 24.00
V = \ 24.00 24.00 24.00 24.00 24.00 24.00 24.00 24.00 24.00 24.00 24.00 24.00
для которого расход топлива составляет JF = 750.83 т.
Далее формируется маршрут, оптимальный по отношению к расходу топлива. Для этого сначала вводится дополнительное ограничение на выбор заданных скоростей, определяемое максимально допустимым расходом топлива на переходе. В качестве величины Jmax возьмем значение, соответствующее экономии в 5% от расхода топлива для оптимального по времени маршрута, что составляет Jjmax = 713.3 т.
Далее для каждой траектории из исходного набора решается задача минимизации времени перехода на сформированном более узком допустимом множестве. В качестве оптимального по отношению к расходу топлива выбирается допустимый маршрут, для которого время перехода наименьшее. В данном случае время перехода для оптимального маршрута JT** = 6.45 суток, расход топлива на переходе J** = 709.3 т. Вектор распределения скоростей v**, соответствующий оптимальному маршруту (Fuel optimal route), имеет вид
17.46 23.25 17.58 20.24 23.25 23.64 23.27 23.27 23.72 23.35 24.00 23.31 21.43 23.33 22.63 23.28 23.96 23.31 23.31 23.33 23.68 23.72 23.25 23.35
Рис. 3. Набор допустимых траекторий (а) и допустимые маршруты в задаче минимизации времени перехода (б)
Этот маршрут также представлен на рис. 3, б.
Из полученных результатов следует, что экономии топлива в 5% соответствует увеличение времени в пути на 3.4%.
6. Заключение. В работе рассмотрен вопрос формирования маршрутов движения судна, оптимальных по отношению ко времени перехода или расходу топлива. Приведена постановка проблемы в форме соответствующей задачи вариационного исчисления и осуществлено ее сведение к задаче конечномерной оптимизации.
Основные трудности в привлечении известных численных методов для поиска решения этих задач порождаются исключительно сложной структурой допустимых множеств, отражающей учет погодных условий и геометрические ограничения в районе плавания. Кроме того, для трансокеанских переходов крайне значимой является высокая размерность конечномерного варианта.
Указанные трудности определяют практически недостижимые вычислительные затраты на поиск решения: особую роль играет существенная ограниченность времени счета, связанная с оперативной обстановкой на маршруте, в частности с периодичностью обновления прогноза погоды. Для выхода из данной ситуации предложены упрощенные варианты постановок соответствующих оптимизационных задач, решения которых можно трактовать как определенное приближение к искомым оптимальным маршрутам.
Разработаны алгоритмы численного решения таких задач, заранее ориентированные на ограниченность вычислительных ресурсов, выделяемых на их реализацию. В любом случае эти алгоритмы приводят к допустимым маршрутам, если они существуют в рамках принятого подхода. При пустоте допустимых множеств указываются пути ослабления ограничений, учитываемых при решении задачи.
Практическая применимость и эффективность предлагаемого подхода проиллюстрированы на конкретном примере.
1. James R. W. Application of wave forecast to marine navigation. Washington, D.C.: US Navy Hydrographic Office, 1957. 78 p.
2. Hagiwara H. Weather routing of (sail-assisted) motor vessels: PhD Thesis. Delft: Technical University of Delft, 1989. 337 p.
3. Bijlsma S. J. A Computational Method for the Solution of Optimal Control Problems in Ship Routing // Navigation. Journal of the Institute of Navigation. 2001. Vol. 48. P. 145-154.
4. Chen H. A dynamic program for minimum cost ship routing under uncertainty: PhD Thesis.
Massachusetts: Massachusetts Institute of Technology, 1978. 163 p.
5. Закатов П. С. Курс высшей геодезии. М.: Недра, 1976. 511 с.
Статья рекомендована к печати проф. Е. И. Веремеем.
Статья принята к печати 25 декабря 2008 г.