Оптимальное управление движением жидкости со свободной поверхностью Текст научной статьи по специальности «Математика»

УДК 517.977, 519.626

Оптимальное управление движением жидкости со свободной поверхностью

1 2 © А. А. Гурченков , А. М. Романенков

1 МГТУ им. Н.Э. Баумана, Москва 105005, Россия.

2 МАТИ РГТУ им. К.Э. Циолковского Москва 109387, Россия.

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

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

Найдем потенциал <р(х, у, 1, поля скоростей идеальной несжимаемой жидкости со свободной поверхностью / (х, у, I) в ограниченном трехмерном пространстве, удовлетворяющий уравнению Лапласа в области о : {(х, у, 1) : 0 < х, у < 1, 0 < 1 < / (х, у, ¿)}, с граничными

условиями Неймана на границах области и двумя нелинейными условиями на свободной поверхности [1]. Отметим, что область, в которой справедливо уравнение Лапласа, не является стационарной (решение не существует) [2].

Предположим, что колебания жидкости слабо возмущенные и свободная поверхность мало отклоняется от равновесной, т. е. граничное условие можно перенести со свободной поверхности на равновесную (1 = 1) [3]. Тогда

А< = 0;

д< д<

дх х=0 дх х

д< д<

ду у=0 ду у

д< =

~д.Г 1=0

= 0;

= 0;

(1)

На свободной поверхности выполняются два нелинейных условия: чисто кинематическое

df ^д/ ду df ду _ ду dt дх дх ду ду dz '

(2)

отражающее тот факт, что частица жидкости, попав на свободную поверхность, навсегда остается на ней, и динамическое

дУ + 1 дt 2

ду дх

2 ^ду^2 Аду V дУ J V &

+ g (f -1)_ 0,

(3)

являющиеся следствием уравнения движения Эйлера.

Применив метод Фурье, можно представить потенциал поля скоростей двойным рядом Фурье:

у(х, у, z, t) _ ^ ^Pnk (t)ch (л^п2 + к2 z) cos;rnx cos яку. (4)

n _0 к _0

Функцию, выражающую форму свободной поверхности, также разложим в ряд Фурье по косинусам:

f (х, у, t) _ 1 + ^ qnk (t) cos япх cos яку.

(5)

n _ 0 к _0

Будем воздействовать на сосуд с колеблющейся жидкостью, ограниченной переменной силой, которую назовем управлением [5-7]:

и {1) = {их (г), иу (г), иг (г)},

где их (г), иу (г), и2 (г) — компоненты силы по соответствующим осям координат, причем иШщ < их (г), иу (г), uz (г) < итах.

Введенное воздействие приведет к изменению уравнения (3):

^ + — + Ф' = 0, дг 2

где Ф' = Ф'(г) = ^(0, 0, + и(г), 5(г)) — потенциальная энергия; 5(г) = {х, у, /(х, у, г)-1}; ,•) — стандартное скалярное произведение в М3 [4].

Окончательно получим динамическое условие с учетом воздействия:

дю 1 — + -

dt 2

где для краткости принять

+ -S2 + ux (t) x + Uy (t) y + (g + uz (t)) (f -1) = 0, (6)

dt 2

& = id^f+f дЮ?+fd£

L dx J L dy J L dz

Будем рассматривать задачу (1), (2), (6). Очевидно, что изменяя ux (t), uy (t), uz (t), получаем различные решения задачи. По теореме

о непрерывной зависимости системы обыкновенных дифференциальных уравнений от начальных условий и правой части решение непрерывно зависит от u (t) .

Введем функционал, характеризующий интенсивность колебаний [8, 9]:

<Х)

J (u (t))= z { [t, u (t)] + fa [t, u (t)]}. (7)

n 2 + k 2 >1

Необходимо подобрать u (t) таким образом, чтобы функционал (6) был минимальным в момент времени T, т. е.

