Введение и постановка задачи

Свойства множеств достижимости в нелинейных системах с интегральными ограничениями исследовались в работах [1,2]. Алгоритмы построения множеств достижимости, основанные на дискретных аппроксимациях, изучались в [3,4]. В работе [5] для систем с интегральными ограничениями на управление доказано, что любое управление, удовлетворяющее интегральным ограничениям и переводящее систему на границу множества достижимости, является локальным решением некоторой задачи оптимального управления.

В настоящей статье результаты из [5] обобщаются на следующий случай: 1) рассматриваются совместные интегральные ограничения на управление и траекторию; 2) изучаются граничные точки проекции множества достижимости на заданное линейное подпространство.

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

Далее будем использовать следующие обозначения. Для вещественной матрицы A через A мы обозначаем транспонированную матрицу, 0 обозначает нулевой вектор подходящей размерности, O нулевую матрицу. Символом I будем обозначать единичную матрицу. Для x, y ∈ R k (x, y) = x y скалярное произведение векторов, x = (x, x)

1 2

евклидова норма в конечномерном пространстве, B(x, r) = {x ∈ R n : x − x r} шар радиуса r > 0 с центром в точке x. Для вещественной прямоугольной k × m матрицы A через A k×m обозначаем норму матрицы, подчиненную евклидовым нормам векторов. Для S ⊂ R n символом ∂S обозначается граница S, ∂g ∂x (x) матрица Якоби отображе-

ния g(x), ∇ = ∂ ∂x1 • • • ∂ пространства суммируемых, суммируемых с квадратом и непрерывных вектор-функций на [t 0 , t 1 ]. Нормы в этих пространствах будем обозначать символами • L1 , • L2 , • C . Мы будем рассматривать управляемые системы с ограничениями вида ẋ(t) = f 1 (t, x(t)) + f 2 (t, x(t))u(t), t 0 t t 1 , x(t 0 ) = x 0 ,(1)где x ∈ R n вектор состояния, u ∈ R r управляющий параметр, f 1 : R n+1 → R n , f 2 : R n+1 → R n×r непрерывные отображения; управление u(t) вектор-функция из L 2 [t 0 , t 1 ]; t 0 , t 1 и начальное состояние x 0 фиксированы. Отображения f 1 : [t 0 , t 1 ]×R n → R n , f 2 : [t 0 , t 1 ]×R n → R n×r непрерывны, непрерывно дифференцируемы по x и удовлетворяют условиям: f 1 (t, x) l 1 (t)(1+ x ),(2)f 2 (t, x) n×r l 2 (t),(3)где x ∈ R n , l 1 (•) ∈ L 1 , l 2 (•) ∈ L 2 , t ∈ [t 0 , t 1 ]. Решением (траекторией) системы (1), отвечающим управлению u(•) ∈ L 2 , будем называть абсолютно непрерывную функцию x : [t 0 , t 1 ] → R n , для которой равенство (1) выполняется для почти всех t ∈ [t 0 , t 1 ]. При условиях (2), (3), для любых u(•) ∈ L 2 существует единственное решение x(t).

Определим интегральный функционал J(u(•)) равенством

J(u) = t1 t0 Q(t, x(t)) + u (t)R(t)u(t) dt,(4)

где

Q(t, x) неотрицательная непрерывная функция, а R(t) положительно определенная матрица, непрерывно зависящая от t. Через U обозначим непустое множество тех u(•) ∈ L 2 , для которых J(u(•)) µ 2 , где µ > 0 заданное число. Определение 1 Множеством достижимости G(t 1 ) системы (1) в момент времени t 1 будем называть совокупность всех концов траекторий x(t 1 ) в R n , отвечающих управлениям из U . Определим линейное отображение P k : R n → R k (k < n) следующим образом: y = P k x,

где матрица P k = I k×k O k×(n−k) . Таким образом, y состоит из первых k координат вектора x.

