Барилович В.А. Основы термогазодинамики двухфазных потоков и их численное решение
Подождите немного. Документ загружается.
В.А. БАРИЛОВИЧ
ОСНОВЫ
ТЕРМОГАЗОДИНАМИКИ
ДВУХФАЗНЫХ ПОТОКОВ И ИХ
ЧИСЛЕННОЕ РЕШЕНИЕ
Учебное пособие
Санкт-Петербург
Издательство Политехнического университета
p
w
п
w
гк
G
п
G
гв
G
хв
w
хк
w
п
Федеральное агентство по образованию
САНКТ-ПЕТЕРБУРГСКИЙ
ГОСУДАРСТВЕННЫЙ ПОЛИТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ
ОБЩИЕ ПРОФЕССИОНАЛЬНЫЕ ДИСЦИПЛИНЫ
В ПОЛИТЕХНИЧЕСКОМ УНИВЕРСИТЕТЕ
В. А. БАРИЛОВИЧ
Основы
термогазодинамики
двухфазных потоков и их
численное решение
Учебное пособие
Санкт-Петербург
Издательство Политехнического университета
2009
УДК 621.1.013
Основы термогазодинамики двухфазных потоков и их численное решение:
Учебное пособие / В.А.Барилович. СПб.: Изд. Политехнического
университета , 2009, 425с.
Рассматривается течение двухфазных потоков, содержащих частицы,
капли и пузырьки в каналах переменного сечения с прямой и криволинейной
осью в абсолютном и относительном движении. Приводятся системы
дифференциальных уравнений, решение которых реализовано в виде программ
на языках Фортран, Бейсик и Паскаль, представлены результаты расчетов.
Данное учебное пособие по сравнению с
первым и вторым (вышли в 2001
и 2005г. в издательстве “Нестор“) содержит значительный собственный
экспериментальный материал по исследованию процессов во вскипающем и
двухфазном потоках в каналах энергетических устройств. Изучаются
кавитационные явления при движении холодной капельной жидкости в канале
переменного сечения. Значительное внимание уделяется кризисным явлениям в
спутных и двухфазных потоках. Рассматриваются
процессы в
парогенерирующих каналах. Приводятся результаты опытов по конденсации
паровой струи в спутном потоке холодной воды. Изложены результаты
охлаждения двухфазным потоком высоконагретых поверхностей.
Учебное пособие расширяет раздел двухфазных течений в общем курсе
гидрогазодинамики и служит дополнительным учебным материалом для
самостоятельной и научной работы студентов бакалавриата и магистратуры
энергомашиностроительного факультета.
Илл.79,
Табл. 8, Библ. 57 назв.
Компьютерная графика в пособии и программы на языке “BASIC”
выполнены доц. каф. ТОТ СПбГПУ Смирновым Ю.А., за что автор выражает
ему глубокую благодарность.
Рецензенты: проф. кафедры "Теоретические основы теплотехники"
СПбГПУ, д.т.н. Э.Л. Китанин;
Зав. кафедрой "Гидроаэродинамика" СПбГПУ проф., д-р физ.- мат. наук
Е.М. Смирнов;
Главный научный сотр. ин-та проблем Машиноведения РАН,
д.т.н. В.В. Фисенко;
Зав. лабораторией ОАО ЦКТИ д.т.н., проф. СПбГПУ Б.С.Фокин.
© Владимир Антонович Барилович
© Санкт-Петербургский государственный
политехнический университет, 2009
3
Введение
Основные разделы курса гидрогазодинамики, читаемого студентам
энергомашиностроительного факультета, посвящены изучению движения
однородных потоков. Такие потоки реализуются в элементах проточных частей
газовых турбин, компрессоров, гидравлических турбин и насосов. При этом,
как правило, считают, что газовые среды подчиняются уравнению состояния
идеального газа, а жидкость - несжимаема.
Такое упрощение в большинстве
случаев оправдано и не приводит к
заметным погрешностям в расчетах. Однако в ряде энергетических машин и
устройств в процессе расширения потока вещество из однофазного переходит в
двухфазное состояние, и расчет в предположении однородной среды становится
неправомерным. Подобные явления имеют место в последних ступенях
паровых турбин, когда в паре образуются капельки
влаги; во входных
устройствах насосов, где в зоне пониженного давления может возникать
“холодное” кипение (из жидкости выделяются пузырьки газа и пара),
приводящее к снижению энергетических и прочностных характеристик
насосов. В холодильных машинах в процессе дросселирования жидкого
хладоагента происходит его вскипание. Аналогичные явления имеют место в
турбинах, работающих на вскипающих потоках,
где в процессе расширения
вначале образуются пузырьки, а затем капли. В инжекторах происходит
конденсация пара на каплях охлаждающей жидкости. В камере сгорания
газотурбинной установки воздушный поток смешивается с каплями топлива. У
ракетного двигателя, работающего на твердом топливе, в сопле движется газ,
несущий твердые частицы. В струйных насосах ускорение жидкой фазы
осуществляется
газовым потоком. Подобные примеры из техники можно было
бы продолжить.
Таким образом, будущему инженеру-энергомашиностроителю
необходимо изучение раздела “Двухфазные потоки”.
Отметим, что успех решения любой научной задачи заключается в
правильном понимании физического явления и создании на этой основе
адекватных физической и математической моделей. В нашем случае
необходимо ясно представлять,
какой режим течения реализуется в тех или
иных условиях.
Существует большое число режимов течения двухфазных потоков,
которые зависят от соотношения массовых расходов компонентов, плотности,
вязкости фаз и коэффициента поверхностного натяжения. К основным режимам
можно отнести пузырьковый, пробковый, расслоенный, кольцевой, дисперсно-
кольцевой и дисперсный, которые в первом приближении могут быть
определены
по диаграммам Бейкера [1], исследовавшего газожидкостные
4
смеси, и Б.К. Козлова [2], изучавшего структуру двухфазного потока в
вертикальных трубах.
В методическом пособии приводятся уравнения сохранения
применительно к гомогенным, дисперсно-кольцевым, дисперсным и
пузырьковым потокам.
В первом разделе рассматривается движение двухфазного потока в
каналах переменного сечения с прямой осью. Во втором – течение двухфазного
потока в криволинейном канале,
движение капель в абсолютной и
относительной системах координат. В третьем – изложены основные
положения математического описания пространственного гетерогенного
потока. Четвертый раздел посвящен экспериментальным исследованиям сопел,
инжекторов и турбины, работающих на перегретой воде. Рассматриваются
процессы в спутных конденсирующихся потоках и охлаждение высоконагретых
поверхностей двухфазными средами.
Численное решение систем обыкновенных дифференциальных
уравнений в
первом и втором разделах осуществляется методами Рунге-Кутта и
Рунге-Кутта-Фельдберга.
В качестве тестового примера, позволяющего оценить точность расчета,
рассматривается расширение идеального газа в сопле Лаваля.
5
Раздел I. Одномерный двухфазный поток
1. Течение идеального газа в сопле Лаваля
В разделе курса гидрогазодинамики “Одномерная газодинамика
стационарных течений” студенты изучают течение идеального газа в каналах и
соплах, выполняют расчет идеального сопла Лаваля с использованием
газодинамических функций. Расчет ведется на основе интегральных
соотношений, включающих уравнения сплошности, энергии и процесса
для
адиабатического канала, уравнения состояния идеального газа и уравнения,
связывающего площадь поперечного сечения канала с продольной
координатой:
wF
ρ
= const, i*= const , p/
ρ
k
= const , p/(
ρ
T)=R= const , F=
ϕ
(x) ,
где k=c
p
/c
v
- показатель изоэнтропы, с
p
- c
v
= R .
Чтобы решить эту задачу методом Рунге-Кутта, необходимо представить
уравнения в дифференциальной форме. Для этого продифференцируем
вышеприведенные уравнения, записав вместо уравнения процесса уравнение
движения (Бернулли) для обратимого течения:
⎪
⎪
⎪
⎪
⎪
⎩
⎪
⎪
⎪
⎪
⎪
⎨
⎧
=+
=−−
=+
−=+
.0
0
0
2
2
p
dp
w
dw
p
w
d
T
dT
p
dp
T
dT
w
dw
Tc
w
F
dFd
w
dw
p
ρ
ρ
ρ
ρ
ρ
(1.1)
Главный определитель такой системы уравнений имеет вид
001
1110
010
1001
2
2
p
w
Tc
w
p
ρ
−−
=Δ
. (1.2)
Раскрывая (1.2), найдем, что
6
1
)1(
1)1(
2
2
2
2
2
22
−=
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
−
−
+−=−+−=Δ M
a
kw
a
wk
p
w
Tc
w
p
ρ
, (1.3)
где
akR
T
= - термодинамическая скорость звука в идеальном газе; М=w/a -
критерий Маха.
Разрешив (1.1) относительно
dp/p , получим
F
dF
M
kM
p
dp
)1(
2
2
−
−=
. (1.4)
Отметим, что выражение (1.4) можно получить, используя формулу Крамера [3]
dp
p
p
=
Δ
Δ
, где Δ
p
- определитель, получающийся из Δ заменой столбца,
составленного из коэффициентов при
dp/p, столбцом, составленным из
свободных членов. Для нашего случая имеем:
F
dF
kM
F
dF
a
kw
F
dF
kRT
kw
F
dF
p
w
p
w
Tc
w
F
dF
p
p
2
2
222
2
2
000
1100
010
101
−=−=−=−=
−−
−
=Δ
ρ
ρ
.
Тогда
F
dF
M
kM
p
dp
p
)1(
2
2
−
−=
Δ
Δ
=
, или
dx
dF
F
p
M
kM
dx
dp
)1
2
(
2
−
−=
. (1.5)
Остальные производные из (1.1) перепишем в виде
dx
dp
wdx
dw
ρ
1
−=
, (1.6)
dx
dw
w
T
Mk
dx
dT
2
)1( −−= , (1.7)
dx
dT
Tdx
dp
pdx
d
ρ
ρ
ρ
−= . (1.8)
Уравнения (1.6)-(1.8) можно также преобразовать к виду:
7
dx
dF
FM
w
dx
dw
)1(
2
−
=
, (1.9)
()
(
)
dx
dF
T
FM
Mk
dx
dT
1
1
2
2
−
−
−=
, (1.10)
dx
dF
F
M
M
dx
d
ρρ
)1(
2
2
−
−=
. (1.11)
Используя понятие логарифмической производной, перепишем систему
уравнений в иной форме:
dx
Fd
M
kM
dx
pd
ln
)1(
ln
2
2
−
−=
, (1.12)
dx
Fd
M
dx
wd ln
)1(
1ln
2
−
=
, (1.13)
dx
Fd
M
Mk
dx
Td ln
)1(
)1(
ln
2
2
−
−
−=
, (1.14)
dx
Fd
M
M
dx
d ln
)1(
ln
2
2
−
−=
ρ
. (1.15)
Полученная система обыкновенных дифференциальных уравнений
позволяет при заданной функции
F=f(x) определить производные dp/dx, dw/dx,
dT/dx, d
ρ
/dx и, используя численный метод, найти распределение p,w,T и
ρ
вдоль координаты
x при заданных граничных условиях : x=0, p=p
o
*, w=w
o
,
T=T
o
*
,
ρ
=
ρ
o
*
(задача Коши).
При изоэнтропийном критическом и сверхкритическом режимах
истечения в минимальном сечении канала
ρ
w=(
ρ
w)
max
, M=1, а dF/dx=0, т.е.
имеет место неопределенность, так как
dp/dx=-0/0 , которую можно раскрыть по
правилу Лопиталя. Для этого перепишем уравнение (1.12) в виде
2
1
1
ln
⎟
⎠
⎞
⎜
⎝
⎛
−
−=
w
a
dx
dF
F
k
dx
pd
(1.16)
и возьмем отдельно производные от числителя и знаменателя. Для числителя
будем иметь:
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
−+−=
⎟
⎠
⎞
⎜
⎝
⎛
−
dx
dF
F
F
dF
dx
Fd
F
k
dx
dF
Fdx
d
k
22
2
111
,
но в критическом сечении
dF/dx=0, тогда числитель примет вид
8
2
2
dx
Fd
F
k
−
. (1.17)
Для знаменателя напишем окончательное выражение
dx
pd
k
k
w
a
dx
d
ln
1
2
2
⎟
⎠
⎞
⎜
⎝
⎛
+
−=
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
−
. (1.18)
Здесь мы не приводим вывод (1.18), ввиду его громоздкости, но для тех,
кто хотел бы его выполнить, подскажем ряд положений. Так, из уравнения
Бернулли, записанного для изоэнтропийного процесса
wdw=-dp/
ρ
, можно
получить
dx
pd
kw
a
dx
pd
kw
kpv
dx
pd
w
vp
dx
wd
lnlnln
ln
2
2
22
−=−=−= .
Дифференцируя уравнение состояния
v=RT/p, найдем, что
dx
pd
dx
Td
dx
vd lnlnln
−=
.
Используя (1.7), можно написать:
dx
wd
a
w
k
dx
Td
ln
)1(
ln
2
2
−−= .
Кроме того, в критическом сечении
w=a и т.д. Подставляя (1.17) и (1.18) в
(1.16), получим
Fk
dxFd
k
dx
pd
)1(
/
ln
22
+
±=
, (1.19)
или
2
2
2
)1(
dx
Fd
Fk
p
k
dx
dp
+
±=
. (1.20)
Уравнение (1.20) определяет наклон кривых давления в критической точке.
Видно, что имеет место двойственность (±), то есть вниз по потоку от
критического сечения возможно развитие процесса по одному из двух
направлений: с ростом давления (
dp/dx>0) и с уменьшением (dp/dx<0), когда
поток становится сверхзвуковым. Численные эксперименты показывают, что
при решении прямой задачи, даже в случае, когда расход равен критическому
&
()GFw
ΚΡ Γ ΚΡ
=
ρ
, процесс протекает при dp/dx>0. Для того чтобы процесс
пошел по нижней ветке кривой (
dp/dx<0), достаточно незначительно увеличить
расход сверх критического. Для рассматриваемого случая (
p
0
=6⋅10
5
Па, T
0
=426
K,
k=1,4 , D
0
=3⋅10
-2
м, D
Г
=8⋅10
-3
м, D
СР
=1,2⋅10
-2
м, х
К
=4⋅10
-2
м, х
МИН
=1⋅10
-3
м,
х
D
= 0,105 м, где D и х – диаметр и протяженность конвергентной,
9
цилиндрической и дивергентной частей сопла Лаваля) отношение
&&
GG
ΚΡ
равнялось 1,007939. Получено оно следующим образом: рассчитывается
критический расход
&
()GFw
ΚΡ Γ ΚΡ
=
ρ
, затем скорость газа в нулевом сечении
при критическом расходе
wG F
000
=
&
//
ΚΡ
ρ
, далее эту скорость увеличивают
до тех пор, пока не выполнится условие dp/dx<0. В данном примере скорость
была увеличена с 17,034 до 17,200 м/с. Таким образом, решение прямой задачи
связано с определенными трудностями, но они, как видим, преодолимы.
Решение задачи при таком подходе не является строгим и возможно лишь
тогда, когда особая
точка находится внутри шага интегрирования.
Следовательно, в этом случае невозможно определить величину минимального
сечения канала и критический расход с высокой точностью. Изменение шага
интегрирования может привести к неустойчивости решения в виде осцилляций
(всплесков) в области горла. Отметим, что при решении прямых задач с целью
снижения трудоемкости вычислительного характера необходимо вначале
реализовать обратную задачу (см. ниже), что позволит в первом приближении
определить характерные поперечные сечения канала.
Если задан (известен) закон распределения статического давления вдоль
сопла p=f(x) (обратная задача), то вышеупомянутая проблема не возникает.
Опытное исследование сопел Лаваля, работающих на вскипающей воде,
показывает, что в качестве такой функции можно использовать зависимость
11
0
0
+
⎟
⎠
⎞
⎜
⎝
⎛
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
−
=
n
cp
l
x
p
p
p
p
, (1.21)
где l- длина сопла; p
ср
- статическое давление на выходном срезе сопла, n=2…6.
В этом случае система дифференциальных уравнений при n = 2 примет вид:
()
⎪
⎪
⎪
⎪
⎪
⎪
⎩
⎪
⎪
⎪
⎪
⎪
⎪
⎨
⎧
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
−=
−
−=
−=
⎥
⎥
⎦
⎤
⎢
⎢
⎣
⎡
+
⎟
⎟
⎠
⎞
⎜
⎜
⎝
⎛
−
−
−=
.
11
,
)1(
,
1
,
11/
)1/(2
0
2
2
0
2
00
dx
dT
Tdx
dp
pdx
d
dx
dw
R
w
k
k
dx
dT
dx
dp
wdx
dw
l
x
ppl
xppp
dx
dp
cp
cp
ρ
ρ
μ
ρ
(1.22)
Задаваясь расходом
&
G и решая систему уравнений (1.22), находим текущие
значения скорости, плотности и температуры.