J (T ) = J [u (T min.

Утверждение 1. Управления не входят в кинематическое условие (2) (оно остается неизменным), а входят линейно в уравнения для определения коэффициентов потенциала поля скоростей.

Доказательство. Первая часть утверждения очевидна. Докажем вторую часть утверждения.

Подставим в условие (6) выражения для потенциала поля скоростей и формы свободной поверхности и выполним стандартные действия для получения уравнения для РПк (t): умножим на cos nnx cos nky и проинтегрируем по квадрату [0, 1] х [0, 1].

Обозначим множество интегрирования как 12 = [0, 1] х [0, 1] и вычислим интеграл, в который входит управление:

jj{ux (t) x + uy (t) y + [g + uz (t)] (f -1)} cos nnx cos nky dxdy =

jj ux (t) x cos nnx cos nkydxdy +

12

jj uy (t) y cos nnx cos nky dxdy +

12

jj[ g + uz (t)] (f - 1)cos nnx cos nky dxdy.

12

+

I 2

+

12

Найдем третий интеграл:

jj[g + uz (t)]( f -1) cos япх cos яку йхйу

[g + uz (t)]jj

^ ^ Чпк (t) cos яп 'x cos як 'у

n'_0 к'_0

cos япх cos яку йхйу _

си си

[ g + uz (t)] ^ ^ qпк (t) jj cos яп 'х cos як 'у cos япх cos яку йхйу.

n ' _ 0 к '_0

Двойной интеграл равен некоторой константе C (п, к), которая зависит от п и к, т. е.

jj [g + uz (t)] (f -1) cos япх cos яку йхйу _ [g + uz (t)] qn (t) C (n, к).

Теперь вычислим первый интеграл:

jj ux (t) х cos япх cos якуйхйу _ ux (t) j cos якуйу jx cos япх dx.

Очевидно, что

jcosяkydy _ j i

0, к * 0;

к _ 0;

i 1 i 1 \ 1 i

jx cos япх dx _ — jx d sin япх _ — х sin япх|0 - j sin япхйх

{ яп{ яп '

(яп )

|1 cos яп -1 2 cos япх|0 _-

(яп )

0, n _ 2к;

2

(яп )

n _ 2к -1;

п _ 0,

где к е N.

i

i

i

i

i

Второй интеграл вычисляется аналогично. Окончательно имеем

JJ {их (t) х + uy (t) y + [g + uz (t)] ( f - 1)J cos nnx cos nky dxdy -

i2

+ [g + uZ (t)]qok (t)C(0,k), n - 0, k * 0(mod2);

^ + [g + uz (t)] qno (t) C (n, 0), n * 0 (mod 2), k - 0;

(8)

[g + uz (t)]qnk (t)C(n, k),

(nk * 0) или (n - 0, k = 0(mod2)) или k - 0, n = 0(mod2).

Утверждение 2. Симметричные колебания невозможно погасить, применяя лишь воздействия по осям х и y.

Доказательство. Управления ux (t) и uy (t) входят лишь в уравнения для определения коэффициентов потенциала поля скоростей с нечетными индексами, а в уравнения с четными индексами входит только управление uz (t) (см. формулу (8)). Из формул (4) и (5) следует, что коэффициенты с нечетными индексами представляют собой амплитуды асимметричных возмущений свободной поверхности и потенциала скоростей, т. е. утверждение очевидно.

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

Получим систему дифференциальных уравнений, ограничившись лишь несколькими членами ряда:

ж2 ж2

q01 = жshжP01 - — q10 P^Wl^- — q11Pwchж-

-Ж q0p02ch^ - ж q02 Poich^

q[0 - жshжpo - ж2q10P20ch2ж - ж2q20pochж -

2 2 - ж^ q0ipi icW^ - ж^ q„P0!ch^

q[1 - ж^Í2sh^f2жP11 - ж2^Р20Л2ж - ж2q20P^W^ -

-ж2 q11P02ch2ж - ж^02 P^W^; ж2

q02- 2жsh2жPo2 +—qolpolchж;

(9)

ж2

q2o = + — qwPwchЖ

р + и ж2 2

Р = "^Г""2^01 "-¡^0^1 + РиРи-ПЗж) 2^; епж -пж ж

р + и ж^ 2

Р0 = -^г-^ю --г-^оРго-ЬЗж + Р01Р11С1) - — их; (9) -пж -пж ж

+ 2 Р'1 = - - d—2-(P"lPloSh2ж + 2 ^)■

2

20 -Ь2ж 4-Т2ж^ 11 10'

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

Отметим, что важен лишь факт возникновения того или иного колебания и не важна его амплитуда (все коэффициенты равны 1).

Утверждение 4. Если изначально возбуждены асимметричные колебания, то с течением времени будут возбуждены и симметричные колебания. Более того, если в начальный момент времени возмущена любая ассиметричная гармоника, то с течением времени будут возмущены другие ассиметричные, а также симметричные гармоники.

Утверждение 5. Если изначально возбуждены только симметричные гармоники, то асимметричные гармоники не будут возбуждены и если в начальный момент времени возбуждена любая симметричная гармоника, то с течением времени все остальные симметричные гармоники будут возбуждены, а ассиметричные гармоники возбуждены не будут.

Сформулируем задачу оптимального управления. Пусть определен функционал

3\и рт)] = Р2 [Г, и рт)] + Р0 [г, и рт)] + р [г, и рт)] + +Р22 \Т, и рт)] + Р0 \Т, и рт)] + q00l \Т, и рт)] + ^ \Т, и рт)] + + ^ \т, и рт)] + ^ [т, и рт)] + q22o \т, и рт)].