Матрицу P k , которая осуществляет переход от x ∈ R n к y ∈ R k , назовем матрицей проектирования на подпространство первых k координат.

Определение 2 Проекцией множества достижимости G 1k (t 1 ) системы (1) на подпространство первых k координат в момент времени t 1 будем называть совокупность векторов y 1 в R k вида y 1 = P k x(t 1 ), где x(t) траектории системы 1, отвечающие управлениям из U . Если G(t 1 ) множество достижимости, то G 1k (t 1 ) = P k G(t 1 ).

Рассмотрим задачу оптимального управления для системы (1).

Задача 1 J(u) → min, u(•) ∈ L 2 , x(t 0 ) = x 0 , P k x(t 1 ) = y 1 ,(5)где y 1 -заданный вектор из R k . Управление u(•) ∈ L 2 , удовлетворяющее ограничению задачи 1, назовем допустимым управлением.

Подчеркнем, что единственным ограничением в задаче 1 выступает терминальное ограничение на траекторию P k x(t 1 ) = y 1 . Схема доказательства приводимых далее результатов во многом повторяет схему из работы [5]. При необходимости мы будем делать соответствующие ссылки.

Определение 3 Управление u(•) ∈ L 2 , удовлетворяющее ограничению задачи 1, доставляет локальный минимум функционалу J(u), если существует ε > 0 такое, что для любого другого v(•) ∈ L 2 , удовлетво- ряющего ограничению задачи 1, такого, что u(•) − v(•) L2 < ε, имеет место неравенство J(u) J(v).

Вспомогательные результаты

Рассмотрим линейную нестационарную систему

ẋ(t) = A(t)x(t) + B(t)u(t), t ∈ [t 0 , t 1 ], x(t 0 ) = x 0 ,(6)где x ∈ R n , u(t) ∈ R r , A(t), B(t) -непрерывные на [t 0 , t 1 ] матричные функции. Определим симметричную матрицу V (t 0 , t) равенствами V (t 0 , t) = P k W (t 0 , t)P k , W (t 0 , t) = t t0 X(t, τ )B(τ )B (τ )X (t, τ )dτ, здесь W (t 0 , t) грамиан управляемости, X(t, τ ) = Φ(t)Φ −1 (τ ) матрица Коши однородной системы, где Φ(t) фундаментальная матрица, удовлетворяющая уравнению Φ(t) = A(t)Φ(t). Условие P k x(t 1 ) = y 1 ∈ R k означает, что траектория x(t) должна быть приведена в момент времени t = t 1 на k-мерную плоскость n-мерного пространства. Если для любой точки y 1 ∈ R k существует управление u(•) ∈ L 2 такое, что y 1 = P k x(t 1 ), то система (6) называется управляемой по выходу y на отрезке [t 0 , t 1 ]. Система (6) вполне управляема по выходу y = P k x на [t 0 , t 1 ] в том и только том случае, когда V (t 0 , t 1 ) положительно определена (см., например, [6, с. 204]). Очевидно, что если u(•) ∈ U , то t1 t0 u (t)R(t)u(t)dt µ 2 .

Так как R(t) положительно определена и непрерывна, то существует δ > 0 такое, что (R(t) − δI) 0 и, следовательно,

t1 t0 u (t)u(t)dt µ 2 δ .

Таким образом, все управления из U принадлежат гильбертову шару

D = u(•) ∈ L 2 : t1 t0 u(t) 2 dt µ 2 δ . Отсюда и из [5, утв. 2], следует компактность множества траекторий системы (1), удовлетворяющих огра- ничению U . Зададим последовательность управлений v m (•), которая отображается в последовательность траекторий x m (•). Линеаризованная вдоль (x m (•), v m (•)) система (1) задается парой матриц (A m (t), B m (t)), где A m (t) = ∂f 1 ∂x (t, x m (t)) + ∂ ∂x [f 2 (t, x m (t)) • v m (t)], B m (t) = f 2 (t, x m (t)). Справедлива Лемма 1 Если u m (•) → u(•) в L 2 и пара (A(t), B(t)), отвечающая управлению u(•), управляема по выходу y, то, начиная с некоторого m, пара (A m (t), B m (t)) будет управляемой на подпространстве R k (по выходу y).

