Акулов Л.А. Установки и системы низкотемпературной техники. Ожижение природного газа и утилизация холода сжиженного природного газа при его регазификации - файл n3.doc

Акулов Л.А. Установки и системы низкотемпературной техники. Ожижение природного газа и утилизация холода сжиженного природного газа при его регазификации
скачать (4119.2 kb.)
Доступные файлы (6):
n1.doc388kb.04.12.2006 15:29скачать
n2.doc62kb.29.11.2006 19:54скачать
n3.doc2333kb.10.11.2006 00:50скачать
n4.doc867kb.13.11.2006 04:02скачать
n5.doc773kb.23.05.2007 14:54скачать
n6.doc3194kb.10.11.2006 03:16скачать

n3.doc

  1   2   3   4   5


Федеральное агентство по образованию

Государственное образовательное учреждение

высшего профессионального образования

САНКТ-ПЕТЕРБУРГСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ


НИЗКОТЕМПЕРАТУРНЫХ И ПИЩЕВЫХ ТЕХНОЛОГИЙ




Учебное пособие

Рекомендовано Учебно-методическим объединением
по университетскому политехническому образованию
в качестве учебного пособия для студентов
высших учебных заведений, обучающихся
по направлению подготовки 140400 «Техническая физика»


Санкт-Петербург 2006


УДК 536.483, 621.59

ББК 31.392
Б 82


Борзенко Е.И., Зайцев А.В.
Установки и системы низкотемпературной техники. Автоматизированный расчет и моделирование процессов криогенных установок и систем: Учеб. пособие. – СПб.: СПбГУНиПТ, 2006. – 232 с.


ISBN 5-89565-146-1
Содержит материал, включенный в программу курсов по направлению «Техническая физика» и специальности «Техника и физика низких температур». Рассматриваются алгоритмы расчета термических, калорических и транспортных свойств чистых криогенных рабочих веществ, реальных веществ и их смесей, алгоритмы расчета фазового равновесия многокомпонентных смесей. Приводятся математические модели контактных устройств и колонных аппаратов узлов ректификации криогенных смесей. Анализируются статистические характеристики ректификационных колонн.

Может быть использовано для разработки программного обеспечения курсовых и дипломных проектов студентов, квалификационных работ бакалавров и диссертаций магистрантов, а также полезно для инженеров и аспирантов, занимающихся моделированием криогенных процессов на ЭВМ.
УДК 536.483, 621.59

ББК 31.392
Рецензенты

Кафедра Э-4 Московского государственного технического университета им. Н.Э. Баумана (зав. кафедрой доктор техн. наук, проф. А.М. Архаров).

Доктор техн. наук, проф. М.М. Пеньков (Военно-космическая академия им. А.Ф. Можайского).
Рекомендовано к изданию редакционно-издательским советом университета
ISBN 5-89565-146-1
© Санкт-Петербургский государственный университет низкотемпературных

и пищевых технологий, 2006

75-летию Санкт-Петербургского

государственного университета

низкотемпературных и пищевых технологий посвящается



ПРЕДИСЛОВИЕ

Развитие криогенной техники и технологии способствует даль­нейшему прогрессу в ряде новейших отраслей науки и техники. Сов­ременные газоразделительные установки и криогенные системы отно­сятся к наукоемким объектам, поэтому на стадии их разработки необ­ходимо применять новейшие методы проектирования и конструирования.

При анализе существующих и вновь разрабатываемых криогенных установок и систем в настоящее время широко применяются системы автоматизированного проектирования (САПР).

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

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

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

1. РАСЧЕТ ПАРАМЕТРОВ СОСТОЯНИЯ

РЕАЛЬНЫХ КРИОГЕННЫХ РАБОЧИХ ВЕЩЕСТВ
1.1. Уравнение состояния реальных веществ
В криогенных системах и установках в качестве рабочих веществ наиболее широкое применение получили гелий, водород, воздух и такие его основные компоненты, как азот, кислород, аргон, неон, криптон, ксенон, и их смеси. Особенностью криогенных систем является то, что температура и давление рабочих веществ изменяются в широких пределах и для описания их свойств, особенно в областях низких температур и высоких давлений, необходимо учитывать характер межмолекулярных взаимодействий, т. е. применять уравнения состояния реального газа.

Простое эмпирическое уравнение состояния реального газа, полученное с учетом сжимаемости 1 кг рабочего вещества, имеет вид
, (1.1)
где p – давление; v – удельный объем; z – коэффициент сжимаемости; R – газовая постоянная; T – температура.

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

