Глава 3. АЛГОРИТМ МОДЕЛИРОВАНИЯ КОМПЛЕКСНОГО

Глава 3. АЛГОРИТМ МОДЕЛИРОВАНИЯ КОМПЛЕКСНОГО

ПРОИЗВОДСТВА И ЕГО РЕАЛИЗАЦИЯ НА ЭВМ.

 

Приведенные ниже алгоритмы расчета позволяют реализовать на ЭВМ методику математического моделирования, основные принципы которой описаны в Главе 2.

 

3.1. Блок переработки природного газа.

 

3.1.1. Риформинг метана.

 

В аппарат первичного риформинга входит два материальных потока: пар и природный газ (см. рис. 8). Расход и состав образующейся газовой смеси определяется покомпонентным сложением с учетом формулы (2.29). Температура на выходе – регулируемый параметр, поэтому ее значение задается до начала расчета.

Из [74] известно, что в первичном риформинге протекают две независимые равновесные реакции (2.13) и (2.14). Они описываются системой из двух уравнений типа (2.27), которая решается итерационным методом [19]. Пусть a – количество СН₄, прореагировавшего по реакции (2.13), b – количество СО прореагировавшего по реакции (2.14) (в нм³/ч). Для первого шага итерации были заданы некоторые начальные значения a и b. Тогда, считая, что реакции (2.15) протекают до конца, операциями покомпонентного сложения и вычитания рассчитали состав равновесной газовой смеси. Затем, вычислив константы равновесия реакций (Кр₁ и Кр₂) при температуре Т_р (формулы (7,8) Приложения III), записали систему уравнений в следующем виде (с учетом формул (2.28-2.29):

 

(3.1)

(3.2)

где V_i – объемы компонентов в равновесной газовой смеси, нм³; V_вых – общий объем равновесной газовой смеси, нм³.

Таким образом, a и b есть решения системы уравнений (3.2-3.3), если c = 1и d = 1. В противном случае необходимо задать новые значения a = a + Da и b = b + Db. Алгоритм моделирования представлен на рис. 15.

Здесь и далее на блок-схемах приняты следующие обозначения:

- прямоугольник – входные и выходные данные или рассчитываемые параметры;

- штриховой прямоугольник – принимаемые значения для итерационных вычислений;

- ромб – условие;

- двойная стрелка - выполнение математических операций в соответствии с формулами, отмеченными над стрелкой;

- простая стрелка – передача начальных или расчетных данных.

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

 

 

 

Рис. 15. Блок-схема моделирования стадии первичного риформинга.

 

В аппарат вторичного риформинга подаются паро-газовая смесь (ПГС) из первичного риформинга и паро-воздушная смесь (ПВС). В I зоне (над катализатором) происходит полное сгорание кислорода с водородом по реакции (2.16).

Расход и состав газовой смеси до и после 1-й зоны определялся покомпонентным сложением и вычитанием согласно стехиометрии реакции (2.16). Во 2-й зоне на катализаторе протекают реакции (2.13,2.14), описываемые системой уравнений (3.1,3.2). Так как реактор адиабатический, то температура на выходе до расчета неизвестна, поэтому необходимо добавить уравнение теплового баланса (2.63), записанное в виде:

Q₁ + Q₂ + Q₃ + Q₄ - (Q₅ + Q₆ +Q₇) = f (3.3)

Тепловые эффекты реакций (2.13,2.14) известны [74]. При определении энтальпий газовых смесей по формуле (2.9) использовали зависимости H = f(t), полученные для каждого компонента по уравнению (2.2).

Пусть e – температура на выходе из реактора, °С. На первом шаге итерации были заданы некоторые значения a,b,e. Таким образом, для стадии вторичного риформинга a,b,e есть решения системы уравнений (3.1-3.3), если c = 1, d = 1 и f = 0. В противном случае необходимо задать новые значения a = a + Da, b = b + Db и e = e + De. Алгоритм моделирования представлен на рис. 16.

 

 

Рис. 16. Блок-схема моделирования стадии вторичного риформинга.

 

3.1.2. Конверсия СО и очистка газа.

 

Конверсию СО по реакции (2.14) проводят в 2 стадии, отличающиеся температурой процесса [81]. На обеих стадиях применяют адиабатические реакторы [74]. Поэтому для определения состава и температуры газа на выходе необходимо решить систему из двух уравнений. Одно из них - уравнение (3.2), второе – уравнение теплового баланса:

Q₁ + Q₂ - (Q₃ + Q₄) = f (3.4)

Итерационные вычисления организованы аналогично стадиям первичного и вторичного риформингов.

Алгоритм моделирования одинаков для каждой стадии конверсии СО и приведен на рис. 17.

Рис. 17. Блок-схема моделирования стадии конверсии СО.

 

За равновесную температуру абсорбции приняли температуру глубокорегенерированного раствора МЭА, подаваемого на верх абсорбера. Равновесное давление СО₂ над раствором определяем по практическим данным [74]. Реакции (2.19,2.20) протекают до конца, поэтому состав газа после метанирования рассчитывали покомпонентным сложением и вычитанием согласно стехиометрии реакций.

Алгоритм моделирования представлен на рис. 18.

Рис. 18. Блок-схема моделирования стадии очистки.

 

3.2. Блок выработки пара высокого давления.

 

3.2.1. Котел-утилизатор.

При расчете котла-утилизатора из уравнения теплового баланса (2.64) определяли количество теплоты, затраченное на испарение воды. Величины Q1 и Q2 вычисляли по формуле (2.9), зависимости энтальпий компонентов газовой смеси от температуры – по уравнению (2.2). Количество образовавшегося пара находили по формуле (2.65). Зависимость Нисп = f(t) в виде полинома получили аппроксимацией литературных данных из [58]. Так как заданной величиной является давление вырабатываемого пара, то температуру насыщения находили из уравнения tнас = f(Рнас), также полученного путем аппроксимации литературных данных из [58].

Алгоритм расчета котла-утилизатора не требует итерационных процедур и приведен на рис. 19.

Рис. 19. Блок-схема моделирования котла-утилизатора.

 

3.2.2. Вспомогательный котел.

 

Количество пара, вырабатываемого во вспомогательном котле, определяли из баланса выработки и потребления пара высокого давления. Из теплового баланса котла рассчитывали необходимый расход природного газа. Для этого уравнение (2.66) записали в виде:

Q₁ + Q₂ + Q₃ - (Q₄ + Q₅ + Q₆) = f (3.5)

Расчет по уравнению (3.5) вели методом итераций. Пусть е – расход природного газа во вспомкотел. Далее исходя из стехиометрии реакций (2.26) и принятого значения коэффициента избытка воздуха рассчитывали материальный баланс процесса горения: количество воздуха, необходимое для горения, и количество и состав дымовых газов. Затем определяли величины, входящие в уравнение (3.5) и проверяли сходимость теплового баланса. Если f = 0 с заданной точностью, то е есть решение уравнения (3.5). Алгоритм расчета показан на рис. 20.

Рис. 20. Блок-схема моделирования вспомогательного котла.

 

3.3. Блок синтеза аммиака.

 

В соответствии с рекомендациями из [66,74] мы предположили, что после стадии очистки в газовой смеси присутствуют только H₂, N₂, CH₄, Ar, т.к. концентрации СО, СО₂ и Н₂О пренебрежимо малы.

Объем и состав потока 0 (см. рис. 14) следуют из расчета стадии очистки. Состав потока 0₁ неизвестен, поэтому приняли начальные значения расходов каждого из компонентов (i = 1,2..5) в нм³/ч. Состав потока 1 рассчитали покомпонентным сложением. Далее по уравнениям (15-18) Приложения III определили концентрацию аммиака в газовой смеси после вторичной конденсации и, тем самым, потоки 2,3.

Пусть a – количество NН₃, прореагировавшего по реакции (2.21). Для первого шага итерации задали некоторое значение a и вычислили константу равновесия К_р (уравнение (14) Приложения III) при температуре Т_р (формула (2.30)), записали уравнение (2.27) в следующем виде (с учетом формул (2.28-2.29):

(3.6)

где V_i – объемы компонентов в равновесной газовой смеси, нм³; V_вых – общий объем равновесной газовой смеси, нм³; Р – давление, МПа.

Если c = 1, то a есть решение уравнения (3.6) и определен поток 4.

Затем по уравнениям Милельса вычислили концентрацию аммиака в газовой смеси после первичной конденсации и потоки 5,6. После этого отвели поток 7, расход которого в модели может изменятся в заданных пределах. Потоки 9,10 рассчитали также по уравнениям Михельса и далее поток 11. Состав потока 8 сравнивали с потоком 0₁ :

e_i = N_01,i/N_8,i (3.7)

где N – расход компонента в нм³/ч; i = 5 для рассматриваемой смеси.

Если e_i ¹ 1, то необходимо задать новые значения .

Массовый расход аммиака определяли по формуле (2.39). Потери аммиака и потребности аммиака на технологические нужды приняли равными 5 % от производительности.

Алгоритм моделирования стадии синтеза аммиака представлен на рис. 21.

Рис. 21. Блок-схема моделирования стадии синтеза аммиака.

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

 

3.4. Блок синтеза метанола.

 

3.4.1. Реактор синтеза метанола.

 

Основным аппаратом процесса синтеза метанола является реактор. Поэтому описание алгоритма моделирования начнем с этого аппарата.

При моделировании реактора синтеза метанола необходимо знать некоторые его конструкционные характеристики. Для этого приняли внутренний диаметр реакционный труб равным 120 мм, материал труб - 12Х18Н12Т в соответствии с рекомендациям [43]. Далее по формулам из [48] рассчитали: толщину труб из условия прочности, шаг между трубами, число труб, внутренний диаметр корпуса. После этого из стандартных рядов [48] выбрали диаметр корпуса и рассчитали его толщину из условия прочности. Описанный расчет выполнен в ''Microsoft Excel'' в виде блока и здесь не приводится. Основные конструкционные характеристики аппарата приведены ниже:

диаметр корпуса 4148 ´ 78 мм

диаметр труб 132 ´ 6 мм

число труб 397

Кроме того, необходимо задать размеры таблеток катализатора и его насыпную плотность. Согласно [8,14,82] в синтезе метанола применяют медьсодержащие катализаторы с размером таблеток (диаметр ´ высота): 3 ´ 3 мм и 5 ´ 5 мм, 4 ´ 6 мм. Для определенности приняли таблетки 4 ´ 6 мм. Из литературы [8, 82] известна насыпная плотность катализатора – 1300 кг/м³.

Прежде, чем начать расчет по слою, были выполнены некоторые предварительные расчеты.

Поверхность F_з (м²) и объем V_з (м³) зерна указанной формы определили по формулам:

F_з = p×d_з×h_з + 2×p×d_з²/4 (3.8)

V_з = p×d_з² ×h_з /4 (3.9)

Порозность слоя зерен чаще всего лежит в пределах 0,35-0,45 [31]. Приняли e = 0,4. Удельная поверхность слоя катализатора a (м²/м³) и эквивалентный диаметр d_э (м) каналов в слое [31]:

a = (1 - e)× F_з/ V_з (3.10)

d_э = 4×e/a (3.11)

Объемный расход V_вх (м³/ч) газовой смеси при условиях входа в реактор:

V_вх = V_вх,0 ×Т_вх / (273×Р_вх) (3.12)

где V_вх,0 – объемный расход в нм³/ч; Р_вх – давление в атм.

Фиктивная средняя линейная скорость w₀ (м/с) потока:

(3.13)

Действительная скорость w (м/с) потока в каналах слоя:

w = w₀/e (3.14)

Далее вычислили расход N_вх (моль/с) и состав n_вх,i (% мол.) газа на входе в одну трубу:

(3.15)

n_вх,i = V_вх,i

где n_вх,i, V_вх,i – мольная и объемная доли компонентов в газовой смеси.

При моделировании реактора необходимо рассчитывать вязкость и теплопроводность газовой смеси, содержащей 8 компонентов, с учетом влияния давления. Выполнение такого расчета на каждом шаге значительно усложняет модель, поэтому мы воспользовались следующим приемом. Приняли температуру на выходе из реактора 265 °С и рассчитали равновесный состав газа на выходе, а затем покомпонентно усреднили состав газа между входом и выходом. Далее для полученного среднего состава газа в пакете MathCAD 2000 по формулам (1-6 Приложения III) рассчитали вязкость и теплопроводность смеси при давлении 70 атм в диапазоне температур 180 – 300 °С (с шагом 20 °С). Полученные расчетные данные аппроксимировали полиномами 1-й степени: m_см = f(t) и l_см = f(t).

Задав шаг по высоте Dl (м), определили геометрические параметры:

F_ст,вн = p× d_тр,в×Dl, м² (3.16)

F_ст,н = p× d_тр,н×Dl, м² (3.17)

, м³ (3.18)

Материальный баланс в слое рассчитывали следующим образом. Для определения скоростей реакций вычисляли парциальные давления реагентов Р_i (атм):

Р_i = n_вх,i×Р_вх (3.19)

Далее:

Dn₁ = r_i×Dv×r_нас (3.20)

Dn₂ = r₂×Dv×10⁶ (3.21)

Далее рассчитывали мольные расходы реагентов (моль/с) на выходе из слоя:

(3.22)

(3.23)

(3.24)

(3.25)

(3.26)

Состав газа (% мол.) на выходе из слоя:

(3.27)

(3.28)

Зависимость давления насыщенного водяного пара Р_нп (атм) от температуры кипения получили аппроксимацией табличных данных из [57]:

Р_нп = 94,2972 – 1,0908×t_в + 3,4893×10^-3×t_в² (3.29)

Зависимость коэффициента теплопроводности материала труб l_ст (Вт/(м×К)) от температуры получили аппроксимацией табличных данных из [43]:

l_ст = 15,00 + 0,011×t + 6.526×10^-6×t² (3.30)

Из уравнений (2.60,2.61) получили выражение для Т_ст,вн:

(3.31)

где B - параметр, определяемый уравнением:

(3.32)

Температуру на выходе из слоя рассчитывали по формуле:

(3.33)

Слой катализатора оказывает гидравлическое сопротивление DР (атм), которое рассчитывали по формуле [48]:

(3.34)

(3.35)

где l_сопр - коэффициент сопротивления. Тогда давление на выходе из слоя:

Р_вых = Р_вх - DР (3.36)

Алгоритм моделирования реактора синтеза приведен на рис. 22.

 

 

Рис. 22. Блок-схема моделирования реактора синтеза метанола.

 

3.4.2. Конденсация метанола и воды.

 

В контуре метанола (см. рис. 8) из газовой смеси за счет охлаждения конденсируется вода. Алгоритм математического описания этой стадии представлен на рис. 23.

 

Рис. 23. Блок-схема моделирования стадии конденсации воды.

 

При конденсации метанола-сырца (см. рис. 8) конденсируется смесь из 2-х компонентов (метанол и вода), поэтому применили итерационный метод расчета. Пусть a и b – количества сконденсировавшихся метанола и воды соответственно (в нм³/ч газа). Для первого шага итерации задали некоторые значения a и b. Тогда (индексы 1 и 2 соответствуют метанолу и воде):

(3.37)

(3.38)

Учитывая, что количества сконденсировавшихся метанола и воды в кмолях равны a/22,4 и b/22,4, получили:

(3.39)

(3.40)

Далее рассчитали коэффициенты активности и записали уравнение (2.31) для метанола и воды в виде:

(3.41)

(3.42)

Таким образом, a и b есть решения системы уравнений (3.41-3.42), если c = 1и d = 1. В противном случае необходимо задать новые значения a = a + Da и b = b + Db. Алгоритм моделирования представлен на рис. 24.

Рис. 24. Блок-схема моделирования конденсации метанола и воды.

 

3.4.3. Контур синтеза метанола.

 

Объем и состав потока 0 (см. рис. 13) были определены из расчета стадии конденсации воды. Приняв для определенности K = 3 в формуле (2.76), находили объем потока 0₁. Состав потока 0₁ неизвестен, поэтому приняли некоторые значения концентрации компонентов N_01,i (в % об.) в потоке 0₁. Состав потока 1 вычисляли покомпонентным сложением. Далее по приведенным выше алгоритмам рассчитывали стадии синтеза и конденсации (потоки 2,3,4). Затем отводили точный объем циркулирующего газа согласно принятой кратности циркуляции, определяя тем самым потоки 5 и 6. В случае несовпадения концентрации компонентов в потоках 6 и 0₁ задавали новый состав потока 0₁ и повторяли расчет до достижения совпадения с заданной точностью e:

e_i = N_01,i/N_6,i (3.43)

где N – объемная доля (в %) компонента в газовой смеси; i = 8 для рассматриваемой смеси.

Массовые расходы метанола и воды в составе метанола-сырца определяли по формуле (2.39). Алгоритм моделирования контура представлен на рис. 25.

Рис. 25. Блок-схема моделирования контура метанола.

 

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

 

3.4.4. Теплообменное оборудование.

 

Тепловой баланс контура определяется суммированием балансов входящих в него теплообменных аппаратов и потоков между ними.

На рис. 26 представлена общая блок-схема моделирования теплообменников. Зависимость свойств теплоносителей от температуры мы описали полиномами различных степеней по справочным данным [43,51,56,74,82] с помощью программы ''MathCAD''.

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

Неизвестный параметр е - это температура на выходе или расход одного из теплоносителей. Задав первое приближение для е, рассчитывали материальный баланс конденсации по алгоритмам, описанным ранее (см. рис. 23,24). Далее вычисляли приход и расход тепла и проверяли сходимость теплового баланса (f = 0). При выполнении условия определяли тепловую нагрузку и среднюю разность температур.

После этого из справочников и каталогов [37,54,78,84] выбирали стандартный теплообменник и проводили проверочный расчет на прочность по методике из [48]. Затем по описанным выше уравнениям (см. раздел 2.4.3) вычисляли коэффициенты теплоотдачи и теплопередачи, необходимую поверхность теплообмена и запас поверхности.

 

Рис. 26. Блок-схема моделирования теплообменника.

3.5. Общий алгоритм моделирования и его реализация

на ЭВМ.

 

При построении общего алгоритма моделирования важно правильно организовать передачу данных между отдельными расчетными блоками, т.к. выходные данные одних блоков являются входными для других. Порядок передачи основных данных (расходы, давления, температуры) показан стрелками на структурной схеме построения модели ХТЭС ''по горизонтали'' (см. рис. 10).

Исходные данные для модели совместного производства (регулируемые параметры модели) приняты согласно [1,13,18,81,82]. Список исходных данных с их конкретными значениями приведен в Приложении I.

Результатами моделирования является набор расчетных данных (контролируемые параметры модели), перечень которых без конкретных значений также приведен в Приложении I. Полученные численные значения и зависимости между ними обсуждаются в Главе 4.

Математическая модель совместного производства реализована на ЭВМ в программе ''Microsoft Excel'' в виде отдельной ''книги'' [79]. Каждый алгоритм выполнен на отдельном ''листе'' ''книги'', все ''листы'' взаимосвязаны между собой. Выделены ''листы'' исходных данных и результатов моделирования.

Итерационные вычисления организованы в надстройке ''Поиск решения'', где поиск нужных значений переменных осуществляется методом Ньютона [79]. Относительная погрешность вычислений задается и определена нами величиной 0,001 %.

Отметим, что в контуре синтеза метанола (см. рис. 25) итеративный расчет проводится по 8 параметрам. На каждом шаге итерации выполняется довольно объемный расчет синтеза метанола по слою и расчет конденсации (в свою очередь тоже итеративный). Тем не менее, благодаря четко составленному алгоритму, надстройка ''Поиск решения'' справляется с поставленной задачей.

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

Дело в том, что многие пользователи, в том числе и автор диссертации, искренне желая применить компьютерное моделирование в своей практической деятельности, сталкиваются с серьезными трудностями при освоении и использовании современных программных средств [39]. Для работы с ними все еще требуются знания, не относящиеся непосредственно к моделированию, а проведение вычислительного эксперимента остается кропотливой и многотрудной работой. Опыт нашей работы показал, что намного проще и быстрее реализовать разработанные алгоритмы с помощью стандартных, совместимых систем типа ''Microsoft Excel'' и ''MathCAD'' , и на их основе создавать уникальные системы расчета и оптимизации ХТЭС.

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

 

Выводы по главе:

 

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

2. Показано, что модель достаточно просто реализуется в ''Microsoft Excel'' без необходимости разработки специальных программ. Особое внимание уделено организации итерационных вычислений, что позволяет рассчитывать материальные и тепловые балансы в достаточно сложных технологических системах с многочисленными рециклами.