Научная статья на тему 'Трехжидкостная формулировка и численный метод решения стационарной задачи теплогидравлики двухфазного потока в дисперсно-кольцевом режиме течения'

Трехжидкостная формулировка и численный метод решения стационарной задачи теплогидравлики двухфазного потока в дисперсно-кольцевом режиме течения Текст научной статьи по специальности «Физика»

CC BY
122
20
i Надоели баннеры? Вы всегда можете отключить рекламу.
Ключевые слова
ДВУХФАЗНОЕ ПАРОВОДЯНОЕ ТЕЧЕНИЕ / ТРЕХЖИДКОСТНАЯ МОДЕЛЬ / ЧИСЛЕННОЕ РЕШЕНИЕ / TWO-PHASE STEAM-WATER FLOW / THREE-FLUID MODEL / NUMERICAL SOLUTION

Аннотация научной статьи по физике, автор научной работы — Авдеев Евгений Эдуардович, Плетнев Александр Александрович, Булович Сергей Валерьевич

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

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

Three-fluid formulation and a numerical method for solving the stationary problem of thermal hydraulics of a two-phase annular dispersed flow

The article presents a one-dimensional three-fluid model for solving the stationary flow problem on a two-phase steam-water annular dispersed stream in a vertical heated channel, based on nine balance equations with an equal phase pressure assumed. The validation of the marching algorithm of the described stationary problem has been carried out by comparison with the published experimental data relating to a two-phase flow in a circular pipe under adiabatic conditions for pressures of 3 – 9 MPa, total flow rates of 500 – 3000 kg /(m^2∙s) and internal diameters of 10 and 20 mm. The total pressure differences in the channel were calculated. Good qualitative agreement with experimental data was obtained. Small quantitative disagreements were found. They were shown to be reduced by the refinement of the closing relations.

Текст научной работы на тему «Трехжидкостная формулировка и численный метод решения стационарной задачи теплогидравлики двухфазного потока в дисперсно-кольцевом режиме течения»

doi: 10.18721/jpm.11311 удк 536.7; 532.5; 519.6

трехжидкостная формулировка и численный метод РЕШЕНИЯ стационарной ЗАДАЧИ ТЕПЛОГИДРАВЛИКИ

двухфазного потока в дисперсно-кольцевом

РЕЖИМЕ ТЕЧЕНИЯ Е.Э. Авдеев, А.А. Плетнев, С.В. Булович

Санкт-Петербургский политехнический университет Петра Великого, Санкт-Петербург, Российская Федерация

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

Ключевые слова: двухфазное пароводяное течение, трехжидкостная модель, численное решение

Ссылка при цитировании: Авдеев Е.Э., Плетнев А.А., Булович С.В. Трехжидкостная формулировка и численный метод решения стационарной задачи теплогидравлики двухфазного потока в дисперсно-кольцевом режиме течения // Научно-технические ведомости СПбГПУ. Физико-математические науки. 2018. Т. 11. № 3. С. 122-132. DOI: 10.18721/JPM.11311

three-fluid formulation and a numerical method

FOR SOLVING THE STATIONARY PROBLEM OF THERMAL HYDRAuLICS OF A TWO-PHASE ANNuLAR DISPERSED FLOW

E.E. Avdeev, A.A. Pletnev, S.V. Bulovich

Peter the Great St. Petersburg Polytechnic University, St. Petersburg, Russian Federation

The article presents a one-dimensional three-fluid model for solving the stationary flow problem on a two-phase steam-water annular dispersed stream in a vertical heated channel, based on nine balance equations with an equal phase pressure assumed. The validation of the marching algorithm of the described stationary problem has been carried out by comparison with the published experimental data relating to a two-phase flow in a circular pipe under adiabatic conditions for pressures of 3 — 9 MPa, total flow rates of 500 — 3000 kg /(mA2s) and internal diameters of 10 and 20 mm. The total pressure differences in the channel were calculated. Good qualitative agreement with experimental data was obtained. Small quantitative disagreements were found. They were shown to be reduced by the refinement of the closing relations.