При обобщении экспериментальных данных часто используют вириальное уравнение Боголюбова–Майера [1], которое базируется на наиболее строгом теоретическом подходе и имеет вид

или в сокращенной форме
, (1.2)
где – вириальные коэффициенты, зависящие от температуры; – плотность.

Число вириальных коэффициентов в уравнении (1.2) выбирают, исходя из требований к точности обобщения экспериментальных данных; тогда уравнение принимает вид
. (1.3)
Каждый из вириальных коэффициентов, в свою очередь, записывают как
,
где si – границы изменения j; bij – коэффициенты разложения; – приведенная температура, здесь Tкр – критическая температура.

Часто уравнение Боголюбова–Майера (1.3) представляют в виде [2]
. (1.4)
По мере накопления экспериментальных данных по термодинамическим свойствам рабочих веществ теплоэнергетических установок и систем в научно-технической литературе появилось множество уравнений состояния, описывающих с определенной полнотой и точностью параметры их состояния в pvT-пространстве. Наряду с уравнением Боголюбова–Майера получили применение уравнения Битте–Бриджмана, Бенедикта–Рубина [1], а в последние годы – уравнения Рабиновича, Сычева, Вассермана и др. [3–6], взаимосогласованные уравнения Клецкого, уравнения состояния гелия-4 Тарана, Мак-Карти и т. д.

1.2. Расчет параметров состояния

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

Единое уравнение состояния криогенных рабочих веществ. При разработке системы расчета термодинамических и теплофизических свойств криогенных рабочих веществ в основу было положено уравнение состояния Сычева–Вассермана, которое имеет единую форму для основных криоагентов и обеспечивает относительно высокую точность расчета pvT-данных [7].

При аппроксимации экспериментальных данных по основным криогенным веществам было получено уравнение состояния (1.4) в следующем виде:
, (1.5)
где – приведенная плотность, здесь ?кр – критическая плотность.

Для повышения точности описания pvT-данных водорода и гелия-4 используют два набора коэффициентов для определенных температурных зон. В связи с тем что характер изменения термических свойств гелия-4 в области температуры от ?-линии до 25 К значительно отличается от характера изменения их в диапазоне температур от 25 до 1500 К, одно уравнение действует в зоне температур от ?-линии до 25 К, а второе – при значениях T =-15…1500 К. Диапазон 15…1500 К частично охватывает каждое из двух уравнений, поэтому переход от одной группы коэффициентов разложения bij к другой осуществляется при температуре 20 К. Однако необходимо отметить, что уравнение состояния в вириальной форме (1.5) не отображает достоверно поведение термодинамических функций в области критической точки.

Данные о давлении насыщенных паров широко применяются при моделировании процессов в элементах криогенных систем. В связи с этим при разработке алгоритма расчета термодинамического поля рабочих веществ могут быть использованы интерполяционные уравнения [8], в дальнейшем приведенные к виду
, (1.6)
где fi – коэффициенты уравнения; – относительная температура.

Калорические уравнения для криогенных рабочих веществ. Калорические уравнения, связывающие энтальпию, энтропию и теплоемкость с давлением и температурой, могут быть получены на основе обобщения экспериментальных данных либо путем интегрирования pvT-зависимости. В работах [3–6, 9] используется второй путь и предлагаются формулы для определения следующих значений:

энтальпии
; (1.7)
энтропии
; (1.8)
изобарной теплоемкости
, (1.9)

где – энтальпия, энтропия и изохорная теплоемкость в идеально-газовом состоянии*.

Энтальпию h0 и энтропию s0 определяют из соотношений
; (1.10)
, (1.11)
где h00 и s00 – энтальпия и энтропия при температуре T0; – теплота сублимации при T1=10 К; – константа.

Значения констант уравнений (1.7)–(1.11) приводятся в таблице.

Изобарную теплоемкость в идеально-газовом состоянии рассчитывают с помощью полинома
, (1.12)
где ? j и ? j – коэффициенты полинома.

При аппроксимации значений изобарной теплоемкости в идеально-газовом состоянии для азота в полиноме (1.12) вместо ?T используется значение приведенной температуры .

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

Параметры криогенных веществ и константы уравнений (1.7)–(1.11) для T0 = 100 K [3–6, 9]


Рабочее

вещество

Газовая

постоянная R, Дж/(кг∙К)

Фактор

акцент-ричности ?

Критические параметры

