Научная статья на тему 'Анализ гибкости оптимального проекта колонны дебутанизации'

Анализ гибкости оптимального проекта колонны дебутанизации Текст научной статьи по специальности «Математика»

CC BY
140
44
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
РЕКТИФИКАЦИОННАЯ КОЛОННА / ОПТИМАЛЬНОЕ ПРОЕКТИРОВАНИЕ / РАБОТОСПОСОБНОСТЬ / ИНДЕКС ГИБКОСТИ / ФУНКЦИЯ ГИБКОСТИ / DISTILLATION COLUMN / OPTIMAL DESIGN / PERFORMANCE / FLEXIBILITY INDEX / FLEXIBILITY FUNCTION

Аннотация научной статьи по математике, автор научной работы — Зиятдинов Н. Н., Лаптева Т. В., Богула Н. Ю.

Дается постановка и решение задачи определения индекса гибкости для спроектированной колонны дебутанизации. Приводятся алгоритмы вычисления индекса гибкости

i Надоели баннеры? Вы всегда можете отключить рекламу.
iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.
i Надоели баннеры? Вы всегда можете отключить рекламу.

We give a formulation of the designed debutanizer flexibility index calculation problem. We present the algorithms of flexibility index calculation.

Текст научной работы на тему «Анализ гибкости оптимального проекта колонны дебутанизации»

Н. Н. Зиятдинов, Т. В. Лаптева, Н. Ю. Богула

АНАЛИЗ ГИБКОСТИ ОПТИМАЛЬНОГО ПРОЕКТА КОЛОННЫ ДЕБУТАНИЗАЦИИ

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

гибкости.

Дается постановка и решение задачи определения индекса гибкости для спроектированной колонны дебутанизации. Приводятся алгоритмы вычисления индекса гибкости.

Keywords: distillation column, optimal design, performance, flexibility index, flexibility function.

We give a formulation of the designed debutanizer flexibility index calculation problem. We present the algorithms of flexibility index calculation.

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

системы сохранять работоспособность и найти оптимальную конструкцию, гарантирующую

сохранение работоспособности, несмотря на изменение внутренних и внешних факторов на этапе функционирования. Выделяя в жизни спроектированной системы этапы проектирования и функционирования, отметим, что на этапе

функционирования могут изменяться условия эксплуатации ректификационной колонны.

В [1] было приведено решение задачи оптимального проектирования ректификационной

колонны для выделения бутановой фракции, с использованием метода, предложенного в работах [23].

Проведем оценку работоспособности

спроектированной колонны. Будем предполагать, что изменения условий эксплуатации происходят в виде изменения состава сырья. Тогда, в число неопределенных параметров в: войдут массовые доли н-бутана (/ = 1), н-пентана (/ = 2), н-гексана (/ = 3), т.е. размерность вектора неопределенных параметров пв= 3.

Тогда определим вид области неопределенности как:

Г(5) = (вв < вN < в?,/ = 1,...,пв,|]в, = 1}

/=1