Доказательство. Запишем соответствующие матрицы управляемости

V (t 0 , t 1 ) = t1 t0 P k M (t)B(t)B(t) M (t) P k dt, V m (t 0 , t 1 ) = t1 t0 P k M m (t)B m (t)B m (t) M m (t) P k dt, где Ṁ (t) = −A (t)M (t), Ṁm (t) = −A m (t)M m (t),(7)M (t 1 ) = M m (t 1 ) = I. Так как A m (•) → A(•) в L 1 [5, док-во леммы 1], то существует суммируемая функция γ(t) такая, что A m (t) n×n γ(t), ∀m, ∀t ∈ [t 0 , t 1 ]. Покажем ограниченность матрицы M m (t), повторяя схему доказа- тельства теоремы 3 [7, с. 9]. Существует K > 0 такое, что M m (t) n×n K ∀m, ∀t ∈ [t 0 , t 1 ]. Перейдем от (7) к интегральному уравнению и построим последовательные приближения: M m0 (t) ≡ I, M mp+1 (t) = I − t t1 A m (s)M mp (s)ds, p = 0, 1, 2, ... Из условия A m (t) n×n γ(t) и того факта, что A m (t)Z n×n A m (t) n×n Z n×n γ(t) Z n×n для любой матрицы Z, получаем существование и непрерывность всех приближений, где M mp+1 (t) − M mp (t) n×n t t1 γ(s) M mp (s) − M mp−1 (s) n×n ds, p = 1, 2, ... Положим ξ = max m,t M m1 (t) − M m0 (t) n×n на отрезке [t 0 , t 1 ], ψ(t) = t t1 γ(s)ds. По индукции доказывается, что на этом отрезке M mp+1 (t) − M mp (t) n×n ξ (ψ(t)) p p! , p = 0,M m0 (t) n×n −ξ • C M mp+1 (t) n×n M m0 (t) n×n +ξ • C для некоторого C > 0, получаем ограниченность M m (t) по норме. Поскольку d dt (M (t) − M m (t)) = −A (t)(M (t) − M m (t)) + (A m (t) − A(t)) M m (t),

то получаем следующую формулу

M (t) − M m (t) = t t1 Y (t, τ )(A m (τ ) − A(τ )) M m (τ )dτ. Здесь Y (t, τ ) есть фундаментальная матрица системы ẋ = −A (t)x, непрерывная по t, τ . Следовательно, последнее равенство означает, что M m (t) ⇒ M (t). Матрица B m (•) равномерно сходится к B(•) [5, док- во леммы 1]. Отсюда V m (t 0 , t 1 ) → V (t 0 , t 1 ), m → ∞. Если V (t 0 , t 1 ) положительно определена, то для достаточно больших m матрица V m (t 0 , t 1 ) будет также положительно-определенной, следовательно, пара (A m (t), B m (t)) управляема на подпространстве R k .

3 Экстремальные свойства граничных точек множества достижимости

Определим отображения F : L 2 → R n и S : R n → R k следующим образом: F (u(•)) = x(t 1 ), Sx(t 1 ) = P k x(t 1 ), где x(t) траектория системы (1), отвечающая u(•). Тогда их композиция H : L 2 → R k есть не что иное, как H(u(•)) = P k F (u(•)) = P k x(t 1 ). Cправедлива Лемма 2 Пусть функции f 1 (t, x), f 2 (t, x) непрерывны, непрерывно дифференцируемы по x и удовлетво- ряют условиям (2) и (3). Тогда функция H непрерывно дифференцируема по Фреше для ∀u(•) ∈ L 2 [t 0 , t 1 ], ее производная Фреше H : L 2 → R k определена равенством H (u(•))δu(•) = P k δx(t 1 ).(8)