Коэффициент сжимае-мости zкр

Т
0
еплота сублима-ции (при T = 0 К) h0, кДж/кг

Энтальпия h00 при T0, кДж/кг

Энтропия s00 при T0, кДж/(кг∙К)

К
0
онстанта s0, кДж/(кг∙К)

Tкр, К

pкр, МПа

?кр, кг/м3

N2

296,8

0,0400

126,2

3,398

313,1

0,291

247,6

103,60

5,6997

0

Ar

208,146

–0,0020

150,65

4,864

531,0

0,290

193,40*

238,60*

3,8673*

3,2369*

H2

4124,2

–0,2801

33,19

1,297

30,11

0,292

377,9937

210,081

30,64318

0

p-H2

4124,2

–0,2491

32,984

1,287

31,43

0,292

377,9937

210,081

30,64318

0

Воздух

287,1

–0,0094

132,5

3,760

316,56

0,291

253,4

100,00

20,0824

0

4He

2077,252

–0,3310

5,19

0,22746

69,64

0,3

14,7404

519,313

25,81691

0

O2

259,835

0,0213

154,581

5,107

436,2

0,292

275,542

90,66

5,4124

0

Kr

99,215

–0,0020

209,40

5,49

912,0

0,291

128,78*

163,54*

1,9539*

1,7319*

Xe

63,322

0,0020

289,74

5,82

1100,0

0,290

114,38*

146,93*

1,2897*

1,1988*

Ne

411,94

–0,0388

44,45

2,721

484,0

0,296

105,08*

120,48*

7,2323*

4,7797*


* Условия определения отмеченных величин – см. в литературе [2].

Подпрограммы расчета термодинамических и калорических пара­метров криогенных рабочих веществ. В программах систем автоматизированного проектирования тепловых процессов в элементах криогенных установок с использованием современных ЭВМ особое место отводится подсистеме, описывающей термические и калорические свойства рабочих веществ.

На базе уравнений (1.5)–(1.13) разработана единая подсистема расчета термодинамических и калорических параметров основных рабочих веществ [10, 11] криогенных установок и систем. Значения параметров этих веществ можно рассчитывать с помощью отдельных подпрограмм-функций, структурно входящих в данную подсистему.

Коэффициент сжимаемости z и давление p реального криогенного рабочего вещества по известной плотности ? и температуре T определяем с помощью подпрограмм, изображенных на рис. 1.1 и 1.2.


Рис. 1.1. Подпрограмма вычисления коэффициента сжимаемости
В подпрограммах на рис. 1.1 и 1.2 приняты следующие идентификаторы: W – приведенная плотность; RO, ROKR – плотность и критическая плотность; TAU – приведенная температура; Т, TKR – температура и критическая температура; Р, PKR – давление и критическое давление; R – специфическая газовая постоянная; B(I,J) – коэффициенты разложения уравнения (1.5).



Рис. 1.2. Подпрограмма вычисления давления
Уравнение состояния (1.5) может быть решено в явном виде только относительно коэффициента сжимаемости и давления. Поэтому для определения плотности с помощью указанного уравнения применяем один из итерационных методов, в частности метод Ньютона, для чего вначале определяем частную производную , численное значение которой находим с помощью вспомогательной подпрограммы (рис. 1.3).


Р
(1.15)
ис. 1.3. Вспомогательная подпрограмма для определения плотности

рабочего вещества
Подпрограмма на рис. 1.3, в свою очередь, используется в основной расчетной подпрограмме CRRO(P,T,N) при определении плотности рабочего вещества по заданным значениям давления и температуры.

Плотность рабочего вещества определяем по подпрограмме, показанной на рис. 1.4.


Рис. 1.4. Подпрограмма расчета плотности рабочего вещества

в зависимости от фазового состояния

Как видим, сходимость решения осуществляется по заданной точности EPS, а число приближений задается величиной KN. Область состояния рабочего вещества задается значением переменной N: для жидкости N = 0 , для парогазовой смеси N =1. В соответствии с этим подпрограмма (см. рис. 1.4) для области жидкости записывается как CRRO(P,T,0), а для парогазовой части pvT-пространства – как CRRO(P,T,1).

Ha базе этой подпрограммы разработана универсальная подпрограмма (рис. 1.5), в которой с помощью логических операторов автоматически определяется область состояния реального газа и выбирается соответствующий алгоритм расчета плотности по давлению и температуре.


Рис. 1.5. Универсальная подпрограмма расчета плотности

