Сер. 10. 2012. Вып. 3
ВЕСТНИК САНКТ-ПЕТЕРБУРГСКОГО УНИВЕРСИТЕТА
ПРИКЛАДНАЯ МАТЕМАТИКА
УДК 517.977
А. В. Аргучинцев, В. П. Поплевко
ОПТИМИЗАЦИЯ ПРОЦЕССА РЕКТИФИКАЦИИ В КОЛОННЕ *)
1. Введение. Процесс ректификации широко применяется в нефтегазоперерабаты-вающей, химической, нефтехимической и других отраслях промышленности. Он предназначен для получения из многокомпонентной смеси продуктов с заданной концентрацией. Этот процесс описывается системами дифференциальных уравнений в частных производных [1]. Для него характерно большое количество параметров, связанных между собой сложными зависимостями. Значительная часть параметров является функциями временной и пространственной координат. Цель управления процессом ректификации - поддержание заданного состава целевого продукта. В данной работе в качестве управлений выбираются потоки готовых продуктов. Предполагается, что управления принадлежат классу гладких функций.
Задача оптимального управления процессом ректификации была сведена к оптимизационной задаче с управляемыми дифференциальными связями на границе. Для численного решения задачи был применен метод улучшения гладких управлений [2]. Расчеты проводились с использованием системы МАТЬАВ 7.0. Результаты численного эксперимента приведены в таблице (см. п. 5.1, 5.2).
2. Постановка задачи. Математическая модель процесса ректификации представляет собой систему уравнений, записанных относительно концентраций компонентов в жидкой и паровой фазах [1]:
Аргучинцев Александр Валерьевич — доктор физико-математических наук, профессор, заведующий кафедрой методов оптимизации Института математики, экономики и информатики Иркутского государственного университета. Количество опубликованных работ: 87. Научные направления: оптимальное управление, системы с распределенными параметрами. E-mail: [email protected].
Поплевко Василиса Павловна — кандидат физико-математических наук, доцент кафедры методов оптимизации Института математики, экономики и информатики Иркутского государственного университета. Количество опубликованных работ: 13. Научные направления: оптимальное управление, системы с распределенными параметрами. E-mail: [email protected].
*) Работа выполнена при частичной финансовой поддержке Российского фонда фундаментальных исследований (грант № 11-01-00713) и федеральной целевой программы «Научные и научно-педагогические кадры инновационной России» на 2009—2013 гг.
© А. В. Аргучинцев, В. П. Поплевко, 2012
d(GxXj) _ d(Lxj) dt ds
kvM - vD + Ф
(1)
d(GyVi) d(VVi) dt ds
+
kvM - Vi) + ^Vi >
N N
¿=1 ¿=1
где Ж - число компонентов в смеси; x¿,y¿ - концентрации г-го компонента в жидкой и паровой фазе соответственно; в - координата вдоль колонны; Ь - время работы колонны; Ох, Оу - количество жидкости и пара в колонне соответственно; Ь = Ь(Ь), V = V(Ь) -потоки жидкости и пара в колонне; кУ1 - коэффициент массопередачи; у* - равновесная концентрация компонента в паре; Фхн, ФУ€ - вводимый поток жидкости и пара в колонну соответственно. Концентрация компонентов, находящихся в испарителе, определяется из уравнения материального баланса
-"77- = -Я-- Уг(во^о) = Уго(«о)- (2)
аЬ Цу
Уравнения материального баланса в конденсаторе
ахн (з1,г) Ь(г) + Б(г)
dt Qx
(yi(si,t) - xi(si,t)), Xi(si,to) = Xio(si). (3)
В граничных условиях (2), (3) присутствуют управляемые потоки D(t), W(t). Здесь Qy, Qx - количество жидкости в испарителе и конденсаторе соответственно.
Решение начально-краевой задачи (1)-(3) понимается в классическом смысле. Цель задачи - достижение заданных параметров 6ц, 02i в конечный момент времени
N tl
J =Е К ii(xi (si,t) - вц)2 + K2i(yi(so,t) - e2if] dt ^ mm, (4)
i=1 to
где К1г, K2i - коэффициенты, определяющие ценность продукта.
Функции управления (потоки D(t), W(t)) удовлетворяют дополнительным ограничениям, задающим баланс потоков сырья и готового продукта в колонне в каждый момент времени:
C(t) = D(t)+W(t), t e T, (5)
или за весь период работы колонны:
j[C(t) - D(t) - W(t)] dt = 0, t e T, (6)
T
в которых C(t) - заданная функция. Предполагается, что W(t), D(t) принадлежат классу гладких функций.
Задача (1)-(6) рассматривается при следующих предположениях:
1) подвод сырья осуществляется только в жидкой фазе
Фуг (S,t) = 0, Фхг (S,t) = Xfi C(t^ (s),
где Xfi - концентрация компонента в жидкости потоков сырья; C(t) - вводимый поток сырья; фхн (s) - функция распределения потока по длине колонны - заданные функции;
2) в ректификационной колонне с увеличением потоков L и V возрастают удерживающие способности Qx и Qy. Будем считать эту зависимость прямо пропорциональной: L/Qx = c1 = const, V/Qy = c2 = const;
3) коэффициент массопередачи в паровой фазе имеет вид ky^ — ky = kV(s,t), k = const;
4) концентрация паровой компоненты в равновесной ситуации находится по формуле у*(в,1) = р(в,£)хг(в,£), где р(в,£) определяется из эксперимента. С учетом всех предположений система (1) преобразуется к виду
dXi^ _CldxiM =Ai{s,t)xi{s,t) + Bi{s,t)yi{s,t) + Fi{s,t),
+ C2MhH = A2(s,t)xi(s,t) + B2{s,t)yi{s,t) + F2{s,t), г = 1JV. dt ds
(7)
Здесь
, , С1С{1)фхМ - CikV{t)p{s,t) - L'{t)
Ai{s,t) =-—-, A2(s,t) = c2kp{s,t),
ci C(t)<f>xi(s)
' ' =-Ut)-' =
*,.,> = = +
Задача (2)—(7) принадлежит классу задач оптимизации гиперболических систем с управляемыми дифференциальными связями на границе, где управления выбираются из класса гладких функций и удовлетворяют поточечным (интегральным) ограничениям. В работе [1] данная задача рассматривалась в классе кусочно-непрерывных функций. Методика [2, 3] позволяет исследовать ее в классе гладких управляющих функций.
3. Необходимые условия оптимальности. Введем следующие функции: H(ф,x,y,s,t) =< ф^^), Ai(s,t)x(s,t) + Bi(s, t)y(s, t) + Fi(s,t) > + + < ф2(s,t),A2(s,t)x(s,t)+B2(s,t)y(s,t)+F2(s,t) >, h^(p^,x,y,t) =<pV(t), m^D{t\y(sbt)-x(gbt)) > -Ki{x{Si,t) -
Qx
hW{pQ\x,y,t) =<p{2\t), V{t)^W{t\x(s0,t)-y(s0,t)) > -K2(y(s0,t) -e2)2.
QV
Потребуем, чтобы функции ф(s,t) и p(t) = (p(i) ,p(2)) удовлетворяли сопряженной задаче:
Ms,h) = 0, =
ci Qv
^ +С2^ = -в2(3,г)ф2,
ф2(з,г1) = о, ф2(*= (8)
c2 Qx
^ = ~ + 2KM'ut) - 0i), PV(ti) = о;
-^- = -д-Ру '(Ч - + 2К2{у{в0^) - 02), рУ '{и) = 0.
Необходимое условие оптимальности в случае ограничений (5) сформулировано в теореме 1.
Теорема 1. Если процесс {Б, W,x,y} является оптимальным в рассматриваемой задаче, тогда всюду на Т выполняется условие
(1),х(в1,г),у(в1 ,ь),ь)Б = 0, ь € Т, н$(р(2)(г),х(во,ь),у(в0,ь),ь)\¥ = о, ь € т,
где р(1)(ь), р(2)(Ь) - решение сопряженной задачи (8).
Необходимое условие оптимальности в случае ограничений (6) сформулировано в теореме 2.
Теорема 2. Если процесс {Б, W,x,y} является оптимальным в рассматриваемой задаче, тогда всюду на Т выполняется условие
оЛШЮ = о,
где р(1)(Ь), р(2)(ь) - решение сопряженной задачи (8).
Доказательства теорем основаны на применении специальных вариаций, обеспечивающих сохранение гладкости управлений и выполнение ограничений [2]: в случае поточечных ограничений (5) вариация:
ПкЕ1к (ь)= Бк (ь + Е1кб£\г)), Wk2k (ь) = Wк (ь + е2кбк\г)), (9)
¿(1) м = (*-*о)(*1 -¿Я1^) ¿(2) ф = (¿-¿оХ^-гЯ2^)
в случае интегральных ограничений (6) -
Бкк (ь) = (1 + е1к5(1\ь))Вк (ь + £ 1 кб£\ь)), Wk2k (Ь) = (1 + ^(Ь)Щк(Ь + е2к6™(Ь)),
(10)
дк к! ' дк К2 '
(Н-^гпвх^т' к (Н-^гпвх^т'
К1 = тх^(ь)\, к2 = тх^](ь)\-
4. Общая схема оптимизационного алгоритма.
1. Выбираются начальные гладкие управления П°(ь), W°(ь), удовлетворяющие поточечным (5) или интегральным ограничениям (6). Пусть с помощью разработанного численного метода на к-й итерации найдены управления Бк(ь), Wк(ь).
2. На управлениях Бк(ь), Wк (ь) находим: (а) решение прямой системы хк(в,ь), ук(в,ь); (б) решение сопряженной системы фк(в,ь), фк(в,ь), р(к^(ь), р(к\ь).
3. На полученных решениях вычисляется значение функционала Jk = Jk (Dk (t),
Jk (t) и шk
Wk (t)) (4), строятся функции ш^ (t) и ш(\t), в случае поточечных ограничений (5)
имеющие вид
ш{к](t) = h™ (p[1)(t),xk (sut),yk (s1 ,t),t)D k (t),
Uk\t) = h($k (Pk\t),xk (so,t),yk(so,t),t)W k (t), а для интегральных ограничений (6) -
d
= jt(hV{(pf\t),xk(s0,tlykM,t))Wk(t)-
Далее проверяется условие оптимальности uk\t) = 0, wk\t) = 0. Если условие оптимальности выполнено, метод заканчивает свою работу.
4. Если данные управления не удовлетворяют условию оптимальности, строятся их гладкие вариации, которые в случае поточечных ограничений (5) задаются по правилу (9), в случае интегральных ограничений (6) - по правилу (10).
Параметры £1k, £2k определяются из численного решения одномерной задачи минимизации
£ik ,e2k : Jk(Dklk (t), Wjk (t))= min Jk(Dk^ (t),W?2 (t))
E1 ,£2E [0,1]
среди значений 1, ■p-,... Случай, когда найденное значение этих параметров очень близко к нулю, соответствует неулучшению функционала на шаге метода.
5. В качестве к + 1-го приближения выбираются
Dk+1(t) = Dkik (t), Wk+1(t) = Wk2h (t),
и итерационный процесс продолжается дальше.
Критерием остановки служит одна из ситуаций (аналогично [2]), полученных на к-й итерации метода:
а) выполнение с заданной точностью необходимого условия оптимальности для функций Dk(t), Wk (t). Так, близость к нулю отвечающих им функций uk\t), (t) в каждой точке t G T можно гарантировать, если справедливы соответствующие неравенства
max |ш (1)(t)| < 10-5, max \ш(2) (t)\ < 10-5; t е т k t е т
б) достижение заданной точности по значению функционала. Поскольку минимальная величина функционала известна и равна J* = 0, то условием остановки может быть, в частности, неравенство
Jk < 10-3;
в) неулучшение величины функционала, полученной на предыдущей (к — 1)-й итерации
Jk — Jk-1 > 10-6.
5. Численный эксперимент. В качестве примера приведем процесс ректификации в колонне К-34, предназначенной, в частности, для сернокислотного алкилирова-ния изобутана бутиленами (разделяемая многокомпонентная смесь сведена к бинарной [1]). В качестве входных данных были выбраны следующие данные одного из режимов работы колонны: Т = [0, 20] - временной промежуток работы колонны; Б = [0, 20] -геометрические размеры колонны; во = 0 - координата испарителя, в! = 20 - координата конденсатора; е\ = 36 м/ч, е2 = 520 м/ч - потоковые коэффициенты жидкости и пара соответственно. Также задавались концентрации компонентов в готовых продуктах, удовлетворяющие начальным условиям, и функции отбора готового продукта в испарителе Ш*(г) и конденсаторе Б* (г). Остальные параметры задачи, такие как поток жидкости в колонне Ь(г), поток пара в колонне V(г) и коэффициенты задачи, определялись по этим функциям.
5.1. Случай поточечных ограничений на управления. Был задан вводимый поток С(г) = 86. Функции отбора готового продукта в испарителе Ш*(г) и конденсаторе Б* (г) были заданы в следующем виде:
Б* (г) = 2еов(г/4) - 3вт(г/2) + 30, Ш*(г) = -2еов(г/4) + 3 вт(г/2) + 56.
Функционал (4) на этом процессе 3* = 0.
Далее решалась оптимизационная задача с начальными управлениями
Б0(г) = -10вшг + 32, Ш0(г) = 10вшг + 54,
удовлетворяющими ограничению (5). Функционал на начальном приближении 30 = 11.0364.
Решение задачи выполнялось описанным выше методом. В таблице представлены результаты расчетов. Метод закончил работу, достигнув заданной точности
Результаты расчетов для случаев поточечных и интегральных ограничений на управления
N 4 О* \¥к{Ь)
Случай поточечных ограничений на управления
1 0 32 32 32 54 54 54
6 2.173 29.011 22.160 29.017 56.989 63.84 56.981
13 5.328 29.099 40.713 29,094 56.901 45.287 56.906
20 8.445 31.618 28,421 31.622 54.382 57.579 54.387
26 11.116 30.126 25.744 30.120 55.874 60.256 55.995
32 13.778 26.383 41.892 26.378 59.617 44,108 59.592
40 17.333 27.197 23.171 27.191 58.803 62.829 58.797
46 20 32.199 32.401 32.194 53.801 53.801 53.799
Случай интегральных ограничений на управления
1 0 40 40 40 61 61 61
7 2.666 31.431 33.882 31.436 58.058 53.918 58.053
15 6.222 26.369 24.004 26.373 55.618 68.974 55.614
19 8 22.323 26.771 22.329 54.469 61.957 54.472
24 10.222 29.517 35.106 29.512 57.139 53.002 57.133
30 12.889 38.182 39.896 38.179 57.318 64.871 57.321
39 16.889 37.974 27.546 37.968 59.084 57.556 59.080
46 20 25.406 25.287 25.401 59.010 58.939 59.005
по величине функционала за 56 итераций (Jк = 0.00084). Отметим, что управления на выходе близки к оптимальным на всей области определения. При этом
max (t)\ = 0.0013, max ^ (t)\ = 0.051.
5.2. Случай интегральных ограничений на управления. Был задан вводимый поток C(t) = 94.07. Функции отбора готового продукта в испарителе W* (t) и конденсаторе D*(t) были заданы в виде
D* (t) = -2 sin(t/4) + 6 cos(t/2) - 6 sin(t) + 34,
W*(t) = cos(t/2) - 4sin(t/6) - sin(t) + 60.
Функционал на этом процессе J* = 0.
Далее решалась оптимизационная задача с начальными управлениями
D0(t) = 8cos(t/2) + 32, W0(t) = -8sint + 61,
удовлетворяющими интегральному ограничению (6). Функционал на начальном приближении J0 = 5.46.
Как видно из таблицы, алгоритм закончил работу, достигнув заданной точности по величине функционала за 123 итерации (Jк = 0.00076). Управления на выходе близки к оптимальным на всей области определения. При этом
max\^k!)(t)\= 0.0031, max\^<k)(t)\=0.0019.
6. Заключение. Проведенные численные эксперименты для оптимизации процесса ректификации показали, что предложенный метод улучшения гладких управляющих воздействий, стесненных поточечными или интегральными ограничениями, может эффективно применяться для решения задач оптимального управления начально-краевыми условиями полулинейных гиперболических систем.
Литература
1. Демиденко Н. Д., Потапов В. И., Шокин Ю. И. Моделирование и оптимизация систем с распределенными параметрами. Новосибирск: Наука, 2006. 551 с.
2. Аргучинцев А. В. Оптимальное управление гиперболическими системами. М.: Физматлит, 2007. 168 с.
3. Аргучинцев А. В., Поплевко В. П. Оптимальное управление граничными условиями гиперболической системы на примере задачи химической ректификации // Труды XV Байкальской междунар. школы-семинара «Методы оптимизации и их приложения», п. Листвянка, оз. Байкал, 23—29 июня 2011 г. Иркутск: РИО ИДСТУ СО РАН, 2011. Т. 3. С. 36-40.
Статья рекомендована к печати проф. Д. А. Овсянниковым. Статья принята к печати 26 апреля 2012 г.