где в" = 0,35 , в2л = 0,3 , в" = 0,35 взяты из характеристик состава сырья задачи проектирования [1]. Параметры в[ = в" - д^,в? = в" + д? ,

/ = 1,...,пв , с учетом условия нормирования состава

£у = 1

определяют

размер

области

неопределенности.

Для оценки работоспособности

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

конструктивные параметры ректификационной

колонны, найденные при решении задачи проектирования включают: число и диаметр тарелок в укрепляющей т1,б1 и исчерпывающей секциях т2,12 колонны; межтарельчатое расстояние ЬГду ; размеры поверхностей теплообмена конденсатора ^сопй и ребойлера БгеЬг. Введем вектор конструктивных параметров колонны

V = {ті, т2 , <Зи с(2, Ь1ту, Бсопа, ЭГеЬг}.

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

на качество целевого продукта:

хС4н > 0,98 ,

С4Н10

на коэффициент захлебывания в колонне: кАО ^ 0,82 ,

на максимальную тепловую нагрузку конденсатора и кипятильника:

Q < Q max

Qcond < Qcond ■ Q < Q max Qrnbr < QtBbr ,

(3)

(4)

Qmax

r

В ограничениях (З), (4) Q™d, у«ГєЬг-предельные тепловые нагрузки для конденсатора и кипятильника, значения которых найдены в результате решения задачи оптимального проектирования колонны дебутанизации.

Задача определения индекса гибкости заключается в нахождении такого размера области неопределенности 8 = {8-,8U,i = 1,...,З}, при котором ректификационная колонна,

характеризуемая конструктивными параметрами v , остается работоспособной, то есть на всей области неопределенности T(8) выполняется требование гибкости:

max min max /, (v,u,0) < Q,

0gT (8) u j=1.4 ‘

где j = 1,...,4 - номера ограничений /j (v,u,0),

представленных неравенствами (1) - (4); u - режим работы колонны, определяется флегмовым числом R и температурой ТгєЬг в кипятильнике колонны,

u = {R,ТгєЬг }.

I=1

Тогда задача определения индекса гибкости ректификационной колонны примет вид:

F (v) = niax X (8L + 8),

sL.sU

i=1

max min max /,(v,u,e) < 0,

вєГ (8) u j=1.4 ‘

ne

где T(8) = в і ві < вN < вU,i = 1,...,Пу,Xв/ = 1} -

i=1

область неопределенности.

Задача (З) представляет собой задачу многоэкстремальной недифференцируемой

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

Для решения задачи (З) будем использовать двухуровневый подход, где на верхнем уровне будет проводиться изменение поисковых переменных 8l ,8u , i = 1,...,З, а на нижнем уровне будет вычисляться значение функции гибкости r(v) = max min max w, (v,u,e)

вєГ (8) u j=1.4 ‘

Тогда задача (З) примет вид:

З

F(v) = max X (6L + 6U)

5iL,5U

i=1

z(v) < 0 .

Для решения задачи (7) на верхнем уровне можно использовать любой поисковый метод, в том числе метод сканирования.

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

Вычисление функции гибкости (6) представляет собой задачу многоэкстремальной недифференцируемой оптимизации, для решения которой неприменимы методы решения задач нелинейного программирования. Поэтому будем использовать для нахождения значения функции ) метод разбиения и границ, предложенный в [4, 5].

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

Верхняя оценка, предложенная в [4, 5], для данной задачи на к -ом шаге метода примет вид:

X u,k(v) = min w

ul,w

max у i(v,ul,9) < w, l = 1,..., Lk, j = 1,...,4,

вегк

где T/k есть подобласть, принадлежащая разбиению

*"К

T (8) = UTi1 , полученному к k -му шагу

метода,

ne

Тк = {в, : в-/к < вN < виік,і = 1,...,пв,^в; = 1}

і=1

. В [5] показано, что %(у) > %и,к (V).

Задача (9) является задачей полубесконечного программирования, для ее решения будем применять метод внешней аппроксимации. На каждой итерации

метода задача (9) будет заменяться на нижнюю оценку решаемой задачи. Для ее формирования

заменим непрерывные области Т)к дискретными к (5)

множествами 5; . В [5] показано, задача нижней оценки задачи (9) имеет вид: (6)

%іик (V) = тіп ю

и ,ю

уі (у,и/ в) < ю, Vв{ч є Бк , / = 1,..., -к,

У = 1,...,4, д = 1,...,Ок,

кк где О/ - число точек в множестве О/ ,

jq

= в і £в, = 1}.