по давлению и температуре
На рис. 1.5 FTS(P) – подпрограмма расчета температуры насыщения по давлению.

Энтальпию реального рабочего вещества определяем по плотности и температуре [см. уравнение (1.7)]. Подпрограмма расчета энтальпии изображена на рис. 1.6.


Рис. 1.6. Подпрограмма расчета энтальпии реального рабочего вещества
Для вычисления энтальпии рабочего вещества используем две вспомогательные подпрограммы: СР0(Т) – для определения изобарной теплоемкости [см. уравнение (1.12)] и Н0(Т) – для определения энтальпии [см. уравнение (1.10)] в идеально-газовом состоянии. Эти подпрограммы показаны на рис. 1.7 и 1.8.


Рис. 1.7. Подпрограмма расчета изобарной теплоемкости


Рис. 1.8. Подпрограмма расчета энтальпии в идеально-газовом состоянии

Для вычисления определенного интеграла в уравнениях (1.10) и (1.11) применяем стандартную процедуру GAUSS. В подпрограммах (см. рис. 1.6–1.8) приняты следующие идентификаторы: H0 – энтальпия в идеально-газовом состоянии; СР0 – изобарная теплоемкость в идеально-газовом состоянии; Н00 – энтальпия при T0 (см. таблицу); Н000 – теплота сублимации при T = 0 К.

Энтропию, согласно подпрограмме на рис. 1.9, рассчитываем подобно энтальпии (см. рис. 1.6) по плотности и температуре. При этом используем две вспомогательные подпрограммы (рис. 1.10 и 1.11): CP0T(T), в которых определяем значение подынтегрального отношения СР0/Т, и S0(T) – для вычисления энтропии в идеально-газовом состоянии [см. уравнение (1.11)].


Рис. 1.9. Подпрограмма расчета энтропии рабочего вещества


Рис. 1.10. Вспомогательная подпрограмма

для вычисления энтропии

Рис. 1.11. Подпрограмма для расчета энтропии

в идеально-газовом состоянии
В подпрограммах (см. рис. 1.9–1.11) использованы следующие идентификаторы: S0 – энтропия в идеально-газовом состоянии; S00 – энтропия при T0; S000 – константа уравнения (1.11).

Подпрограммы вычисления теплоемкости при постоянном давлении или объеме по плотности и температуре (рис. 1.12 и 1.13) разработаны на базе соотношения (1.9).


Рис. 1.12. Подпрограмма расчета теплоемкости

при постоянном давлении

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


Рис. 1.13. Подпрограмма расчета теплоемкости

при постоянном объеме

Рис. 1.14. Подпрограмма расчета давления насыщенных паров

На рис. 1.14 F(I) – коэффициент аппроксимации fi в уравнении (1.6).

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


Рис. 1.15. Подпрограмма расчета температуры насыщения
При анализе тепловых процессов в парогенерирующих элементах криогенных систем часто возникает необходимость в расчете теплоты парообразования. В разработанной системе обеспечения вычислительного эксперимента эта термодинамическая величина может быть найдена по изображенной на рис. 1.16 подпрограмме CRR(T).

Рис. 1.16. Подпрограмма расчета теплоты парообразования
При решении ряда задач гидродинамики важно знать скорость звука. Для ее вычисления используем подпрограмму, разработанную на базе соотношений, приведенных в [5] (рис. 1.17). В подпрограмме на рис. 1.17 VSO – скорость звука
.


Рис. 1.17. Подпрограмма расчета скорости звука
При моделировании тепловых и гидродинамических процессов в каналах проточных трактов элементов систем криогенного обеспечения объектов криостатирования часто возникает необходимость в определении калорических параметров рабочих веществ и скорости звука по известным значениям давления и температуры. Для этого на базе подпрограмм (см. рис. 1.5, 1.6, 1.9, 1.12, 1.13 и 1.17) разработаны соответствующие расчетные подпрограммы (рис. 1.18–1.22).

Рис. 1.18. Подпрограмма расчета энтальпии

по известным значениям давления и температуры
Базовые подпрограммы-функции (см. рис. 1.1–1.22) структурно входят в основное тело программы подсистемы расчета термических и калорических параметров криогенных рабочих веществ. Массивы коэффициентов расчетных уравнений и характерные параметры для каждого криоагента (см. таблицу) описываются в подпрограммах SUBROUTINE N2 для азота, SUBROUTINE О2 для кислорода и в аналогичных подпрограммах для других газов операторами DIMENSION и DATA, а обмен информацией осуществляется с помощью соответствующего оператора СОМMON. Подпрограммы-функции использованы при разработке программы расчета параметров реальных тепловых и гидродинамических процессов в элементах криогенных систем и установок.