Здесь δx(t) решение линеаризованной вдоль (u(t), x(t)) системы (1), отвечающее управлению δu(t) и нулевому начальному условию. Если линеаризованная система вполне управляема по первым k координатам, то Im H (u(•)) = R k .

Доказательство. В [5, л. 2] показана дифференцируемость отображения F по Фреше, производная которого F (u(•))δu(•) = δx(t 1 ). Линейное отображение S также дифференцируемо и S (x(t 1 ))δx(t 1 ) = P k . Таким образом, отображение H дифференцируемо по Фреше как композиция дифференцируемых отображений и имеет место равенство

H (u(•))δu(•) = S (x(t 1 )) • F (u(•)) = P k δx(t 1 ). Лемма 3 Функционал J(u(•)) непрерывен в любой точке u 0 (•) ∈ L 2 .

Доказательство. Как уже говорилось, под решением (траекторией) системы ( 1), отвечающим управле-

нию u(•) ∈ L 2 , понимаем абсолютно непрерывную функцию x : [t 0 , t 1 ] → R n , для которой равенство (1) выполняется для почти всех t ∈ [t 0 , t 1 ]. Зададим управление u 0 (•) из L 2 , которому соответствует траектория x 0 (•). Тогда cправедлива оценка для любых u(•) ∈ L 2 : J(u(•)) − J(u 0 (•)) = = t1 t0 [Q(t, x(t)) − Q(t, x 0 (t))] + (u (t) − u 0 (t))R(t)u(t) + u 0 (t)R(t)(u(t) − u 0 (t)) dt max t Q(t, x(t)) − Q(t, x 0 (t)) •|t 1 − t 0 |+ u(•) − u 0 (•) L2 R(•) C u(•) L2 + u 0 (•) L2 max t Q(t, x(t)) − Q(t, x 0 (t)) •|t 1 − t 0 | + C• u(•) − u 0 (•) L2 , C > 0.(9)

Функция Q(t, x(t)) непрерывна на компактном множестве траекторий системы (1), поэтому равномерно непрерывна на нём. Пусть ∆u

(t) = u(t) − u 0 (t), ∆x(t) = x(t) − x 0 (t), а ∆J(u(•)) = J(u(•)) − J(u 0 (•)). Из оценки ∆x(•) C = O ( ∆u(•) L2 ), которая была получена в доказательстве леммы 2 из [5], вытекает: для любого δ > 0 из ∆u(•) L2 < δ следует ∆x(t) < kδ, k > 0. Тогда, учитывая равномерную непрерывность Q(t, x(t)), всегда найдется такое δ = δ(ε) > 0, что для любого ε > 0 из ∆u(•) L2 < δ будет следовать |∆J(u(•))| < ε. Таким образом, выполнено условие непрерывности исходного функционала. Теорема 1 Пусть: 1) y 1 ∈ ∂G 1k (t 1 ); 2) u(•) ∈ U управление, переводящее систему (1) из x(t 0 ) = x 0 в точку x(t 1 ) такую, что P k x(t 1 ) = y 1 , x(t) соответствующая траектория; 3) Линеаризованная вдоль (x(t), u(t)) система (1) управляема по y = P k x на [t 0 , t 1 ]. Тогда управление u(•) ∈ L 2 доставляет локальный минимум функционалу J(u(•)) при ограничениях x(t 0 ) = x 0 , P k x(t 1 ) = y 1 и величина J(u(•)) = µ 2 .

Доказательство. Проведем доказательство от противного. Пусть найдется u(•), переводящее систему из x(t 0 ) = x 0 в P k x(t 1 ) = y(t 1 ) ∈ ∂G 1k (t 1 ), которое не является локально оптимальным в задаче 1, иначе говоря, для любого p существует допустимое управление u p (•) такое, что