Key words: two-phase steam-water flow, three-fluid model, numerical solution

Citation: E.E. Avdeev, A.A. Pletnev, S.V. Bulovich, Three-fluid formulation and a numerical method for solving the stationary problem of thermal hydraulics of a two-phase annular dispersed flow, St. Petersburg Polytechnical State University Journal. Physics and Mathematics. 11 (3) (2018) 122—132. DOI: 10.18721/ JPM.11311

Введение

Для описания состояния многофазных потоков в теплогидравлическом приближении наиболее распространен подход, основанный на модели взаимопроникающих континуумов. Известные теплогидравли-ческие расчетные коды (КОРСАР, TRAC, RELAP) базируются на так называемом двухжидкостном приближении. Оно относится к случаю, когда многофазный поток представлен двумя «жидкостями»: паром и жидкостью. Соответствующая система основывается на шести дифференциальных уравнениях: баланса массы, импульса и энергии (для каждой из жидкостей). Система замыкается уравнением баланса фаз, уравнением термодинамического состояния (для каждой из жидкостей) и эмпирическими или полуэмпирическими соотношениями (для каждой из жидкостей). Они описывают массообмен, теплообмен и обмен импульсами между фазами и фаз со стенкой канала (называются замыкающими соотношениями).

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

Двухжидкостное приближение было первым в хронологии развития теплоги-дравлических расчетных кодов. Оно до сих пор применяется во многих кодах и часто встречается в литературе. В целом двухжид-костный подход не устарел, так как хорошо описывает все режимы течения, где двухфазный поток было бы правильнее разделить на две любые составляющие. К таким режимам относятся, например, пузырьковый, снарядный, дисперсный, расслоенный. Для всех указанных режимов течения двухжидкостной подход дает приемлемые результаты; и такие описания будут различаться лишь определением «жидкостей» в каждом случае и, соответственно, замыкающими соотношениями.

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

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

Проблема неточного описания дисперсно-кольцевого режима решается заменой двухжидкостного подхода трехжид-костным. В нем каждая из жидкостей, т. е. пар, капли и жидкая пленка, описывается собственными уравнениями баланса массы, импульса и энергии. Трехжидкостная модель не является новшеством. Она разрабатывалась и ранее в ряде отечественных и зарубежных публикаций [1 — 7]. Переход от описания двух жидкостей к трем не несет с собой алгоритмических новаций, но позволяет построить более полную физическую модель, неравновесную по скоростям и температурам рассматриваемых жидкостей.

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

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

В данной работе представлен маршевый алгоритм для численного решения стационарной задачи течения дисперсно-кольцевого режима в одномерном приближении с использованием трехжидкостного подхода. Известно, что при рассмотрении систем в двух- и трехжидкостном приближении в нестационарной постановке, в случае равновесного значения давления (математические модели с общем давлением в жидкостях), при определенном значении режимных параметров утрачивают эволюционные свойства [8]. Восстановление корректности задачи Коши может быть достигнуто различными способами, но все эти приемы «возмущают» исходную систему уравнений. Поэтому, если в рассматриваемой задаче возможен стационарный режим течения, то решение, полученное маршевым методом интегрирования, является «эталонным» для эволюционных задач и тем самым позволяет количественно оценить искажения, вносимые в решения за счет применяемых приемов регуляризации задачи.

Математическая модель

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

Решаемые дифференциальные уравнения записываются в следующем виде.

1) Уравнения баланса массы.

Для пара:

дх

(Aav pu) = mdv пш + mp П,

где x, м, — продольная координата; A, м2, — площадь поперечного сечения канала; av — объемная доля пара; pv, кг/м3, — его плотность; uv, м/с, — его скорость; mdv, mp, кг/(м2-с), — массовые источники генерации пара от капель и пленки соответственно (положительны при конденсации и отрицательны при испарении); n,d, nf м, — периметры поверхности теплообмена капель и пленки с паром, соответственно; индексы v, d, f относятся к пару (vapor), каплям (drops) и к жидкой пленке (film) соответственно. Далее будет использован еще индекс w, который относится к стенке канала (wall).

