Учет неопределенности при проектировании оптимальных химико-технологических систем Текст научной статьи по специальности «Математика»

УДК 66.01

Г. М. Островский, Н. Н. Зиятдинов, Т. В. Лаптева,

И. Д. Первухин

УЧЕТ НЕОПРЕДЕЛЕННОСТИ ПРИ ПРОЕКТИРОВАНИИ ОПТИМАЛЬНЫХ ХИМИКО-ТЕХНОЛОГИЧЕСКИХ СИСТЕМ

Ключевые слова: оптимизация, неопределённость, оптимальное проектирование, химико-технологические

системы.

Проектирование химико-технологических систем (ХТС) обычно выполняется при частичной неопределенности исходной физико-химической, технологической и экономической информации. Вследствие этого возникает необходимость учета условия гибкости проектируемой ХТС. Особенностью этого условия является наличие maxminmax процедуры, которая делает задачу оптимального проектирования задачей недифференцируемой многоэкстремальной оптимизации. В статье предлагаются два подхода к решению задачи, основанные на методах «разбиений и границ» и «внешней аппроксимации». Решение задачи оптимального проектирования продемонстрировано для системы реактор-теплообменник.

Key words: optimization, uncertainty, optimal design problem, chemical-technological systems.

Designing of chemical-technological systems is usually taking place in case of partial uncertainty in source physical-chemical, technological and economical information. Therefore arise the need of taking account for flexibility condition of the system being designed. The feature of this condition is the presence of maxminmax procedure, which transforms the optimization problem into nondifferentiable multiextremal problem. We propose two approaches to solve this problem based on the "outer approximation" and "partition and bounds" methods. The "reactor - heat exchanger" system is used to demonstrate the solution of optimal design problem.

При проектировании ХТС стремятся выполнить некоторые требования, например, ХТС должна работать без аварийных ситуаций, должна быть экологически безопасной и т. д. Для краткости, ХТС, соответствующую этим требованиям, будем называть работоспособной. Задачу построения работоспособной ХТС приходится, как правило, решать в условиях неполной физико-химической и технологической информации. Источниками неопределенности являются: неточность коэффициентов в математических моделях, связанная с неточностью эксперимента; неточность химических и физических закономерностей, положенных в основу математических моделей; изменение части коэффициентов в математических моделях во время эксплуатации ХТС; изменение внешних условий функционирования ХТС во время её эксплуатации.

При проектировании основным вопросом является определение оптимальных размеров аппаратов и их режимов. В общем случае задача примет вид:

min f(d, z, х,в), (1)

d ,z,x

ф (d, z, х,в) = 0, i = 1,..., n, (2)

(d, z, х,в) < 0, j = 1,..., m, (3)

где (2) - уравнения материального и теплового баланса ХТС (математическая модель ХТС), (3) - математическая формулировка проектных требований, х - n -вектор переменных состояния, d е D - вектор конструктивных переменных, D - область допустимых значений конструктивных переменных, z - вектор технологических управляющих переменных. Неполнота исходных данных выражается в том, что часть параметров математических моделей, пусть это вектор параметров в, на этапе проектирования известны неточно, но

известно, что они принадлежат некоторой области неопределенности T = (в: в < в < в ).

Если выразить переменные х из (2), то задачу (1)-(3) можно свести к виду

min f (d, z,G), (4)

d,z

(d,z,0) < 0, j = 1,...,m. (5)

Задача (4)-(5) не может быть решена, поскольку параметры в принимают любые значения из области Т. Обычно поступали следующим образом: брали некоторое значение вы, для которого решали задачу (4)-(5). Поскольку нужно было удовлетворять ограничения (5) при любых значениях в е T , то далее на основе опыта специалистов вносили поправку для полученных значений конструктивных переменных. Полученное таким образом значение конструктивных переменных приводило либо к значительному увеличению конструкции, либо могло не обеспечить выполнения условий (5) для всех возможных значений в еТ .

Учёт неопределённости при проектировании ТС может иметь различную форму. Например, неопределённым переменным на основе опыта и интуиции присваиваются некоторые «номинальные» (обычно - средние) значения. Решается задача оптимизации с этими значениями в традиционной постановке, в результате чего определяются номинальные оптимальные величины параметров оборудования (длины и диаметра реактора, поверхностей теплообмена в теплообменниках, числа тарелок в ректификационных колоннах и т.п.). После этого с учётом знаний о процессе, опираясь на имеющийся опыт, вводят так называемые «коэффициенты запаса» и принимают для проектирования величины параметров оборудования, полученные как произведение номинальных оптимальных величин и запасов. С примером такого подхода можно ознакомиться в [1]. Другой вариант - учитывать неопределённость в математической постановке задачи оптимального проектирования.