і=1

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Рассмотрим построение множества 3^. Пусть {и*,ю*} есть решение задачи (10). Тогда множество 3^ = {в^ :уі (у,и/ ,в ^ю *}. Точки

ввд будем называть критическими точками.

Алгоритм вычисления верхней оценки на к -ом шаге метода имеет следующий вид [5].

Алгоритм 1. (7)

Шаг 1. Задание стартовых значений: номер итерации г = 1, разбиение области неопределенности Т/к; Я/ = 0 ; Б,0 = 0 ,

/ = 1,...,-к ; иГ, ю0. к - номер итерации внешнего метода разбиения и границ.

Шаг 2. Построение множеств Я/,

/ = 1,...,-к, критических точек решением 4-к задач:

Т j = max /j (v,u;r-1,в), j = 1,...,4.

I = 1,..., Lk

В результате получаем решение {Т„; ej}.

Если Т. - wr1 > 0, I = 1,...,Lk , j = 1,...,4,

то соответствующая точка множество R .

ej

заносится в

Lk

Шаг 3. Если множество и R/ =0, то

=1

переходим на шаг 1, иначе

Шаг 4. Построение множеств Б к по

правилу 3; = 3/к-1 и Я/ , / = 1,...,-к . (9)

Шаг 5. Решение задачи (10). Пусть {х-и'г ,иг} - решение задачи (10).

Присвоение юг = х-и,г (у).

Шаг 6. г = г +1. Перейти на шаг 2.

Шаг 7. Построить множество номеров областей Vй , в которых присутствуют активные

ограничения Vй = {/: Зу(і = 1,...,4): Ту - юг-1 > 0}.

В работе [5] показано, что итерационная процедура алгоритма действительно сходится к

верхней оценки функции гибкости Хи (V).

i=1

В качестве нижней оценки функции гибкости используется решение задачи:

XL'k (v) = min w

uql ,w

/j (v,uq ,eq ) < w

V0q £ Sk, j = 1,...,4, q = 1...Qk ,

где Sk - множество точек eq £ T(S),

Пв

в4 = {в, : ^в; = 1}, Gk - количество точек в

i=1

множестве Sk. В [5] показано, что x(v) > XL,k (v).

Поскольку множество Sk содержит конечное количество точек в, задача (11) является обычной задачей нелинейного программирования, для ее решения можно использовать хорошо зарекомендовавшие себя методы нелинейной оптимизации, например, метод последовательного квадратичного программирования.

Для реализации алгоритма разбиения и границ используется разбиение области неопределенности на подобласти. Разбиению подвергаются не все подобласти. На к -ом шаге выбираются для разбиения такие подобласти T f для которых на последней итерации алгоритма вычисления верхней оценки хотя бы одно ограничение /j (v,u[ ,в) стало активным. То

есть, подобласть подвергается дроблению, если выполняется следующее условие:

3j £ J max у j(v,urt,B) = wr

9eTk

где urj, wr - решение задачи (9), полученное на шаге 5 алгоритма вычисления верхней оценки.

Разбиение областей, для которых условие (6) выполняется, приведёт, вообще говоря, к удалению соответствующих активных ограничений, что в свою очередь улучшит значение верхней оценки функции гибкости.

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

/,k(v) > 0,

ZU'k(v) < 0

Выполнение условия (13) означает, что найдена такая точка в , для которой невозможно подобрать управление u , чтобы нижняя оценка функции гибкости была неположительной. Поскольку верно ^(v) >xL (v) , значение самой функции гибкости будет заведомо неотрицательным и рассматриваемая конструкция v будет неработоспособной на рассматриваемой области неопределенности T(S).

При выполнении условия (14), поскольку х(м) >хиу (V), конструкция V является работоспособной на рассматриваемой облШск неопределенности Т(д).

Алгоритм метода разбиения и границ примет следующий вид.

Алгоритм 2.

Шаг 1. Положим к = 1. Выбираем

начальное разбиение области Т(д) на подобласти

ТI , начальные значения и°, V/0, I = 1,...,Ьк, малые

величины 81 < 0 , 82 < 0 .

Шаг 2. Вычисление верхней оценки

функции гибкости, решением задачи (9), используя

алгоритм 3.1. Пусть (хш,г;иг) - решение этой

задачи, Ук - множество активных областей на

последней итерации алгоритма 3.1.

Шаг 3. Вычисление нижней оценки

функции гибкости, решением задачи (13). В качестве множества точек Эк используется набор критических точек, полученных при решении задачи вычисления верхней оценки функции гибкости. _ Ц ^к

Э,к = иэ,к, вк = у Ок .

1=1 ы

Шаг 4. Если выполняется одно из условий

(13), (14), то решение найдено. При выполнении

условия (13) х(У) = Х1~’к&), при выполнении

условия (14) х^) = хи,к^). Останов. (12)