Рис. 1.19. Подпрограмма расчета энтропии

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


Рис. 1.20. Подпрограмма расчета теплоемкости при постоянном давлении

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


Рис. 1.21. Подпрограмма расчета теплоемкости при постоянном объеме

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


Рис. 1.22. Подпрограмма расчета скорости звука

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

Для расчета параметров в двухфазной области применяем процедуру FAZA2(P), которая базируется на уравнении состояния (1.5) (рис. 1.23).

Рис. 1.23. Подпрограмма расчета параметров в двухфазной области
В начале решения находим температуру кипения, применяя подпрограмму CRPS(T) (см. рис. 1.14). Затем для локализации корней уравнения состояния, лежащих на левой и правой пограничных кривых области двухфазного состояния, в процедуре CRRO(P,T,N) принимаем соответствующие значения для N (см. рис. 1.4). Далее вычисляем энтальпию, энтропию жидкости и пара и теплоту парообразования.

П
Рис. 1.24. Схема расчета параметров в конце процесса дросселирования:

1–2 для p2 > pкр ;

1?–2? для p2 < pкр и ;

1?–2? для p2 = pнас и
ри расчете криоген­ных систем возникает необ­ходимость в определении па­раметров в конце процесса дросселирования (рис. 1.24). При этом последний может заканчиваться либо в области газовой фазы (процесс 1–2), тогда изменение темпера­туры перед дросселем или изменение отношения дав­лений будет вызывать изме­нение температуры оконча­ния процесса, либо в двух­фазной области (процесс 1–2), тогда температура конца процесса остается неизменной, равной темпе­ратуре насыщения пара при конечном давлении, а изменение параметров будет приводить к из­менению степени сухости рабочего тела в конце процесса дроссели­рования.

В соответствии с этим предусмотрены три процедуры.

  1. Процедура DROSG(P1,Т1,Р2), с помощью которой по начальным давлению и температуре и конечному давлению находим конечную температуру газовой фазы. В данной процедуре производится определение энтальпии в точке, характеризующей начало процесса дросселирования; далее методом половинного деления решается уравнение hн – hк = 0, из которого определяется температура конца дросселирования. Программа этой процедуры показана на рис. 1.25.

  2. Процедура DROS(P1,Т1,Р2) служит для расчета параметров газа в конце процесса дросселирования в двухфазной области (рис. 1.26). В этой процедуре по конечному давлению рассчитывается температура кипения, которая и является температурой конца процесса дросселирования.

По процедуре FAZA2(P) находим параметры двухфазной области, степень сухости и соотношение между пото­ками пара на входе в дроссель и на выходе из него.


Рис. 1.25. Подпрограмма расчета процесса дросселирования,

заканчивающегося в газовой фазе


Рис. 1.26. Подпрограмма расчета процесса дросселирования,

заканчивающегося в двухфазной области


  1. Процедура DROSL(P1,T1,P2) (рис. 1.27) служит для вычисления параметров в конце процесса дросселирования насыщенной жидкости, так как в установках часто осуществляется данный процесс (процесс 1–2 на рис. 1.23), что приводит к изменению энергетического состояния криоагента после дросселя.



Рис. 1.27. Подпрограмма расчета процесса дросселирования

насыщенной жидкости




На рис. 1.28 показаны процессы расширения рабочего вещества в парожидкостном детандере (1–2) и изоэнтропного сжатия (1–2). В обоих случаях стоит задача определения параметров состояния рабочего вещества в конце процессов для случая, когда заданы начальные давление рн и температура Тн, конечное давление рк и изоэнтропный КПД s.

В
Рис. 1.28. Схема расчета параметров в конце процессов расширения в детандере (1–2) и изоэнтропного сжатия (1–2)

первом случае по известным параметрам рн и Тн находим энтальпию hн и энтропию sн, а по значению рк с помощью проце­дуры FAZA2(PK) – температуру кипения Ткип, энтальпию hL и hG, энтропию sL и sG на пограничных кривых области двухфазного состояния. Далее по известным значениям sн, sL, sG определяем положение точки 2s, а по заданному значению КПД (CPD) – истинное состояние рабочего вещества в конце расширения (рис. 1.29).