Рассмотрим подходы к решению задачи оптимального проектирования ХТС, учитывающие неопределенность в математической постановке задачи. В процессе проектирования необходимо создавать гибкие ХТС - такие, для которых на этапе функционирования с помощью управляющих переменных можно удовлетворить все проектные ограничения, несмотря на изменение внутренних и внешних факторов. Целесообразно до решения задачи оптимального проектирования гибкой ХТС предварительно оценить возможность существования ХТС, удовлетворяющей всем наложенным ограничениям на заданной области неопределенности. Таким образом, актуальным является решение следующих задач:

1. Оценка гибкости имеющейся ХТС.

2. Оценка существования решения задачи (4)-(5) на всем диапазоне изменения значений параметров в е Т .

3. Определение оптимальных значений конструктивных переменных, гарантирующих гибкость ХТС на заданной области неопределенности.

При решении задачи оптимального проектирования будем выполнять следующие предположения:

1. В жизни ХТС имеются два этапа - проектирования и функционирования. На этапе проектирования известны диапазоны изменения параметров в . На этапе функционирования можно получить точные значения в .

2. В задаче поисковые переменные можно разделить на два типа: конструктивные d и управляющие z . На этапе функционирования переменные d остаются постоянными, а управляющие z могут изменяться. В этом случае речь идет о двухэтапной задаче (ДЭЗО). Если считать переменные d и z равноправными, то мы говорим, что найденные на этапе проектирования переменные z не могут изменяться на этапе функционирования при изменении параметров в . Это может привести к удорожанию всей конструкции. Тогда речь идет об одноэтапной постановке задачи оптимизации. В дальнейшем будем рассматривать ДЭЗО.

3. На этапе функционирования будет проводиться оптимизация ХТС на основе математических моделей с уточненными значениями параметров в .

4. На этапе функционирования ХТС ограничения (5) должны безусловно выполняться при любых значениях неопределенных параметров ве T .

Задача оценки гибкости ХТС рассматривалась в [2]. Здесь приведем только условие гибкости ХТС. ХТС является гибкой, если

z(d) = maxminmaxу,(d,zß), %(d) < 0, J = {j: j = 1,m} (6)

веТ z jeJ 1

где x(d) - функция гибкости ХТС [3].

Запишем логическое условие, выполнение которого гарантирует существование гибкой ХТС. Необходимо проверить, существует ли такое значение конструктивных переменных d, что для любого значения ве T будут найдены управляющие переменные z , позволяющие удовлетворять все ограничения (5).

Используя утверждения, приведенные в [4], можно записать вышеуказанное логическое условие в виде

E = minmaxminmax w.(d,z,в) = min r(d) < 0. (7)

deD веТ z jeJ 1 deD

Задача (7) является задачей недифференцируемой многоэкстремальной оптимизации. Рассмотрим подход, который обеспечит получение значения величины E на основе вычисления верхней EU и нижней EL оценок этой величины, а также разбиения области Т на подобласти.

Вводя дополнительную переменную u, приведем задачу (7) к виду

E = min u (8)

dsD

maxminmax у, (d, zfi) < u (9)

веТ z jeJ 1

Найдем верхнюю оценку EU . Поменяем местами операции максимизации и минимизации, получим задачу вида

EU = min u, (10)

d,u

minmaxmaxw:(d,zfi) < u. (11)

z веТ jeJ 1

Можно показать, что EU будет верхней оценкой величины E. В [3] показано, что задачу (10)—(11) можно свести к виду

EU = min u (12)

d ,z,u

maxy, (d, z,в) < u, j = 1,m. (13)

веТ ■'