Для капель:

д

дХ (AadiPdUd ) = -mdvПШ - nf (Sd - Se ), (2)

где ad, pd — объемная доля и плотность капель, соответственно; ud, м/с, — их скорость; Se, Sd, кг/(м2-с), — скорости уноса и осаждения капель на поверхности жидкой пленки, соответственно.

Для жидкой пленки:

д

— (Aa f p fUf ) = -mfv nid + nf (Sd - Se X (3)

где af, pf — объемная доля и плотность жидкой пленки, соответственно; и, м/с, — скорость ее потока.

2) Уравнения баланса импульса.

Для пара:

д _ — дР

— (AavPvuv2) +avA д = mdvПШ (udl - uv ) +

дх дх (4)

+ mfv nf (ufi - uv ) - nf Tvf - nd Tvd + Aav Pvgx ,

где P, Па, — давление; ud,, u,, м/с, — скорость межфазной поверхности (пар — капли и пар — пленка соответственно); т, , Tvd, кг/(с2-м) — сдвиговые напряжения (пар — жидкая пленка и пар — капли соответственно); g, м/с2 — проекция вектора тяжести на ось х; верхняя черта указывает на усредненное значение величины.

Используемая скорость межфазной поверхности ud между каплями и паром считается равной скорости капель ud. Для скоро-

д

сти межфазной поверхности между жидкой пленкой и паром используется зависимость из кода CATHARE [9]:

а„

а

ufi =

а + а

-uf +

f

f

а + а,

-и,.

Для капель:

• _ — дР

-^(АасРdu2) + аА ^Т = -mvПИ (udi - uv ) +

(5)

^^ а и и ' и

ох ох

+ ПиТи + АааРи&х - ПГ (*аиа - 5еи/).

Для жидкой пленки:

д ,, 2Ч _ -лдР — (Аа , р гыг) + а ,А-=

дх ■>■>■> / дх

= -тРП/ (и/ - и) - Пк/т, + Пт, + (6) + Аа/Р/8х + П/ (¿Х - 5еиг),

где П, м — периметр поверхности теплообмена между стенкой канала и пленкой; х,, кг/(с2-м), — сдвиговое напряжение между жидкой пленкой и стенкой канала. 3) Уравнения баланса энергии. Для пара:

д_