u(•) − u p (•) L2 < 1 p , J(u p (•)) < J(u(•)), либо J(u(•)) < µ 2 . Тогда существует последовательность u p (•) → u(•), p → ∞ в L 2 такая, что P k x p (t 1 , u p (•)) = y(t 1 ), а J(u p (•)) < µ 2 .

Далее, выберем p настолько большим, чтобы пара (A p (t), B p (t)) была управляемой по y (см. лемму 1) и зафиксируем это p. Обозначим

δ = µ 2 − J(u p (•)) > 0. Тогда J(u p (•)) < µ 2 − δ/2.

Из леммы 2 получаем, что линеаризованная вдоль

(u p (•), x p (•)) система управляема по y, это эквивалентно условию Im H (u p (•)) = R k . Тогда по теореме Грейвса [8, с. 105] для некоторого m > 0 и всех достаточно малых r выполняется включение B(y 1 , mr) ⊂ F (B(u p (•), r)). Таким образом, B(u p (•), r) ⊂ U ⇒ F (B(u p (•), r)) ⊂ F (U ) ⊂ G 1k (t 1 ), отсюда следует B(y 1 , mr) ⊂ G 1k (t 1 ), что противоречит условию y 1 ∈ ∂G 1k (t 1 ).

Выпишем необходимые условия оптимальности в форме принципа максимума для задачи 1, предполагая, что Q(t, x) непрерывно дифференцируема.

Функция Понтрягина здесь имеет вид

H(p, t, x, u) = −p 0 Q(t, x) + u R(t)u + p (f 1 (t, x) + f 2 (t, x)u), p 0 0.

Локально оптимальное управление удовлетворяет принципу максимума (см., например, [

): существуют (p 0 , p(•)) = 0 такие, что H(p(t), t, x(t), u(t)) = max v∈R r H(p9, §5, т. 2](t) = − ∂ ∂x H(p(t), t, x(t), u(t)) = −A (t)p(t) + p 0 ∇ Q(t, x(t)), p j (t 1 ) = 0, j = k + 1, ..., n.(t), t, x(t), v). ṗ

Равенства (10), в силу частичного закрепления концов траекторий, есть условия трансверсальности на правом конце. Через (A(t), B(t)), как и ранее, обозначаем матрицы линеаризованной вдоль (x(t), u(t)) системы. Если эта линеаризованная на (x(t), u(t)) система вполне управляема по выходу y, то p 0 = 0.

Предположим, что p 0 = 0. Тогда p(•) = 0. Максимум по v линейной по управлению функции

H(t, x(t), p(t), v) = p (t)f 1 (t, x(t)) + p (t)f 2 (t, x(t))v = p (t)f 1 (t, x(t)) + p (t)B(t)v существует, если p (t)B(t) ≡ 0, где p(t) удовлетворяет системе ṗ(t) = −A (t)p(t). С учетом (10) p(t 1 ) имеет вид p(t 1 ) = l k 0 , где l k ∈ R k , l k = 0. Представим p(t) в виде p(t) = l k , 0 X(t 1 , t), где X(τ, t) - фундаментальная матрица системы ẋ = A(t)x. Тогда l k , 0 X(t 1 , t)B(t) ≡ 0, или l k P k X(t 1 , t)B(t) ≡ 0.

Поэтому

l k P k x(t 1 ) = l k P k t1 t0 X(t 1 , t)B(t)u(t)dt = 0 ∀u(•) ∈ L 2 .

Иначе говоря l k ⊥ P k x(t 1 ). Из последнего и P k x(t 1 ) = R k (следует из управляемости на R k ) получаем l k = 0, что противоречит начальному предположению p(t 1 ) = l k 0 = 0. Cправедливо принять p 0 = 1 2 . Тогда из принципа максимума получим u(t) = R −1 (t)f 2 (t, x(t))p(t).

Замыкая исходную систему данным управлением, имеем