Пусть имеется разбиение области Т на N подобластей Т(, l = 1,N таких, что Т П Т2 П... П ^ = 0 и Т U Т2 U... U Т^! = Т. Тогда можно показать, что задача

EU= min u, (14)

d ,zl ,u

max у/,(d, Zi,0) < и, j = 1, m, l = 1, Nr, (15)

где r - номер итерации в алгоритме вычисления значения E, есть улучшенная верхняя оценка задачи (7).

Задача (14)-(15) представляет собой задачу полубесконечного программирования. Ее можно решать методом внешней аппроксимации (ВА) [5].

Алгоритм 1. Метод ВА для решения задачи (14)-(15).

Шаг 1. Задать: номер итерации k = 1; разбиение области T| (l = 1, Nr); значения d0, zf, и0; множества критических точек S°, соответствующие областям Tf. Здесь r - номер итерации метода, используемого для вычисления значения величины E. Шаг 2. Построение множеств S/k . Решить m ■ Nr задач вида

Ч.. = max{dk \гк-\в), j = 1,m, l = 1, N.

' ,j беГ. j ' '

Если для некоторых (j ,l ) выполняется условие T,. - u - > 0, то соответствующие критические точки . j. заносим в множество R. Построить множества S/k , Sk = Sf-1 U R/k, l = 1, . Перейти на шаг 3.

Если не было выявлено ни одной критической точки, то решение задачи найдено. Построить множество Vr номеров областей, в которых присутствуют активные ограничения, то есть, если для некоторых (j , l ) верно . - uk-1 = 0 , то заносим такой номер области l в

множество Vk. Stop. Шаг 3. Решить задачу

EULk = min u (16)

d ,z, ,u

¥j (d, z, ,O!q ) < u, (17)

j = 1m, I = Wr, Oq g Sk, q = ,

k ok i—r jk k

где Pi - количество критических точек в множестве S, . Пусть d , zt есть решение задачи (16)-(17).

Шаг 4. Присвоить k = k +1. Перейти на шаг 2. Запишем задачу

EL'r = min u (18)

d ,Zq,U

(d, Zq, Oq ) < U , (19)

j = im, Oq g S, q = 1P. В качестве множества S будем использовать множество критических точек, полученное на последней итерации Алгоритма 1. Задача (18)-(19) дает нижнюю оценку задачи (7), поскольку ее допустимая область значений поисковых переменных включает допустимую область задачи (8).

Для вычисления значения E будем использовать алгоритм метода РГ [3], опирающийся на сравнение оценок и дробление области неопределенности.

Алгоритм 2. Метод РГ для решения задачи (7)

Шаг 1. Задать: номер итерации r = 0, разбиение области T°, I = 1,Nr, малое число s1 > 0;

число s2 > 0, множество Q0 =0 разбиваемых областей, значения d0 , z° . r = r +1.

Шаг 2. Вычислить Ef,r (алгоритм 1, задача (16)-(17)). Если выполняется Ef,r < 0, то

очевидно, что и E < 0, тогда ХТС существует. Stop. Иначе перейти на шаг 3.

Шаг 3. Найти значение нижней оценки E!^,r решением задачи (18)-(19).

Шаг 4. Если верно Ef,r - EM < s1, решение найдено. Stop. Иначе перейти на шаг 5.

Шаг 5. Для УТ[ gV проверяем ее размер ô(T[ ). Если верно ô(T[ ) > s2, УТ[ GVr, то

s

заносим такую область Т[ в Qr. Если Qr =0, то решение найдено. Stop. Иначе s2 = -2- и переходим к шагу 6.

Шаг 6. Каждую область Т[ g Qr делим на две подобласти. Получаем новый набор подобластей Т[+1, I = 1, Nr+1 . r = r +1. Перейти на шаг 2.

Ранее мы предположили, что на этапе функционирования конструктивные параметры d известны и неизменны, значения параметров O g Т могут быть точно измерены, значения переменных z вычисляются на этапе функционирования решением задачи оптимизации

f\d,0) = min f (d, z,d)

z

y(d,z,0) < 0, j = \m. В качестве критерия задачи ДЭЗО будем использовать величину E (f '(d,9)) = \ f '(d,0)u(0)d0, где ß(0) - функция совместной плотности вероятности

Т

переменных в, E(f (d,9)) - математическое ожидание величины f (d,9) на области T - дает

среднее значение f (d, в) на всем этапе функционирования.

Учитывая требование работоспособности ХТС, запишем задачу ДЭЗО в виде:

f = min E(f *(d,e)) (20)

d

x(d) < 0. (21)

Перейдем к дискретному варианту записи задачи (20)-(21) [3], а также учтем вид функции гибкости

f1 = rmin 2 wf (d, z, ,в) (22)

dz и

y (d, z, ,в) < 0, j = Im, Ig I, (23)

max ОД в) < 0, (24)

веТ

где I - множество номеров аппроксимационных точек, по которым построен дискретный вариант критерия, w■ - веса, а функция h(d, в) = min max у/ (d, z,e).

z jeJ 1

При решении задачи (22)-(24) можно предварительно решить задачу существования решения. Тогда полученное решение и разбиение области неопределенности можно взять как стартовую точку при решении задачи (22)-(24).

Для решения задачи (22)-(24) можно использовать метод ВА (Алгоритм 1). При этом на шаге 2 алгоритма нужно решать задачу H (d) = max h(d^) для формирования множества

веТ

критических точек Sk. Это задача поиска значения функции гибкости [2]. На шаге 3

проводится решение варианта задачи (22)-(24), в котором условия (24), с учетом разбиения области Т на подобласти Т,, заменены конечным количеством условий

y (dk, zqв) < 0, p = й j = im, I e I, в,q e Sk (25)

Рассмотрим другой подход к решению задачи (22)-(24). В [3] было предложено решать задачу методом РГ. Для использования этого метода нужно уметь вычислять верхнюю и нижнюю оценки задачи.

Верхняя оценка для задачи (22)-(24) с учетом разбиения области Т на подобласти Т;

(I = 1,Nr) на r -ом шаге метода РГ может быть записана в виде [3 ]:

fiU r = min 2 wff (d, z, Д) (26)

d z -z lie _

y(d, z, в) < 0, j = im, I e I, (27)

maxy(d,z,,в) < 0, I = \Nr. (28)

веТ,

На основе полученных множеств критических точек можно вычислить нижнюю оценку задачи (22)-(24)

fiL,r = min 2 wf (d, z, ,в,) (29)

d,zi ,zq iel

y/j (d, zi ,в) < 0, ig I, I = 1, Nr, j = 1 m, (30)

y(d,zlq,вq) < 0, e S, . (31)

Алгоритм решения задачи ДЭЗО методом РГ аналогичен Алгоритму 2.

Рассмотрим применение предложенных подходов на примере ХТС реактор-теплообменник [6]. В реакторе идеального смешения объемом V происходит экзотермическая реакция первого порядка вида А ^ В . Теплообменник служит для поддержания температуры внутри реактора ниже заданной.

Рис. 1 - Технологическая схема процесса

Р0,Т0,СА0 - расход (м3/ч), температура (К) и концентрация реагента А (кмоль/м3) в потоке в реактор, соответственно; V,Т1,СА1 - объем реактора (м3), температура (К) и концентрация

33

реагента А (кмоль/м ) в реакторе, соответственно; ^2,Т2 - расход (м /ч) и температура рециркуляционного потока (К), соответственно; ^,Тм1,Тм2 - расход (м /ч), входная и выходная температуры потока охлаждающей воды (К), соответственно. Математическая модель задачи имеет вид [3]:

VkQCA0 ехр(-Е / RT1)

conv =

Ъ + VkoCAo ехр(-Е / RT1), Оме С - ОМЕ

=-^-, =

Ср (Т1 Т2) Срм/ (Тм2 Т\н 1)

= 2(-М- 2ГрСр(Т - Т,) -

'2 = А и ' 1 + 'м 2 + ' '

О = Аи(Т1 Тм2) , (Т Т )

0МЕ =-2-+ ( ' 2 - ' ) ,

где М - теплота реакции, (кДж/кмоль); Ср и Срт - теплоемкости рецикловой смеси и

охлаждающей воды, (кДж/(кг К)), соответственно; А - поверхность теплообмена в

2 2 теплообменнике, (м ); и - коэффициент теплопередачи, (кДж/(м ч К)), kR -скорость

реакции, (м3/(кмоль ч)).

В задаче V и Д - конструктивные переменные, а Т1 и Тм2 - управляющие

переменные. Переменные состояния - СА1, Т2, и . Неопределенные параметры

в=[Я-0,Т0 ,Тт 1,kR, и ]. Размер области неопределенности задается как

Т(у) = 00 -у80! )<0 (1 + у80)}, и характеризуется параметром у, величиной 8

отклонения параметра 0 от его номинального значения в? при у = 1.

Таблица 1 - Отклонение неопределенных параметров от номинального значения

Параметр 70 7Ш1 к* и

Номинал 45.36 393 300 К 9.81 1635.34

5 0.1 0.02 0.03 0.1 0.1

Критерий оптимальности имеет вид

^ = 691.2407 + 873.6Л06 +1.76^ + 7.056Г2. Ограничения задачи имеют вид

V > 0, А > 0, 0.9 < оопу < 1, 72 - 71 < 0, 7т 1 - 7т 2 < 0, Г„1 - 72 +11.1 < 0, 7,2 - 71 +11.1 < 0, 311 < 71 < 389, 311 < 72 < 389, 301 < 7т2 < 355 . Рассматривалась двухэтапная задача поиска оптимальных значений конструктивных переменных V и А( при наличии неопределенности в исходных данных. Задача была решена

на основе предложенных подходов без задания предварительного разбиения области неопределенности для различных у. Результаты решения задачи приведены в табл. 3. При решении задачи методом РГ решение было получено не для всех значений параметра у. Поэтому для получения решения ДЭЗО методом РГ предварительно была решена задача (7). Найденные значения конструктивных параметров V и А(, а также значение функции гибкости

для ХТС приведены в табл. 2. На основе решения задачи (7) была решена ДЭЗО методом РГ. Результаты приведены в табл. 3.

Таблица 2 - Результаты решения задачи (7)

у V А Еи

1,00 12,75 14,97 -0,0426

1,25 12,67 14,19 -0,0391

1,50 12,58 12,22 -0,0362

1,75 13,16 15,09 -0,0352

2,00 13,42 13,96 -0,0306

2,50 26,00 12,65 0,7307

Таблица 3 - Результаты решения задачи ДЭЗО методом РГ

у Без п редразбиения С предразбиением

V А Время, сек V А Время, сек

1.00 10937 6.7 8.9 210.47 10937 6.7 8.9 210.34

1.25 11323 7.1 9.6 221.16 11338 7.12 9.52 212.91

1.5 - - - - 11766 7.5 10.4 1442.03

1.75 - - - - 12293 7.9 11.3 2734.72

Таблица 4 - Результаты решения задачи ДЭЗО методом ВА без предразбиения

у V А Время, сек

1.00 10937 6.7 8.9 209.61

1.25 11323 7.1 9.6 212.14

1.5 11766 7.5 10.4 1284.59

1.75 12293 7.9 11.3 1961.0

Сравнивая результаты решения ДЭЗО методами РГ и ВА нетрудно заметить:

- решение задачи ДЭЗО методом РГ лучше проводить при задании на первой итерации предварительного разбиения области неопределенности. Для получения такого разбиения можно использовать решение задачи существования решения (7);

- в отличие от метода РГ, метод ВА позволяет решить задачу проектирования без предварительного решения задачи (7);

- при решении задачи ДЭЗО подход, опирающийся на метод ВА при вычислении %(d) методом РГ, имеет преимущество в вычислительных затратах по сравнению с подходом, решающим задачу ДЭЗО методом РГ.

Литература

1. Зиятдинов, Н.Н. Поиск энергосберегающих режимов работы установки разделения изоамилен-изопреновой фракции производства изопрена / Н.Н. Зиятдинов, Д.А. Рыжов, Т.В. Лаптева, В.А. Курбатов // Вестник Казан. технол. ун-та. - 2009. - № 6. - С. 249-258.

2. Островский, Г.М. Оценка гибкости химико-технологических систем / Г.М. Островский, Н.Н. Зиятдинов, Т.В. Лаптева, И.Д. Первухин // Теоретические Основы Химических Технологий. -

2007. - том 41, № 3. - с. 249-261.

3. Островский, Г.М. Методы оптимизации химико-технологических процессов: учебное пособие. / Г.М. Островский, Ю.М. Волин, Н.Н. Зиятдинов. - М.: КДУ, 2008. - 419 с.

4. Островский, Г.М. Технические системы в условиях неопределенности: анализ гибкости и оптимизация: учебное пособие. / Г.М. Островский, Ю.М. Волин. - М.: БИНОМ. Лаборатория знаний,

2008. - 319 с.

5. Maine, P.Q. An Outer Approximation Algorithm for Computer-Aided Design Problem. / P.Q. Maine, E. Polak, R. Traham // J. Optim. Theory Applies. - 1979. - v.28. - p. 3

6. Halemane, K.P. Optimal Process Design under Uncertainty. / K.P. Halemane, I.E. Grossmann // AIChE Journal. - 1983. - v.29. - p. 425-433.

© Г. М. Островский - д-р техн. наук, проф. каф. системотехники КГТУ, ostralex@yandex.ru; Н. Н. Зиятдинов - д-р техн. наук, проф., зав. каф. системотехники КГТУ, nnziat@yandex.ru; Т. В. Лаптева - канд. техн. наук, доцент той же кафедры, tanlapteva@yandex.ru; И. Д. Первухин -ведущий программист ЦНИТ КГТУ, pervuhin@kstu.ru.