дх (

(АауpvuvHv ) = (HTC)VdПи(Tsat - Tv) +

u2 ^

h + —

Vv 2 y

mvПи + (HTC)f nf (Tsat - Tv ) 4(7)

( u2 Л

J fi hv + —

v 2

v y

mv nf '

где Н, Дж/кг, — полная удельная энтальпия пара, Н = к + 0,5и2 (к — удельная энтальпия); (НТС)^ (НТС), Вт/(м2К -коэффициенты теплоотдачи от пара к межфазной поверхности с каплями и от пара к межфазной поверхности с жидкой пленкой, соответственно; Та(, Т, К, — температура насыщения и пара, соответственно.

Для капель:

_д_

дх

(ЛайpdudHd ) = (HTC)dVПи(Taat - Td ) -

f II 2 ^ h + ^

(8)

mdv Пи -nf (SdHd - SHf ),

где Ни, Н Дж/кг, — полная удельная энтальпия капель и жидкой пленки, Н/ = каг + 0,5и2 (ки/ — удельная энтальпия);

(НТС)Л, Вт/(м2-К), — коэффициенты теплоотдачи от капель к межфазной поверхности с паром; Ти, К, — температура капель. Для жидкой пленки:

_д_

дх

( Аа f р fufHf ) = (HTC) v Пг (Taat - Tf ) -

( u2 Л hf + f

f 2

V У

mVnf + nf (SdHd - SeHf ) +

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

(9)

qwfnfw 1

IwfiП fw ,

где (НТС/ Вт/(м^К), — коэффициент теплоотдачи от жидкой пленки к межфазной поверхности с паром; Т, К, — температура жидкой пленки; д, Вт/м2, — тепловой поток от стенки канала к жидкой пленке; д,, Вт/м2, — часть теплового потока от стенки канала, идущая непосредственно на генерацию пара.

Система также дополняется уравнением баланса фаз

Ъак =1 (10)

и уравнениями термодинамического состояния вида

(11) (12)

Р* = Р* (Р ,T ), e* = ek (Р ,T ),

где ек, Дж/кг, — внутренняя удельная энергия к-ой фазы.

Замыкающие соотношения

Для замыкания системы уравнений| (1) — (10) требуется набор замыкающих соотношений, описывающий процессы обмена массы, импульса и энергии как между фазами, так и со стенкой канала.

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

mdy = -[(HTC)dv (Taat - Td ) + + (HTC)vd (Taat -T )]/(h - hd ).

(13)

df

df

В дисперсно-кольцевом режиме течения теплообмен непосредственно со стенкой канала будет иметь только жидкая пленка. Соответственно, в правую часть ее урав-

нения баланса энергии добавляется член Ц/ П , описывающий поток тепла, получаемый жидкой пленкой от стенки канала. В уравнение баланса энергии также добавляется член, характеризующий часть теплового потока от стенки канала, идущая непосредственно на генерацию пара: П //, где ^ = V?/-; V = 0 -1.

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

тр = -[(НТС) д (Т^ - Тг) +

/у V аа1

+ (НТС)у/ (Таа( -Ту ) + Ц/, ]/(Ау - кг).

(14)

Для обоих источников массы (13), (14) слагаемые вида (НТС)^ ТаЛ — Тк) представляют собой тепловой поток от жидкой фазы к межфазной поверхности с паром, а (НТС)ук( Т аа1 — Ту) — тепловой поток от пара к межфазной поверхности с жидкой фазой.

Используемая модель фазового перехода допускает наличие разных по величине коэффициентов теплоотдачи с двух сторон от межфазной поверхности, однако в качестве упрощения используется следующее предположение. Считаем, что коэффициенты теплоотдачи по обе стороны от межфазной поверхности одинаковы и равны лимитирующему показателю — коэффициенту теплоотдачи со стороны пара. Тогда коэффициенты теплоотдачи от пара к межфазной поверхности с каплями и от капель к межфазной поверхности с паром равны и представляются в виде корреляции, описывающей обтекание сферических частиц потоком газа [10]:

(НТС)у, = (НТС)Л =

= Б- (2 + 0,5Ке»'5РгГ), (15) Бй

где число Рейнольдса для капель следует выражению

Р V

Яей =

К - иЛВс1

М*

(Б, м, — средний диаметр капель; , Пас, — коэффициент динамической вязкости пара); , Вт/(м-К), — коэффициент теплопроводности пара.

Коэффициенты теплоотдачи от пара к межфазной поверхности с жидкой пленкой и от жидкой пленки к межфазной поверхности с паром рассчитываются по корреляции для однофазной конвекции [11], применительно к парокапельному ядру:

(НТС)у/ = (НТС) /у = К. ______^ (16)

■ (0,023 Ке0'8 РГу0'4),

(Б - 25)

где Б, м, — внутренний диаметр канала; 5, м, — средняя толщина пленки; Ргу — число Прандтля для пара; Яе у — число Рей-нольдса для пара,

|«у - 1 (Б - 25)

Яе,. =

М*

Выражения для источников массы, описывающих гидродинамический унос и осаждение на поверхности жидкой пленки (£ и Sd), заимствованы из работы [2]. Автор работы считает, что доминантным механизмом осаждения капель на поверхность пленки является турбулентная диффузия:

Бл = 9 • 10-3«

( С V

Яе-0'2 8еу-2/3С, (17)

где С, кг/м3, — концентрация капель,

С =

ай Рй а + а,

__+__^ ^-у 1 ^

Ру«у Рй

8су — число Шмидта для пара (в качестве упрощения оно приравнивается единице);

здесь Яеу =

К

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

^ = 1,07

а

\к , если Яеу > 105;

/ \ 0,4

Р /

(18)

|Л:а[2,136^(Кеу) - 9,68], если < 105, где ^ = 0,575 + 21,73 • 10352 - 38,3 • 10653 +

+ 55,68 • 10984.

Представленные в уравнениях баланса импульса сдвиговые напряжения (х, — между стенкой канала и жидкой пленкой; х, — между паром и жидкой пленкой; х^ — между паром и каплями) записыва-

ются в следующем виде:

- rf f I |.

Twf - rJwf 2 uf \uf\;

V - f у-uf)|uv -uf|;

vd

- rfvd у (uv - ud1 uv - ud| ,

(19)

(20) (21)

где С/кк — соответствующие коэффициенты трения (силы сопротивления).

Коэффициент трения между паром и жидкой пленкой рассчитывается по модифицированной корреляции Уоллиса [12]:

f - М| f 1 + 300 !„ vf Re°'25 l D)

(22)

где

Re -

Pv |uv| (D - 28)

^v

Коэффициент аэродинамического сопротивления капель рассчитывается по зависимости [13]:

(

f -

n , 24 0,4 +-+

Л

Re,

VRe

(23)

где

Re d - max

0,1;

Pv Iuv - ud|Dd Л ^v

Коэффициент трения между жидкой пленкой и стенкой канала выражается следующим образом [14]:

rfjw - max

( 16 ; 0,079 ^

Re ;Re0,25 V R Dhf Re Dhf )

(24)

где

Re Dhf -

uf p^a ff Ну

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

При записи замыкающих соотношений следует также уделить внимание формулировке используемых геометрических харак-

теристик (представлены далее); к ним относятся:

периметры межфазных взаимодействий;

периметр взаимодействия жидкой пленки со стенкой канала;

средняя толщина жидкой пленки;

средний диаметр капель.

Периметр межфазной поверхности между паром и каплями следует выражению

6а.,

Пы -

DJ

(25)

Периметр межфазной поверхности между паром и жидкой пленкой:

nf -n(D - 28), (26)

где 8 - 0,5D(1 - ^ 1 - af) (средняя толщина жидкой пленки).

Периметр поверхности между жидкой пленкой и стенкой канала:

nwf - nD. (27)

Диаметр капли рассчитывается по методике, принятой в работе [14]:

Dd - max(8,4 • 10-5;min[Dd1; Dd2]), (28)

Dd1 - 7,96 • 10-

CT

PvJv

Re,

/ у/з 'Pd

/ \2/3 Ц.

^d

,(29)

где jv — приведенная скорость,

uvPv^a

Jv -

PvA

- uvav;

число Рейнольдса

Re -

Pv|Jv|D.

Dd2 - 0,254L

Hv

-0,13 Wev +

V'6 + (0,13 We„ )2

(30)

где Ь — характерный размер (взят как внутренний диаметр канала Б);

рЛ Ь

Wev -

CT

Численный метод

При переходе к конечно-разностной за-

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

(W л ГГ к,т ( АакРкик ^

\¥ = wk , = Аак Рки2 ; (31)

w V к ,е , Аак РкикНк

объединены. Поскольку вектор потоков не является линейным, с целью достижения сходимости процесса разбиваем искомые величины на следующем шаге на сумму известной величины Wn+1's с предыдущей итерации и их приращений на искомой

итерации АЖ"

что приводит рассма-

О =

к,е

триваемую систему уравнений к виду:

(32)

+1,5 + а ^

к, р

к, р

- WJn

к, р

где W — вектор потоков; Q — вектор правых частей; индекс к определяет принадлежность к соответствующей «жидкости» (пар, капля, жидкая пленка).

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

Ах

= ; р = т, е; (36)

Р" + 1'5 + а АI-

+ Аwkn;1'5+1 - W",

Ах

+ АР" + 1,5+1 _ рп ^

Ах )

(37)

_ +1,5 .

= ¿к ,1 ;

"+1,5 , А И+1,5+1

а + Аа-.

) = 1-

(38)

wи+1 - wn

-^ = ¿тр1; р = т, е;

Ах

Wkn+1 - Wkn. _ -( Р"+1 - Р" А ' ' ■ + а^4

Ах

Ах

+

(33)

п+и (34)

К а?1) = 1,

(35)

где индекс р обозначает принадлежность к определенному уравнению баланса (т — массы, I — импульса, е — энергии).

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

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

/ = К ,Тк, ик, Р )Т.

В результате система сводится к матричной записи вида

Мп+1.. А^п+1,.+1 = вn+1. *, (39)

которая решается методом Гаусса; Бп+М —

Таблица

Коэффициенты матрицы преобразования М"+1'*

к

к

Уравнение баланса Аа к АТк Аик АР

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

массы Ари Ааи дТ Аар Ааи — дР

импульса А ри2 Ааи2 — дТ Аар2и Ааи2 — + Аа дР

энергии АриН АаиЕ + дТ де +Аари- дТ А ар(( + и2) Ааи(Е — + р — + 1 | 1, дР дР )

фаз 1 0 0 0

вектор правых частей; Мп+1'5 — матрица преобразования вектора потоков к вектору примитивных переменных.

Коэффициенты матрицы Мп+1-5 можно получить путем векторного дифференцирования вектора потоков по вектору примитивных переменных:

У+1'5

8+1 = I I Д/"+1' '+1 = ^ ^ )

= ЫП+1'8 Д/П+1'5+1

Матрица Мп+15 имеет повторяющуюся блочную структуру с размером блоков 3 х 3 (см. таблицу).-

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

Валидация модели

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

адиабатических круглых трубах с внутренним диаметром 10 и 20 мм, с диапазоном давлений 3 — 9 МПа и общих расходов 500 — 3000 кг/(м2^с).

С помощью модели, представленной в настоящей работе, было рассчитано в общей сложности 90 экспериментальных точек. Сравнение рассчитанных зависимостей падения полного давления ар/аХ от относительного расхода пара с измеренными представлено на рис. 1. — это отношение абсолютного расхода пара к сумме расходов всех трех «жидкостей».

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

Рис. 1. Экспериментальные (символы) и расчетные (линии) зависимости падения полного давления в канале от относительного расхода пара для внутренних диаметров канала Б = 10 мм (а) и 20 мм (Ь), при различных значениях общего расхода, кг/(м2с): 500 (1), 750 (2), 1000 (3), 2000 (4), 3000 (5). Давление Р = 7 МПа

Рис. 2. Связь рассчитанного падения полного давления в канале с измеренными значениями для всех использованных экспериментальных данных работы [15]

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

Рис. 3. Зависимости, аналогичные представленным на рис. 1, которые получены для различных коэффициентов шероховатости: 1,00 (1), 1,20 (2), 1,50 (3), 1,85 (4); Р = 3 МПа, Б = 10 мм

Однако можно показать, что учет дополнительных факторов, например шероховатости трубы (который усложнит выражение для коэффициента трения между пленкой и стенкой канала), позволяет приблизить результаты расчета к экспериментально измеренным значениям. Результат такого уточнения расчетов представлен на рис. 3. Для простоты шероховатость трубы учитывается в коэффициенте трения вариацией множителя, а расчеты проведены только для одной серии экспериментов (давление — 3 МПа, общий расход — 1000 кг/(м2-с), внутренний диаметр канала — 10 мм). Данные рис. 3 наглядно демонстрируют, что путем увеличения коэффициента шероховатости можно приблизить расчетное значение падения полного давления к экспериментально измеренному.

Заключение

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

Выполнено начальное тестирование разработанной вычислительной модели посредством сопоставления результатов расчетов с экспериментальными данными Дж. Вюрца [15]. Сравнение показало, что представленная в настоящей работе трех-жидкостная модель хорошо описывает рассматриваемую серию экспериментов с парожидкостным потоком и качественно полученные зависимости полностью повторяют экспериментальные. При этом существуют количественные расхождения, не превышающие 20 %. Однако показано, что указанные расхождения можно снижать путем корректировки замыкающих соотношений.

список литературы

1. Sami S.M. An improved numerical model for annular two-phase flow with liquid entrainment // Int. Comm. Heat Mass Transfer. 1988. Vol. 15. No. 1. Pp. 281-292.

2. Sugowara S. Droplet deposition and entrainment modeling based on the three-fluid model // Nuclear Engineering and Design. 1990. Vol. 122. No. 1-3. Pp. 62-84.

3. Hoyer N. Calculations of dry-out and postdry-out heat transfer for tube geometry // International Journal of Multiphase Flow. 1997. Vol. 24. No. 2. Pp. 319-334.

4. Алипченков B.M., Зайчик Л.И., Зейгар-ник Ю.А., Соловьев С.Л., Стоник О.Г. Развитие трехжидкостной модели двухфазного потока для дисперсно-кольцевого режима течения в каналах. Размер капель// Теплофизика высоких температур. 2002. Т. 40. № 4. C. 641-651.

5. Jayanti S., Valette M. Prediction of dry-out and postdry-out heat transfer at high pressures using one-dimensional three-fluid model // International Journal of Heat and Mass Transfer. 2004. Vol. 47. No. 22. Pp. 4895-4910.

6. Stevanovic V., Stanojevic V., Radic M.D. Three-fluid model predictions of pressure changes in condensing vertical tubes // International Journal of Heat and Mass Transfer. 2008. Vol. 51. No. 15-16. Pp. 3736-3744.

7. Safari H., Dalir N. Effect of virtual mass force on prediction of pressure changes in condensing tubes // Thermal Science. 2012. Vol. 16. No. 2.

Pp. 613-622.

8. Bulovich S.V., Smirnov E.M. Experience in using a numerical scheme with artificial viscosity at solving the Riemann problem for a multi-fluid model of multiphase flow // AIP Conference Proceedings. 2018. Vol. 1959. Pp. 050007-1-050007-8.

9. Bestion D. The physical closure laws in the CATHARE code // Nuclear Engineering and Design. 1990. Vol. 124. No. 3. Pp. 229-245.

10. Броунштейн Б.И., Бишбейн Г.А. Гидродинамика, массо- и теплообмен в дисперсных системах. Ленинград: Химия, 1977. 280 с.

11. Кутателадзе С.С., Боришанский В.М. Справочник по теплопередаче. М.: Госэнерго-издат, 1958. 414 с.

12. Wallis G.B. One-dimensional two-phase flow. New York: McGraw-Hill, 1969. 408 p.

13. Нигматулин Р.И. Динамика многофазных сред. В 2-х ч. Ч. 1. М.: Наука. Гл. ред. физ.-мат. лит-ры, 1987. 464 с.

14. Юдов Ю.В., Волкова С.Н., Слуцкий А.А. КОРСАР/В1.1. Теплогидравлический расчетный код. Методика расчета замыкающих соотношений и отдельных физических явлений контурной теплогидравлики. Сосновый Бор, 2001. 147 с.

15. Wurtz J. An experimental and theoretical investigation of annular steam-water flow in tubes and annuli at 30 to 90 bar. Roskilde, Denmark: Technical University of Denmark. Ris. Report. 1978. No. 372. 141 p.

Статья поступила в редакцию 13.07.2018, принята к публикации 07.08.2018.

сведения об авторах

АВДЕЕВ Евгений Эдуардович — студент Института прикладной математики и механики Санкт-Петербургского политехнического университета Петра Великого, Санкт-Петербург, Российская Федерация.

195251, Российская Федерация, г. Санкт-Петербург, Политехническая ул., 29 avdeev-evgeni@yandex.ru

ПЛЕТНЕВ Александр Александрович — кандидат технических наук, доцент кафедры гидродинамики, горения и теплообмена Санкт-Петербургского политехнического университета Петра Великого, Санкт-Петербург, Российская Федерация.

195251, Российская Федерация, г. Санкт-Петербург, Политехническая ул., 29 aapletnev@yandex.ru

БУЛОВИЧ Сергей Валерьевич — кандидат физико-математических наук, доцент кафедры гидродинамики, горения и теплообмена Санкт-Петербургского политехнического университета Петра Великого, Санкт-Петербург, Российская Федерация.

195251, Российская Федерация, г. Санкт-Петербург, Политехническая ул., 29 bulovic@yandex.ru

references

[1] S.M. Sami, An improved numerical model for annular two-phase flow with liquid entrainment, Int. Comm. Heat mass transfer. 15 (1) (1988) 281-292.

[2] S. Sugowara, Droplet deposition and entrainment modeling based on the three-fluid model, Nuclear Engineering and Design. 122 (1-3) (1990) 62-84.

[3] N. Hoyer, Calculations of dry-out and postdry-out heat transfer for tube geometry, Int. J. Multiphase Flow. 24 (2) (1997) 319-334.

[4] V.M. Alipchenkov, L.I. Zaichik, Yu.A. Zeigarnik, et al., The development of a three-fluid model of two-phase flow for a dispersed-annular mode of flow in channels: size of droplets, Heat and Mass Transfer and Physical Gasdynamics. 40 (4) (2002) 641-651.

[5] S. Jayanti, M. Valette, Prediction of dry-out and postdry-out heat transfer at high pressures using one-dimensional three-fluid model, Int. J. Heat and Mass Transfer. 47 (22) (2004) 4895-4910.

[6] V. Stevanovic, M. Stanojevic, D. Radic, Three-fluid model predictions of pressure changes in condensing vertical tubes, Int. J. Heat and Mass Transfer. 51 (15-16) (2008) 3736-3744.

[7] H. Saffari, N. Dalir, Effect of virtual mass force on prediction of pressure changes in condensing tubes, Thermal Science. 16 (2) (2012) 613-622.

[8] S.V. Bulovich, E.M. Smirnov, Experience in using a numerical scheme with artificial viscosity at solving the Riemann problem for a multi-fluid model

Received 13.07.2018, accepted 07.08.2018.

of multiphase flow, AIP Conference Proceedings. 1959 (2018) 050007-1-050007-8.

[9] D. Bestion, The physical closure laws in the CATHARE code, Nuclear Engineering and Design. 124 (3) (1990) 229-245.

[10] B.I. Brounshteyn, G.A. Bishbeyn, Gidrodinamika, masso- i teploobmen v dispersnykh sistemakh [Hydrodynamics, mass and heat exchange in the dispersion systems], Leningrad, Chemistry PH, 1977.

[11] S.S. Kutateladze, V.M. Borishanskiy,

Spravochnik po teploperedache [Manual of heat transfer], The State Energetic PH, Moscow, 1958.

[12] G.B. Wallis, One-dimensional two-phase flow, McGraw-Hill, NewYork, 1969.

[13] R.I. Nigmatulin, Dinamika mnogofaznykh sred v 2 ch. Ch.1. [Multiphase mixture dynamics, in 2 parts, P. 1], Moscow, Nauka, 1987.

[14] Yu.V. Yudov, S.N. Volkova, A.A. Slutskiy, KORSAR/V1.1 teplogidravlicheskiy raschetnyy kod, metodika rascheta zamykayushchikh sootnosheniy i otdelnykh fizicheskikh yavleniy konturnoy teplogidravliki [KORSAR/V1.1 thermal and hydraulic calculation code, design procedure of constitutive relationships and separate physical phenomena of contour heathydraulics ], Sosnoviy Bor, 2001.

[15] J. Würtz, An experimental and theoretical investigation of annular steam-water flow in tubes and annuli at 30 to 90 bar, No. 372, Technical university of Denmark, Roskilde, 1978, 141 p.

the authors

AVDEEV Evgeniy E.

Peter the Great St. Petersburg Polytechnic University

29 Politechnicheskaya St., St. Petersburg, 195251, Russian Federation

avdeev-evgeni@yandex.ru

PLETNEV Alexander A.

Peter the Great St. Petersburg Polytechnic University

29 Politechnicheskaya St., St. Petersburg, 195251, Russian Federation

aapletnev@yandex.ru

bulovich Sergey V.

Peter the Great St. Petersburg Polytechnic University

29 Politechnicheskaya St., St. Petersburg, 195251, Russian Federation

bulovic@yandex.ru

© Санкт-Петербургский политехнический университет Петра Великого, 2018

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