УДК 66.011
Г. М. Островский, Т. В. Лаптева, Н. Н. Зиятдинов, А. С. Сильвестрова
ВЫЧИСЛЕНИЕ ОГРАНИЧЕНИЙ ПРИ ОПТИМИЗАЦИИ ХИМИКО-ТЕХНОЛОГИЧЕСКИХ
СИСТЕМ С УЧЕТОМ НЕОПРЕДЕЛЕННОСТИ
Ключевые слова: оптимизация с учетом неопределенности, вероятностные ограничения, жесткие ограничения, оптимизация химико-технологических систем
Проектирование оптимальных химико-технологических систем при учете неопределенности исходной физической, химической и экономической информации должно обеспечивать получение конструкции, обеспечивающей выполнение предъявляемых к ней требований с заданной вероятностью либо жестко за весь период эксплуатации. В задачах проектирования оптимальных технологических систем такие требования представлены в виде жестких или вероятностных ограничений. Расчет таких ограничений сопряжен со значительными вычислительными затратами, возникающими из необходимости вычисления многомерных интегралов и решения задач многоэкстремальной недиффернцируемой либо полубесконечной оптимизации. В статье рассматриваются способы снижения вычислительных затрат на вычисление ограничений в задачах оптимального проектирования с учетом неопределенности в исходной информации.
Key words: optimization under uncertainty, chance constraints, hard constraints, process system optimization
Effective chemical process design under physical, chemical, economic uncertainty must get a flexible process, which provides satisfaction of design requirements for the process operation period. This requirements take the form of probabilistic or hard constraints. The calculation of probabilistic constraints requires the calculation of multidimensional integrals at the each optimization procedure step. The calculation of hard constraints related to the solution of the multiextremal undifferentiated optimization problem or semi-infinite one. We develop the approach that reduces CPU-time for solving of the optimization under uncertainty problem.
Введение
Учет неопределенности в исходной информации при проектировании оптимальных химико-технологических систем (ХТС) проявляется как в формализации целевой функции, так и в виде ограничений соответствующей задачи оптимизации
т\п ^си, 6) (1)
с1,ъеН
д^си6) < 0 , ] = 1,...,т, (2)
где ^ (с1,ъ,6) - некоторая функция характеризующая эффективность работы ХТС за рассматриваемый период функционирования, С - пс -вектор конструктивных параметров, ъ - пъ -вектор управляющих переменных, 6 - п6 -вектор неопределенных параметров, о которых известно, что они принадлежат некоторой области Т, Н -область изменения конструктивных и управляющих параметров.
Очевидно, что решение задачи (1) невозможно, поскольку на этапе проектирования неизвестны точные значения параметров 6 , изменяющихся на этапе функционирования независимо от нашего желания.
При формализации задачи проектирования оптимальных ХТС с учетом неопределенности в исходной информации необходимо учитывать, что, в зависимости от полноты информации на этапах проектирования и функционирования ХТС, проектировщик может сталкиваться с двумя основными типами параметрической
неопределенности в исходной информации [1]:
а) интервальная неопределенность -предполагается, что на этапе проектирования доступна информация только о диапазонах изменения значений параметров 6; е Т, I = 1,...,П6 , формирующих неопределенность в исходной информации;
б) вероятностная неопределенность предполагает, что на этапе проектирования известны закон распределения и его характеристики для неопределенных параметров 6 .
Тип неопределенности во многом определяет формализацию решаемой задачи оптимизации:
а) Случай интервальной неопределенности. В качестве функции цели можно использовать либо наихудшее значение оценки работы ХТС за период функционирования
тахОДъ, 6),
6еТ
либо некоторую оценку среднего значения функции 1Хс1,ъ, 6) за это период, полученную на основе имеющихся экспертных сведений
^ • ОДъ, 61),
1е1
где I - множество номеров точек 61, заданных на основе эмпирических сведений, а - весовые коэффициенты, соответствующие эмпирическим точкам. Следует отметить, что, используя данный критерий, мы полностью полагаемся на эмпирические знания эксперта, которые не всегда удачны [2].
При наличии интервальной
неопределенности ограничения задачи
проектирования могут принять вид как отдельных жестких ограничений
gj(d,z(0),0) < 0, j = 1,...,m, V0e T,
(3)
так и могут быть представлены требованием гибкости проектируемой ХТС, представляемым ограничением на значение функции гибкости , которое не должно быть положительным
х(й) = maxminmaxgj(d,z(8),8) < 0 , ] = 1,...,m. (4) 0eT z ] J
Вычисление значения ограничения (3) связано с решение задачи полубесконечного программирования, для чего можно использовать метод внешней аппроксимации [3]. Для получения значений ограничений вида (4) необходимо решать задачу многоэкстремальной недифференцируемой оптимизации, для чего нами в [4] был предложен алгоритм разбиения и границ.
б) Случай вероятностной неопределенности. В качестве функции цели можно использовать математическое ожидание значения оценки работы ХТС за период функционирования
E0 [f(d, z, 0),T] = J f(d, z, 0)p(0)d0,
где р(0) - плотность распределения параметров 0 .
При наличии вероятностной
неопределенности ограничения задачи
проектирования могут принять вид отдельных вероятностных ограничений
Pr{gj(d,z,0)<0}>аj, j = 1,...,
m,
(5)
где Pr{gj (d, z, 0) < 0} - вероятность удовлетворения ограничения gj(d,z, 0) < 0
Pr{gj(d,z,0) < 0} = f p(0)d0>a j
Qj = {0 :gj(d,z,0) < 0, 0e T}.
(6)
(7)
Использование вероятностных ограничений в задачах оптимизации с учетом неопределенности в исходной информации впервые было рассмотрено в работе [5]. Учет вероятностных ограничений при проектировании в области химической технологии предложен Неппоп К с соавторами, РеЛоу 8.Б. с соавторами [6, 7].
Очевидно, что жесткие ограничения (3) или (4) дадут более затратную конструкцию, однако обеспечат выполнение всех проектных ограничений при изменении условий эксплуатации. Вероятностные ограничения (5) при некоторых условиях эксплуатации не будут удовлетворяться,
однако дадут более экономную конструкцию. При этом, то, при каких условиях функционирования вероятностные ограничения не будут выполняться, зависит и от требуемого уровня вероятности a j.
Как мы видим, вычисление ограничений вида (5) потребует проведения операции многомерного интегрирования на каждом шаге процедуры оптимизации, что даже при использовании современных методов численного интегрирования приведет к значительным вычислительным затратам.
Известны три пути преодоления таких сложностей. Первый подход основан на использовании модификации Гауссовых квадратур. Acevedo и Pistikopoulos [8] рекомендуют три альтернативные схемы. Для нормально распределенных параметров Bernardo с соавторами [9] предложили специальную квадратурную формулу для вычисления многомерного интеграла, которая значительно снижает число требуемых точек (узлов) вычисления функции. Известны работы [10; 11], в которых обсуждаются методы приближенных вычислений многомерных интегралов.
Второй подход основан на использовании методик дискретизации (Монте-Карло, Латинский гиперкуб, дискретизирующие последовательности Hammersley (HSS)) [12], [13]. Авторы показали, что среди всех подходов дискретизации метод HSS наиболее эффективен. Поскольку число аппроксимационных точек, требуемых для получения заданной точности в подходе HSS, не растет с размерностью вектора неопределенных параметров 0i , он применим в случаях большого числа неопределенных параметров. К сожалению, даже подход HSS требует несколько сотен аппроксимационных точек для получения приемлемой точности.
Третий подход состоит в преобразовании вероятностных ограничений в детерминированные. В работе Marañas [14] предложено преобразование для случая нелинейной зависимости ограничений от поисковых переменных и линейной зависимости от неопределенных параметров. Группа ученых под руководством Wozny развили стратегию решения нелинейной одноэтапной задачи оптимизации с вероятностными ограничениями для различных случаев. В работе [15] Wendt с соавторами предложили метод решения стационарной одноэтапной задачи оптимизации с вероятностными ограничениями для случая монотонной зависимости между выходами ограничений и входными значениями неопределенных параметров. Это свойство позволяет определить местоположение вероятностной области в области изменения неопределенных параметров. В результате такого определения вероятность выполнения выходных ограничений может быть вычислена интегрированием совместной функции плотности неопределенных параметров.
В статье будет предложен подход к решению одноэтапной задачи оптимизации (ОЭЗО)
[1], основанный на преобразовании вероятностных ограничений в детерминированные для случая, когда неопределенные параметры статистически независимы и имеют нормальное распределение.
Подходы к вычислению вероятностных ограничений
ОЭЗО с вероятностными ограничениями имеет вид [16]
тт Е0[1^, 0)] d,zeH
0)<0}>аj, ] = 1,...,т .
(8)
(9)
В [16] мы предложили использовать линеаризации д^у,0,01) для функций ограничений дj z, 0), получаемые из их разложения в ряд Тейлора по неопределенным параметрам в точке 01 е Т с отбросом нелинейной части
^(у, 0,01) = gj(y, 01)
, ^^^0|) ,0. 01
i=1
50i
i-0i).
(10)
Мы показали, что в этом случае вероятностное ограничение (9) может быть заменено ограничением вида
У (0^ - Е0Л) I1 - gj(y, 0й)
-^->Ф-"Ц), (11)
1
П0 ^
^ и Э0i i=1
где Ф(а)
функция нормированного,
л-1
распределения, Ф (а) - функция, обратная к -1
Ф(а), (Ф(Ф (а)) = а). Значения функции Ф(а) могут быть найдены из таблиц распределения [17]. Тогда вычисление значения ограничения (11) не будет требовать операции интегрирования.
Подставив (11) в ОЭЗО (89) получим задачу
тт Е0[1^, 0)] d,zеH
(12)
У (^ - вдЦ- - ^0и) i=1 i
1
У 40^
■ л ^
■>Ф-1(аj).
(13)
i=1
Задача (12) детерминированного
является задачей нелинейного
программирования и может быть решена хорошо известными методами.
Очевидно, если функции gj(d,z, 0)
ограничений (9) выпуклы по неопределенным параметрам 0 , то задача (12) даст нижнюю оценку критерия ОЭЗО (8). Для улучшения аппроксимации (10) нами в [16] предложена итерационная процедура, основанная на разбиении области неопределенности. В этом случае, при достаточном разбиении получаемые подобласти области неопределенности будут малы и, в большинстве случаев, решив (12), мы получим нижнюю оценку критерия (9).
Поскольку предложенный подход ограничен требованием на выпуклость функций ограничений по неопределенным параметрам, нами в [18] был предложен подход, основанный на аппроксимации областей выполнения ограничений (9). Задача (8) была сведена нами к задаче
тт Е0 [%, 0)]
уеН,Таj е^аj
тах дj(y,0) < 0, j = 1,...,т .
0еТа j
Рг{0 е Таj} >аj , Таj с Т , j = 1,...,т ,
(14)
(15)
(16)
где Та j - множества областей Та j,
удовлетворяющих (9). Очевидно, что решение задачи (14) потребует поиска формы, размера и местонахождения областей Та j , что является очень
трудоемкой операцией. Мы предложили в [19]
ч
аппроксимировать области Та j многомерными
прямоугольниками
Т. = {0i: 0^ ^ <0^ = 1.....П0},
(17)
а j
0L,j 0U,j
включив их границы 0| , 0| в число поисковых
переменных задачи. В [19] мы показали, что ограничение (16) можно заменить ограничениями
П0 „ j
П[Ф(0^) -Ф(8¡L,j)] > аj i=1
М 0U,j <0и
0|- <0
где
0^ < 0" , i = 1.....П0 ,
(18)
(19)
(20)
0Г-E[0i] eU,j У -E[0i] i ' i ^ '
2
E[0i] - матожидание, (стi) - дисперсия неопределенного параметра 0i, i = 1,...,П0 .
Тогда можно переписать задачу (8) в виде
min E0 [f(y, 6)]
max g¡(y,0) < 0, j = 1,...,m.
6eTaj
ne
^[O(0.U,j) _0(e.L,j)] > aj, j = 1,...,m, i=1
(21)
(22)
(23)
| д^си 6)р(6)С6>а j. (28)
Ра j
Рассмотрим более подробно вычисление левой части ограничения (28). Можно записать
|gj(d,z,6)р(6)<С6 = I1,j • I2,j... • 1п6-1,j • In6,j, (29) Ра j
eL <eL,j, eU,j<eU, i = i,...,ne, j=i,...,m. (24)
Задача (21) является задачей полубесконечного программирования, для ее решения нами в [20] была разработана модификация метода внешней аппроксимации [3].
Поскольку мы сузили вид областей Та j,
задача (21) будет давать верхнюю оценку критерия задачи (8). При этом предложенный подход преобразования ограничений не требует выпуклости функций ограничений по неопределенным параметрам.
Для улучшения получаемой оценки нами в [20] была предложена итерационная процедура улучшения аппроксимации областей Та j
многомерными прямоугольниками Таj , за счет
разбиения Таj на подобласти Ту, I = 1,...,^.
Однако, при увеличении количества подобластей Ту возрастает размерность задачи (21).
Поэтому нами предложено
аппроксимировать области Та у криволинейными
областями Рау, полученными из областей Тау
путем замены плоскости Бп6 у, представляющей одну из сторон многомерного прямоугольника Тау , гиперплоскостью ^6,у. Будем получать гиперплоскости Бп6,у, выразив неопределенный параметр 6п6 из уравнения
где
gj(d,z,e) = 0 .
(25)
Тогда гиперплоскость S^ j будет описана уравнением
ene =Фj(d,z, e1,..., ene-1).
(26)
Ii,j = JP(ei)dei, i = 1,...,ne-1,
Sni,j
Ine,j = jp(ene)dene .
(30)
(31)
S
ne,j
Для одномерных интегралов (30) верно Ii,j = ®(0jU,j)-ф(~и), i = 1,...,ne-1,
поэтому для вычисления (29) потребуется вычисление только одного одномерного интеграла
(31).
Для улучшения аппроксимации области Тау криволинейными областями Рау, будем
разбивать области Рау на подобласти Ру,
I = 1,...,Иу. Тогда задача (8) может быть записана в виде
min Ee[f(d,z, e)] d,zeH
Ni
j f ne-1
z
l=1
где
n
i=1
®(0jU,j,1) -®(eL,jJ)
4JJ
> a j.
(32)
(33)
eL,j,i = eL;j;'-E[ei] eu,j,i eU,jJ-E[ei]
i
CTi
i
CTi
(34)
eL,j,1, ep^1 - границы области R|j.
Полученная задача является задачей нелинейного программирования.
Вычислительный эксперимент
Будем аппроксимировать ее плоскостью, проходящей через ее центр. Тогда мы можем переписать задачу (8) в виде min Ee[f(d,z, e)] (27)
d,zeH
Рассматривается задача проектирования системы реакторов, представленной на рисунке 1 [21]. В реакторах 1 и 2 протекают реакции превращения вещества А в вещество В и, далее, в вещество С.
Математические модели реакторов имеет
вид:
Реактор 1
СА1 + к1СдМ = 1; СБ1 + СА1 + кзСв11^1 = 1 к1 = к10е"Е^КТ1 ; кз = к2ое"Е^КТ1;
Реактор 2
СА2 - СА1 + к2СА2^2 = 0 ;
СБ2 - СБ1 + СА2 - СА1 + к4СБ2^ = 0 :
к2 = кюе-Е^КТ"2 ;
к4 = к2ое-Е^КТ2.
А ^ В ^ С
А ^ В ^ С
к
к
С = 1 А0 1 V1 С А1 С т 2 С А2 С В 2
Ш 0
Рис. 1 - Система из двух реакторов
В качестве неопределенных параметров выбраны 0 = {Е1; Е2; кю; к2о}. Характеристики неопределенных параметров:
Е1: 0М =6665,948; 5 =200; Е2: 0М =7985,248; 5 =240; к10 : 0М =0,715; 5 =0,0215;
к20 : 0М =0,182; 5 =0,0055.
Неопределенные параметры
предполагаются независимыми и подчиняющимися нормальному закону распределения. Область неопределенности имеет вид многомерного
параллелепипеда {0^ - к - 5 < 0 < 0^ + к -5}.
Целевая функция представляет собой капитальные затраты
*=лМ". (35)
Ограничения задачи имеют вид:
0 < СА1 <1 0 < СБ1 < 1 0 < СА2 < 1
0 < СБ2 <1
Св2 ^ СБ . (36)
Все ограничения, кроме (36), жесткие, ограничение (36) - вероятностное. В качестве поисковых переменных выбраны: конструктивные переменные d - объемы реакторов 1 и 2 - V], У2; управляющие переменные ъ - температуры в реакторах 1 и 2 - Т|, Т2. Диапазоны изменения поисковых переменных:
0 < V < 16
601,4 < Т1 < 661,53
0 < V2 < 16
541,26 < Т2 < 601,4 .
Для этой ХТС была решена ОЭЗО с использованием подхода, основанный на получении верхней оценки решением задачи (21), а также подход, основанный на преобразовании вероятностных ограничений к виду (33).
В таблице 1 приведены требуемое значение СБ концентрации целевого вещества Б на выходе
ХТС и значений вероятности выполнения ограничения (36), для которых решалась задача, а также полученные значения критерия * и
затраченное на получение решения время 1, сек индексы 1 - верхней оценке ОЭЗО, 2 - подход,
(индекс А соответствует решению авторов [21]; использующий задачу (32)).
Таблица 1 - Результаты решения задачи оптимизации, верхняя оценка
* Св a fA f1 t1 f2 t2
0,50 0,90 3,62 2,96 96,4 2,71 0,22
0,50 0,95 3,67 3,07 58,3 2,75 0,26
0,52 0,90 3,9 3,58 1 1 248 3,03 0,24
0,52 0,95 3,96 3,83 1 2224 3,13 0,43
Анализируя полученные результаты, заметим, что использование предложенных подходов дает решение лучшее, чем полученное авторами [21]. При этом заметно преимущество подхода, использующего задачу (32), перед подходом, дающим верхнюю оценку, как в найденном значении критерия, так и во времени решения задачи. Отметим, что достигнуто сокращение времени в «50^15000 раз.
Заключение
В статье рассмотрены подходы к решению задач проектирования гибких ХТС, представленных в форме ОЭЗО с вероятностными ограничениями. Сложность решения таких задач состоит в необходимости вычисления многомерных интегралов на каждой итерации метода при прямом решении.
Использование предложенных в статье подходов для решения ОЭЗО показало их преимущество перед известными способами решения ОЭЗО, при этом предложен новый подход, значительно сокращающий время решения задачи, что подтверждено результатами вычислительного эксперимента.
Литература
1. Г.М. Островский, Н.Н. Зиятдинов, Т.В. Лаптева. Оптимизация технических систем. М.: КНОРУС, 2012.
2. Т.В. Лаптева, Н.Н. Зиятдинов, Д.Д. Первухин. Вестник Казан. технол. ун-та, 12, 216-219, (2012).
3. E.Polak, D.Q.Mayne. IEEE Transactions on Automatic Control, 21, 2, 184-193, (1976).
4. Г.М. Островский, Н.Н. Зиятдинов, Т.В. Лаптева, И.Д. Первухин. Вестник Казан. технол. ун-та, б, 199-206, (2011).
5. A.Charnes, W.W. Cooper. Management Science, 6, 73-79, (1959).
6. R.Henrion, A.Moller. Comput. Math. Appl., 45, 247-262,
(2003).
7. S.B.Petkov, C.Maranas. Ind. Eng. Chem. Res., 36, 48644881, (1997).
8. J. Acevedo, E.N. Pistikopoulos. Comput. Chem. Eng., 22, 647-671, (1998).
9. F.P. Bernardo, E.N. Pistikopoulos, P.M. Saraiva. Ind. Eng. Chem. Res., 38, 3056-3068, (1999).
10. N. Pintaric, Z. Kravanja. Chem. Eng. Sci., 28, 6-7, 1087,
(2004).
11. J.Wey, M.J. Realff. Comp. Chem. Eng., 28, 3, 319, (2004).
12. F.P. Bernardo, P.M. Saraiva. AIChE J., 44, 2007-2117, (1998).
13. U.M. Diwekar, J.R. Kalagnanam. AIChE J., 43, 440-447, (1997).
14. C.D. Maranas. AIChE J., 43, 1250-1264, (1997).
15. P.Li, H.Arelano-Garcia, G.Wozny. Comp. Chem. Eng., 32, 25-45, (2008).
16. Т.В. Лаптева, Г.М. Островский, Н.Н. Зиятдинов, Д.Д. Первухин. Вестник Казан. технол. ун-та, 14, 7, 218-224, (2011).
17. Г. Крамер. Математические методы статистики. М.: РХД, 2003.
18. Т.В. Лаптева, Н.Н. Зиятдинов, Д.Д. Первухин, Г.М. Островский. Вестник Казан. технол. ун-та, 14, 9, 281-287, (2011).
19. Т.В. Лаптева, Н.Н. Зиятдинов, Г.М Островский, И.В. Зайцев. Докл. АН, 435, 4, 497-500, (2010).
20. Г. М. Островский, Т. В. Лаптева, Н. Н. Зиятдинов. ТОХТ, 48, 5, 1-11, (2014).
21. M. Wendt, P. Li, G.Wozny. Ind. Eng. Chem. Res., 41, 3621-3629, (2002).
© Г. М. Островский - д.т.н, проф. каф. системотехники КНИТУ, ostralex@yandex.ru; Т. В. Лаптева - к.т.н., доц. той же кафедры, tanlapteva@yandex.ru; Н. Н. Зиятдинов - д.т.н., проф., зав. каф. системотехники КНИТУ, nnziat@yandex.ru; А. С. Сильвестрова - ассистент, аспирант той же кафедры sensoriumsa@mail.ru.
© G. M. Ostrovsky - full professor, professor of Process system engineering department, KNRTU, ostralex@yandex.ru; T. V. Lapteva - PhD, ass. professor of Process system engineering department, KNRTU, tanlapteva@yandex.ru; N. N. Ziatdinov -full professor, professor of Process system engineering department, KNRTU, nnziat@yandex.ru; A. S. Silvestrova - assistant, graduate student of Process system engineering department, KNRTU, sensoriumsa@mail.ru.