Рис. 1.29. Подпрограмма расчета параметров после расширения в детандере
Для определения состояния рабочего вещества в конце изо­энтропного сжатия применяем процедуру SCONST(PH,TH,PK,K) – рис. 1.30.


Рис. 1.30. Подпрограмма расчета параметров после изоэнтропного сжатия

По этой процедуре находим значения энтальпии и энтропии в точке 1 и определяем температуру . Поскольку показатель k в зависимости от области состояния может отличаться от принятого первоначального значения, найденное значение температуры Тs является первым приближением, которое в последующем уточняется в результате решения методом половинного деления уравнения sнsк = 0.

Часто при проведении вычислительных экспериментов возникает необходимость расчета температуры рабочего вещества по известным значениям давления и энтальпии, для этого может быть использована процедура TH(P,H,N) – рис. 1.31.

Рис. 1.31. Подпрограмма расчета температуры рабочего вещества по известным значениям давления и энтальпии
При моделировании тепловых и гидродинамических процессов необходимо иметь информацию о транспортных свойствах рабочих веществ – вязкости и теплопроводности. Теплопроводность может быть найдена по полиному, полученному при аппроксимации экспериментальных данных,
. (1.14)
Реализация решения полинома (1.14) приводится в процедуре FLA(RO,T) – рис. 1.32.


Рис. 1.32. Подпрограмма расчета теплопроводности

На рис. 1.32 FLA – теплопроводность, Вт/(мК).

Динамическую вязкость рассчитываем по следующему уравнению:
. (1.15)

Слагаемые уравнения (1.15) вычисляем по формулам
; (1.16)
. (1.17)
Для определения динамической вязкости по уравнениям (1.15)–(1.17) в программе предусмотрена процедура FMU(RO,T) – рис. 1.33.


Рис. 1.33. Подпрограмма расчета динамической вязкости
На рис. 1.33 FMU – динамическая вязкость, Пас.

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

Уравнение Редлиха–Квонга для смесей в модификации Соаве имеет вид [12]
, (1.18)
где zm – коэффициент сжимаемости для смеси; m – молярный объем смеси, м3/кмоль; A и B – простые числа (A = 0,4274802327; B = 0,086640350);

, (1.19)
здесь yj – молярная доля j-го компонента;
; (1.20)

Mj – молярная масса компонента смеси.

Коэффициент Fm учитывает параметры бинарных взаимодействий компонентов смеси и, согласно модификации Соаве, определяется как
, (1.21)
где i, j (i-й и j-й) – парные компоненты смеси; – параметр бинарного взаимодействия компонентов; Fi, Fj – коэффициенты компонентов, которые определяются по уравнению
, (1.22)
здесь – фактор ацентричности молекул вещества.

На базе алгоритма (1.1) и выражений (1.18)–(1.22) разработаны блок-схемы расчета коэффициента сжимаемости, давления и плотности смеси.

На рис. 1.34 показана блок-схема расчета коэффициентов уравнения (1.18).

Рис. 1.34. Блок-схема алгоритма расчета коэффициентов уравнения (1.18)

Структурно в блок-схеме реализованы три цикла, позволяющих определить суммы, зависящие от числа компонентов N и состава смеси,

.
В блоке 1 осуществляется ввод исходной информации, характеризующей состав и свойства компонентов смеси.

В блоке 2 рассчитываются индивидуальные коэффициенты Fi (F(I)) по уравнению (1.22), bi (В(I)) по (1.20), а также суммарный коэффициент bm (ВМ) по (1.19).

В блоках 4 и 5, которые структурно входят в циклы по I и по J, рассчитываются знаменатель и числитель уравнения (1.21), а в блоке 8 завершается его решение определением коэффициента Fm (FM)*, учитывающего параметры бинарных взаимодействий компонентов смеси.

Блоки 3, 6 и 7 организуют расчетные циклы по I и J по числу компонентов смеси от 1 до N.

Вывод и передачу рассчитанных величин в соответствующие подпрограммы и процедуры осуществляет блок 9.

На основании приведенной блок-схемы разработана и реализована подпрограмма MIXTURE.

Для расчета давления смеси [см. уравнение (1.1)] по известной плотности и температуре может быть использована блок-схема (рис. 1.35).

В

Рис. 1.35. Блок-схема алгоритма расчета давления
блоке 1 вводятся значения плотности, температуры, универсальная газовая постоянная, константы уравнения Редлиха–Квонга и коэффициенты, предварительно рассчитанные в подпрограмме MIXTURE.