= р,

Необходимо найти вектор-функцию йор1 (г), при которой

3 (( (Т )) = т т 3 \и (Т)],

где и (Т) = {их (Т), иу (Т), иг (Т)}: их (г), иу (г), иг (г) е ¿2 ([0; Т]); итп < < их (г),иу (г),иг (г) < итах, удовлетворяющую системе дифференциальных уравнений (9) и начальным условиям Коши

ГЧ01 (0 Ч10 (0) Чп (0)

Р0 (0)

Р20 (0)

где Р е м10.

Данная задача является задачей оптимального управления со свободным концом и с фиксированным временем [4]. Ее решение будем искать, основываясь на принципе максимума Понтрягина [10, 11].

Составим, согласно этому принципу, функцию Гамильтона для нашей задачи:

н ((Р, Ч), V, и (г)) =

Г я2 я1

= I я8ЬяР01 -—ч10РисЬ>/2я- —ч11Р0еЬя-

-я Ч01Р02сЬ2я-я Ч02Р01сЬя 1^1 +

+1 явЬяР0 - я2ч10Р20сЬ2я - я2ч20Р10сЬя -

(10)

(

- — Ч01Р11С^>/2я-^ ЧпР0Ая^2 +

-я ч11Р20сЬ2я-я ч20

- я2ч11Р02сЬ2я - я2ч02 РисЬ>/2я) V +

> Г

^4 +

я

2я8Ь2яР02 + —ч01Р01сЬя

я

Л

2явЬ2яР20 н—Ч10Р0сЬя

V5 +

__._

-Пж 01 -Пж 12

—(-ВДС + PmPm-h3ж\- —^

р+иq -ж_ГРР -П3ж+1РРС иАи

, ql0 , I р0Р20-П3^ + Р01р1С1 I 2 их -Пж -Пж I 2 у ж

11 + 11 +

Г

Р + uz

Г р + Uz ж2

--qm--

-П2ж 4-П2ж

ж

рР01Рю8П2ж + 2РцР20С2)) + (10)

v Г

рр1§П2^^2 -Р021 )]^9 +

ж

Р + Uz

-П2ж ^20 4-П2ж

рр8П2ж72 - р0)

1

где 1 = 1 Pt), (/ = 1.. .10) — сопряженные функции.

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

, дН 2 л 1 /-> 1 =--= ж Р02-П1жу/1 +

01

ж „ , ГГ ж^, Р + и

+—р1-П^2ж12 р01-Пж^4+—— ц;

2 2 -Пж

11 =

дН

= ж2 Р20-П2жц +

10

ж ^ . пт ж _ , р + и +—Рп-Ы2ж^1 - — РР0-Пж^5 + ——^ц; 2 2 -Пж

дН ж2-Пж / _ г, ч 1 = = —^—рр10^1 + Р01^2 ) +

11

(11)

+п2 рр20+Р02) -П2ж 12 + -Пж и^ 1»;

14 =■

15=-

дН

02

дН

^20

= ж2 р Р01-Пж1 + Р11-П\[2жу/3 ) ■ ж2 рР10-Пжу/2 + Р11-Пу[2жщ3 ) +

р+и

-П2ж

р+и

-П2ж

1;

1;

V = ■

дН

дР0

- (Ья - я2ч02сЬя) V

01

я2 я2

+ — Ч11сЬя2 Ч01сЬя V4 -

я я Г 1