ẋ(t) = f 1 (t, x(t)) + f 2 (t, x(t))R −1 (t)f 2 (t, x(t))p(t), x(t 0 ) = x 0 , ṗ(t) = − ∂f1 ∂x (t, x(t)) + D (t, x(t), v) p(t) + 1 2 ∇ Q(t, x(t)), p j (t 1 ) = 0, j = k + 1, ..., n. (11) где обозначено D(t, x, v) = ∂ ∂x (f 2 (t, x)v), v = R −1 (t)f 2 (t, x(t))p(t).

4 Численное моделирование

Описание алгоритма

Соотношение (11) можно положить в основу следующего алгоритма построения границы множества достижимости.

Выбирая p(t 0 ) = 0 и интегрируя систему (11) (учитывая условия трансверсальности), мы получим управление и траекторию, удовлетворяющие принципу максимума. Перебирая p(t 0 ) из регулярной сетки, аппроксимирующей область {p(t 0 ) ∈ R n : |p i (t 0 )| c i , i = 1, 2, ..., n}, интегрируя систему (11) и отбирая те траектории, для которых |J(u(•)) − µ 2 | не превосходит малого ε 1 > 0, мы получим аппроксимацию части границы множества достижимости, образованную точками y 1 = P k x(t 1 ). Плюс ко всему, при отсеивании траекторий, нужно проверять выполнение условий трансверсальности на правом конце: |p j (t 1 )| < ε 2 , j = k + 1, ..., n, для малого ε 2 > 0. Для достаточно больших c i и малых ε 1 и ε 2 аппроксимироваться будет вся граница, если на каждой из возможных траекторий выполняется условие полной управляемости линеаризованной системы. Данный алгоритм, использованный в [5] при построении множеств достижимости, требует достаточно больших вычислительных затрат и априорной оценки области начальных данных для p(t 0 ).

Здесь мы используем другой алгоритм, обеспечивающий более экономную схему вычислений. Приведем его описание для двумерного случая и P k = I. Будем выбирать направления на единичной сфере: {p * ∈ R 2 : p * = cos α sin α , α ∈ [0; 2π]}. На каждом луче {p ∈ R 2 : p = kp * , k > 0} мы ищем значение k, для которого достигается минимум J(u(•)) − µ 2 2 , на полуинтервале [0; +∞). Для этой цели можно использовать стандартную процедуру одномерной оптимизации. Найденное значение k определяет наш вектор p = k p * , выбирая который в качестве начального условия для системы принципа максимума, мы получим точку x(t 1 ), лежащую на границе множества достижимости. В трехмерном случае рассуждения похожи, если использовать, например, сферические координаты.

Пример 1

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

ẋ1 = ax 1 − cx 1 x 2 , a > 0, c > 0, ẋ2 = −bx 1 + dx 1 x 2 , b > 0, d > 0,(12)

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

Обозначим объемы закупок электроэнергии в секторе свободной торговли через x 1 , цену электроэнергии через x 2 . Спрос тем быстрее уменьшается, чем больше проданной электроэнергии и чем больше ее цена. Иными словами, цена тем быстрее снижает закупки, чем более насыщен рынок электроэнергии по этой цене (x 1 x 2 ). Поэтому, если объемы закупок электроэнергии ненулевые, то объемы закупок электроэнергии меняются по закону ẋ1 = ax 1 −cx 1 x 2 (a > 0, c > 0). С другой стороны, прибыль, получаемая поставщиками от продажи электроэнергии, стимулирует увеличение цены. Прибыль пропорциональна количеству проданного товара по его цене (x 1 x 2 ). Поэтому имеем ẋ2 = −bx 2 + dx 1 x 2 (b > 0, d > 0). Модель "цена-объемы покупок" построена.

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

Обозначим за управление u неопределенные факторы, влияющие на объем покупок и наложим на x 1 (объем закупок) и u квадратичные интегральные ограничения: ẋ1 = ax 1 − cx 1 x 2 + u, a > 0, c > 0, ẋ2 = −bx 1 + dx 1 x 2 , b > 0, d > 0, Результаты численного моделирования представлены на рис. 3 и 4.

Благодарности

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