В блоке 2 по уравнению (1.18) определяется значение коэффициента сжимаемости смеси, а в блоке 3 решается уравнение состояния реального вещества относительно давления.

Вычисление плотности с использованием уравнения (1.18) производится в соответствии с блок-схемой (рис. 1.36).

Рис. 1.36. Блок-схема алгоритма расчета плотности смеси

Уравнение Редлиха–Квонга может быть решено в явном виде только относительно коэффициента сжимаемости и давления. Поэтому для  определения плотности применяется один из итерационных методов, в частности метод Ньютона, которым решается уравнение
. (1.23)
Для локализации корней уравнения состояния принимаем первоначальные значения плотности пара и жидкости, которые идентифицируются целочисленной переменной М: для жидкости М = 0, для пара М = 1 (см. рис. 1.36, блок 2).

В блоках 4 и 6 обеспечивается реализация метода Ньютона, для чего при заданных возможных значениях плотности при постоянной температуре с помощью уравнения (1.23) находим частную производную , уточненное значение плотности и степень сходимости DELTA решения уравнения (1.23).

После сравнения (в блоке 7) DELTA с заданной точностью решения EPS принимаем решение либо об обращении в блок 3 для уточнения исходного значения плотности и повторения расчетов, либо, при достижении необходимой сходимости, о передаче информации в блок 8 вывода информации, на чем решение задачи заканчивается.

2. СТАТИЧЕСКИЕ ХАРАКТЕРИСТИКИ

ЭЛЕМЕНТОВ КРИОГЕННЫХ СИСТЕМ
2.1. Моделирование статических нерасчетных режимов

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

Криогенную систему можно представить в виде конструкции, состоящей из следующих ступеней: ступени подготовки рабочего вещества (СПТ), на которой производятся его сжатие с отводом теплоты в окружающую среду и очистка от примесей; ступени предварительного охлаждения (СПО); ступени окончательного охлаждения (СОО), на которой ведется процесс охлаждения путем преобразования внутренней энергии рабочего вещества при дросселировании, расширении в детандере, выхлопе и т. п.; ступени использования эффекта охлаждения (СИО), на которой обеспечивается вывод из системы субстанции либо в виде теплоты, либо в виде конденсированного криопродукта. Структурно ступени СИО и СОО могут частично или полностью совпадать. Часто в структурную схему системы включается дополнительная криогенная установка, предназначенная для получения внешнего криогенного продукта, используемого в СПО для отвода теплоты от потока с высоким давлением [1,13, 14].

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

К радиационным теплообменным аппаратам прежде всего следует отнести элементы объекта криостатирования с радиационным теплоподводом и внутристеночным выделением теплоты, участки криотрубопроводов и другие элементы криогенных систем.

Современные системы криогенного обеспечения рассчитывают с учетом последних достижений криогенной техники [13, 15]. Они способны работать при оптимальных параметрах с высокими энергетическими показателями. Однако при эксплуатации могут возникать ситуации, которые переводят отдельные элементы или всю систему криогенного обеспечения на нерасчетные режимы работы. Для выявления степени влияния отдельных возмущающих факторов на параметры потоков рабочих сред и оболочки каналов может быть использован метод математического моделирования [10, 16, 17].

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

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

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

Если математические модели теплообменных аппаратов рассматриваются как одномерные в направлении оси x, тогда система уравнений, описывающая исходный равновесный режим или новое устойчивое состояние рекуперативного теплообменника, имеет вид:

уравнения сплошности
; (2.1)
уравнения энергии
; (2.2)
уравнения состояния

(2.3)

уравнения движения
, (2.4)
где G , G  – расход вещества прямого и обратного потоков; h, h – энтальпия потоков; ?, ? – коэффициенты теплоотдачи потоков; Fi, Fi – площадь теплообменной поверхности наружной и внутренней стенок оболочки канала; T , T  – температура прямого и обратного потоков; ,  – температура наружной и внутренней поверхностей оболочки каналов; p, p – давление прямого и обратного потоков; ?, ? – плотность потоков; – потери давления в потоках.

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

Математическая модель противоточного рекуперативного теплообменного аппарата. На рис. 2.1 показана физическая модель рассматриваемого теплообменного аппарата.


а

б

Рис. 2.1. Физическая модель рекуперативного теплообменного аппарата:

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