-Р02сЬЗя6 Р1^7С1 + сья2 [ ^в^- 2Р0^9

(11)

V? = ■

дН

дР

- (Ья - я2д02сЬя) V2

10

я2 я2

+ "2 ЧпсЬя^ - ^ ЧюсЬя -

20 сЬЗя я2 п ^ х л , п и2

-я Р2 "' "'

20 сЬя V7 2сЬя

Начальными данными для системы (11) является следующий вектор:

Г-2Ч01 (Т -2Ч10 (Т )

-2Ч11 (Т)

V(T )= -2Р01 (Т) -2Р10 (Т) -2 Р1 (Т )

-2Р20 (Т)

где v(T)- вектор из М10.

Чтобы определить начальные данные для системы (11), необходимо решить систему (9). Решение этой системы найдем методом Рунге — Кутты четвертого порядка. Назовем этот процесс прямым проходом. При этом необходимо хранить данные со всех слоев, так как на их основе будет строиться оптимальное управление. Для решения системы (11) также используем метод Рунге — Кутты, но при этом знак шага изменим на минус. Этот процесс назовем обратным

проходом. При обратном проходе необходимо выбирать новое управление. На каждой итерации управление следует выбирать так, чтобы функция Гамильтона была максимальной по управлению, т. е.

m ax H = - min

U(t) umin , uy, u

-(y/7ux + ¥Uy )

f

¥б001 + ¥7^0 , w%4ll , ¥9^2 + ¥10020

Л

chn

ch2^

u, + D

где .О — слагаемые функции Гамильтона, которые не зависят от управлений.

В силу того, что управления линейно входят в функцию Гамильтона (9), поиск оптимального управления достаточно прост:

u = <

u„ =

uy =

u, = <

| Umin , если ¥7 > 0;

tUmax, еСли ¥7 ^ 0;

f Umin, еСли ¥б > 0;

tUmax, еСли ¥б ^ 0;

Umin, еСли

Umax , еСли

¥б001 + ¥7010. + ¥011 + ¥9002 + ¥щ020 1 > ().

ch^

¥б001 +¥7010

ch^

chW2

ch2^

, ¥8011 + ¥9002 +¥10020 |< 0

chW2

ch2^

(12)

Ниже приведены численные тесты и результаты вычислений для различных возмущений. На рис. 1 рассмотрены симметричные возмущения. В начальный момент времени д02 ^ 0. Численный эксперимент

показывает, что асимметричные гармоники не возбуждаются.

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

На рис. 3 показана асимметричная гармоника д01 ^ 0. Видно, что возникают и симметричные возмущения (д02 ^ 0).

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

^ Численное решение задачи о движении жидкости со свободной поверхностью (п=3) нив

Для расчета поверхности Для расчета поля скоростей Аквариум в разрезе Х=0 иг=01к=01И>

ч01= 0 Р01= 0 №

410= 0 Р10= 0

ч11= 0 Р11= 0

ф2= -0,197756994759908 Р02= 9,38961028497538Е-05

д20= 0 Р20= 0 1к

Поверхность 1= 0,80224300524009 Скорость *= 0.0251403774539665

д= 10 0,001 Значение и= 19,55

Отключить управления? 0 по осих [7 по оси у О по оси г

11так= 5 1Ып= -5 Цу

В точке [7] Считать Интервал таймера 1 | Принять

Х= 0 Фиксировать при расчете !§) X 0 У

0 I Итерация | | Обнулить ] [ Выход |

409

Рис. 1

Для расчета поверхности д01= О

я10» о ч11= о

задачи о движении жидкости со свободной

Дли расчета поля скоростей

Р01 = Р10= Р11 =

поверхностью (п-3)

Аквариум в разрезе Х=0

д02= 0.0114-312846535436 д20= 0

Р02= 0,0012219813294139 Р20=

Поверхность I = 1.01143128465354 Скорость у= 0.327181543544178

д= 10 Ь= 0,001 Отключить управлении? [V] п

В точке [7] Считать

Значение ^ 0,066

«ну □

11пип= -5 Интервал таймера

Фиксировать при расчете ® X © У ..Запустшь 3 [ Итерация | [ Обнулить""] | Выход

673 1,636; 27,59; 6,210; 65,32;

Ив

aj Численное решение задачи о движении жидкости со свободной гтовернностьго (п=3} ЕЗЖв

Для расчета поверхности Дгр расчета поля скоростей Атвариум в разрезе Uz=BLk=0Uy=D

q01= 0,29632178176635 Р01- DDCÍ44691921S90529 Uz

q10= В Р10- В

q11- В Р11= В

qB2= -D.BB13732DS3E7BB41 РВ2= 4.57575554D5D753E-B7

q2B= В Р20- 0 Ux

Поверхность f- Скорость V-

g= 10 h- 0,001 Значение J- 4410

Отключить управления? П по ОСИ X ПО ОСИ Y И| по ОСИ Z

Umax- В,5 Lknin— -0,5 иу

Вто'ке ОСгтать Интервал таймера 1 ; Приять

Х- 0 Фиксировать при расчете ® X © Y

Y= [' [ Запустить ¡ [ Итерация | | Обнулить | Выход

35

Рис. 3

Di1 Чис/генное решение задачи о движение жидкости со свободной поверхностью í"=Tl

Для рамета поверхности Для расчета поля скоростей Аквариум в разрезе Х=В

401- -0.1581315В66Ш65

q10= -0.0609437S3B642337

Р01-Р10=

qll= -D.D17040BS61SS0767

-D.M134S4467521645 адШ14028Е142158

Pll= -D.00161S2S247B233S

qB2= -D.D12553B3723425B4

q2D- B,DB255424S5757B37B

PD2= D,DBB1D5S5SBB4S67346

P2D- <777ZM5B)t»a

Поверхность f= 0,75388497260871; Скорость v= -S.721752T41S1S13E-g= 10 h= 0 001 Значение J= 1559

Отклони управления? П по осих по оси у Цпоосит ÛTTÏ1-

-0.5

Umax- 0,5

В точке |V] Считать

х. ;

Y= 0 | Запустить | Итерация | Обнулить ] Вьиод

Интервал таймера |1 ; Принять Фиксировать при расчете ® X © Y

UT =-0.5 Lk-0.5L^-0.5

ЛИТЕРАТУРА

1. Черноусько Ф.Л. О движении тела с полостью, частично заполненной вязкой жидкостью. ПММ, 1966, т. 30, вып. 6, с. 476-494.

2. Красовский Н.Н. Теория управления движением. Москва, Наука, 1968, 434 с.

3. Гурченков А. А. Динамика завихренной жидкости в полости вращающегося тела. Москва, Физматлит, 2012, с. 221 с.

4. Крейн С.Г., Моисеев Н. Н. О колебаниях твердого тела, содержащего жидкость со свободной поверхностью. ПММ, 1957, т. 21, вып. 2, с. 97-114.

5. Гурченков А.А., Есенков А. С., Цурков В.И. Управление движением ротора с полостью, содержащей идеальную жидкость. Изв. РАН ТСУ, 2006, № 3, с. 82-89.

6. Гурченков А.А., Есенков А.С., Цурков В.И. Управление движением ротора с полостью, содержащей вязкую жидкость. Автоматика и телемеханика.

2007, № 2, с. 81-94.

7. Дьяченко В.П. О колебаниях гироскопа, частично заполненного жидкостью. ДАН УССР. Сер. А, 1971, № 10, с. 56-64.

8. Гурченков А.А., Носов М.В., Иванов И.М. Оптимальное управление движением волчка с жидким наполнением. Тез. докл. ХУП Всероссийской конференции «Теоретические основы и конструирование численных алгоритмов и решение задач математической физики». Абрау-Дюрсо. 15-21 сентября 2008 г. Новосибирск, НИИ гидродинамики им. М.А. Лаврентьева, СО РАН,

2008, с. 143-144.

9. Гурченков А.А., Есенков А.С., Цурков В.И. Управление движением ротора с полостью, содержащей идеальную жидкость. Изв. РАН ТСУ, 2006, № 1, с. 135-142.

10. Гурченков А. А., Носов М. В., Цурков В. И. Управление вращающимися твердыми телами с жидкостью. Москва, Физматлит, 2011, 202 с.

11. Gurchenkov А.А., Nosov М.У., Tsurkov V.I. Control of Fluid-Containing Rotating Rigid Bodies. CRC Press, 2013, 147 p.

Статья поступила в редакцию 21.02.2013

Ссылку на эту статью просим оформлять следующим образом:

А.А. Гурченков, А.М. Романенков. Оптимальное управление движением жидкости со свободной поверхностью. Инженерный журнал: наука и инновации, 2013, вып. 2. URL: http://engjournal.ru/catalog/appmath/hidden/613.html

Гурченков Анатолий Андреевич — профессор кафедры «Высшая математика», д-р. физ.-мат. наук; автор более 100 научных работ, 10 из которых — монографии; сфера научных интересов: управление вращательными твердыми телами с низким наполнением, устойчивость динамических систем с жидкостью. e-mail: challenge2005@mail.ru

Романенков Александр Михайлович — ассистент кафедры «Прикладная математика и информационные технологии»; сфера научных интересов: обыкновенные дифференциальные уравнения, уравнения в частных производных, оптимальное управление, управление нестационарными процессами. e-mail: romanaleks@gmail.ru