Если условия (13), (14) не выполняются, то проверить условие:

|хи,к^) -ХЬ,к^)| <81 .

Если условие (15) выполняется, то решение найдено, XV) = хи,к (V), останов.

Шаг 5. Если выполняется условие: г(Тк)<82, I = 1,..., 1к, утк еУк

то значение х(У) = хик^), останов. В противном случае, заносим область, для которой нарушается условие (16), в Ук .

Шаг 6. Если выполняется условие Ук = 0 ,

то х^У) = хи,к (V), останов.

Шаг 7. Разбиение подобластей. Разбить каждую подобласть Т к е Ук на две подобЩЪти

так, что Т к = Т/1 и Т/2 . Образовать (нйвое

множество подобластей Тк+1 из старого множества Тк, заменяя подобласть Т]к новыми подобластями Тк и Т/2.

Шаг 8. Положить к = к +1 и перейти к

шагу 2.

На шаге 7 подобласть Т ]к е Ук , I = 1,..., 1-к разбивается на две подобласти Т/к1 и Т/к2 :

Тк = {в : в є Тк в < с3}, Т2 = {в : в є Т/к в < с3}, где в3 - переменная разбиения, с3 - точка разбиения.

Верхний уровень алгоритма решения задачи (15) заключается в определении предельно допустимых значений переменных

8-,тах,8и,тах,і = 1,...,3,, задающих размер области неопределенности

пв

Т(8) = {в, : в- - 8-,тах < вN < ви + 51и тах, і = 1,пв, ^ в = 1} , на

і=1

которой спроектированная колонна остается работоспособной. Определение значений параметров 8-,тах, 8и,тах, і = 1,...,3, проводилось

последовательным увеличением значений параметров 8-,8и, і = 1,...,3, на постоянную величину.

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

Результаты

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

работоспособности формируется из индивидуальных диапазонов по следующему правилу Таблица 1 - Индивидуальные диапазоны

Правые границы индивидуальных диапазонов на массовые доли н-бутана и н-пентана были определены нарушением ограничения (4) на фактор захлебывания в колонне, т. е. при превышении указанных правых границ области неопределенности было выявлено захлебывание в колонне.

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

Литература

1. Н.Н. Зиятдинов. Н.Ю. Богула, Т.В. Лаптева, Г.М. Островский, Вестник Казанского технологического университета, 5, 118-123 (2011).

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

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Богула, Доклады Академии Наук, 431, 6, 768-771

(2010).

3. Г.М. Островский, Н. Н. Зиятдинов, Т. В. Лаптева, Н.Ю. Богула, Теоретические основы химической технологии, 45, 1, 88-97 (2010).

4. Г.М. Островский, Н.Н. Зиятдинов, Т.В. Лаптева. Оптимизация технических систем, КНОРУС, Москва, 2012. 432 с.

5. Первухин И.Д. Автореф. дисс. канд. тех. наук, КГТУ, Казань, 2011.

изменения неопределенных параметров

Неопределенный параметр вN - 8-,тах вN +8и,тах

Массовая доля н-бутана 0,11 0,87

Массовая доля н-пентана 0,05 0,75

Массовая доля н-гексана 0,05 0,65

© Н. Н. Зиятдинов - д-р техн. наук, проф., зав. каф. системотехники КНИТУ, nnziat@yandex.ru; Т. В. Лаптева - канд. техн. наук, доц. той же кафедры; Н. Ю. Богула - канд. техн. наук, асс. той же кафедры, nellybog@gmail.com.

i Надоели баннеры? Вы всегда можете отключить рекламу.