1 – прямой поток; 2 – оболочка канала; 3 – обратный поток
Рекуперативный теплообменник состоит из трех подсистем: оболочки канала и двух рабочих сред – прямого G  и обратного G  потоков, которые имеют на входе параметры Для расчета распределения температур вдоль ординаты x может быть применен метод элементарных балансов [10, 18]. Для решения поставленной задачи поверхность теплообменного аппарата разбивается на n участков Fi = F/n (рис. 2.2) и для i-го участка каждой подсистемы записываются уравнения теплового баланса:

(2.5)

где  – теплоемкость в i-м сечении прямого и обратного потоков соответственно; – коэффициенты теплоотдачи в i-м сечении потоков; – площадь поверхности теплообмена элементарного участка внутренней и наружной поверхностей оболочки канала.



Рис. 2.2. Расчетная схема к определению статического распределения термических параметров по высоте

рекуперативного теплообменного аппарата
Выражения для расчета температуры каждой подсистемы в i-м сечении могут быть записаны следующим образом:
(2.6)
Для решения на ЭВМ система уравнений (2.6) может быть представлена в виде
(2.7)
где – коэффициенты, характеризующие условия теплообмена на границах подсистем прямого и обратного потоков и оболочки канала.

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

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

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

, (2.8)

где n, n– необходимое расчетное число участков прямого и обратного потоков соответственно.

Истинное число участков должно быть больше или равно большему из чисел n и n.

Блок-схема реализации математической модели рекуперативного теплообменного аппарата изображена на рис. 2.3. С помощью блока 1 производится ввод исходных данных. Блок 2 осуществляет обращение к программе расчета термодинамического поля рабочих веществ. В блоке 3 вычисляется начальное распределение температуры потоков рабочих сред. После задания начальных условий счета по соответствующим процедурам в блоке 4 рассчитываются теплоемкость и другие теплофизические параметры рабочих веществ, необходимые для расчета коэффициентов теплоотдачи. Локальные коэффициенты теплоотдачи потери давления pi вычисляются по соответствующим процедурам, вид которых выбирается исходя из принятой конструкции теплообменного аппарата, геометрии каналов и условий теплообмена. В блоке 5 происходит определение коэффициентов уравнений типа (2.7), а в блоке 6 – расчет методом итерации с заданной точностью EPS температуры в i-м сечении каждой подсистемы и давлений потоков. В блоке 7 происходит наращивание цикла расчета и передача информации в блоки 4 и 5. После завершения всего цикла результаты выводятся на печать или информация передается соответствующим подпрограммам математической модели криогенной системы (блок 8).

На рис. 2.4 отображены результаты расчета и экспериментального исследования распределения температуры потоков гелия при высоком и низком давлении в теплообменнике нижней ступени охлаждения сателлитной гелиевой установки с избыточным обратным потоком. Показано распределение температур в теплообменном аппарате типа «труба в трубе» при следующих параметрах потоков: G  = 0,0192 кг/с; = (– )/= 0,08 при p1= 2,0 MПa, p2 = 0,13 МПа. Необходимо отметить вполне удовлетворительное подтверждение расчетных зависимостей экспериментальными данными [18].



Рис. 2.3. Блок-схема алгоритма расчета рекуперативного

теплообменного аппарата



Рис. 2.4. Распределение температуры потока в теплообменном аппарате:

расчет; – экспериментальные точки
Математическая модель рекуперативного теплообменного аппарата с кипением криоагента в межтрубном пространстве. На рис. 2.5 показана физическая модель указанного аппарата.


Рис. 2.5. Физическая модель рекуперативного теплообменного

аппарата с кипением криоагента:

а – расчетная схема; б – распределение температуры по длине аппарата;

1 – прямой поток; 2 – оболочка канала;

3 – кипящая жидкость; 4 – насыщенный пар
Рассматриваемый теплообменник относится к испарительным аппаратам и состоит из трех подсистем: прямого потока, оболочки канала и кипящего криоагента.

Тепловой баланс i-го участка для подсистем Прямой поток и Оболочка канала может быть записан в виде (рис. 2.6)
. (2.9)


Рис. 2.6. Расчетная схема к определению
статического распределения параметров
в подсистемах испарительного теплообменного аппарата

Значение коэффициента теплоотдачи будет определяться теплофизическими свойствами кипящего криоагента, режимом процесса кипения и условиями его осуществления. Для пузырькового режима кипения в большом объеме коэффициент Зная, что данное уравнение можно записать следующим образом [19, 20]:
. (2.10)
  1   2   3   4   5


Учебный материал
© bib.convdocs.org
При копировании укажите ссылку.
обратиться к администрации