Приложения дифференциальных уравнений
УДК 532. 534. 539 А.К. Кубанова
МОДЕЛИРОВАНИЕ НЕСТАЦИОНАРНОГО ТЕЧЕНИЯ ГАЗА В ДВУХФАЗНОЙ СРЕДЕ
Приведено решение нелинейной задачи об истечении газа из пористой среды за поршень, который движется с постоянной скоростью ¥0 , для двухфазной среды, в которой газ и твердые частицы закреплены и образуют скелет.
Рассматривается решение нелинейной задачи об истечении газа из пористой среды за поршень, который движется с постоянной скоростью У0. Среда двухфазная - газ и твердые частицы закреплены, образуют скелет.
Выведем уравнение движения газа, считая это движение одномерным; движущийся газ -совершенным, не обладающим теплопроводностью. Исследование будем проводить в координатах Лагранжа. Ось координат направим по движению поршня.
Таким образом, движущийся за поршнем газ будет находиться на отрицательной половине оси х. Введем перемещение сечения газа V = V(х, Ґ). В дальнейшем мы будем использовать
также обозначения —— = —х; —— = —. Физическая картина движения имеет вид, представлен-дх Э/
ный на рис. 1.
X = 0
1
о N * N 1 ^
х = —а0і
х
Р и с. 1. Схема течения газа
Область х = У0Ґ - область "чистого" газа, прилегающая к поршню (на рис. 1 обозначена цифрой 1). В этой области уравнение движения газа имеет вид [1]:
=-
—х
(1 + их)г+1 “ ’
где а0 - скорость звука в невозмущенном газе.
Область х = —а0, - область распространения возмущения в двухфазной среде (на рис. 1. обозначена цифрой 2).
Получим уравнение движения для области (2):
и = -—^ - Ш,, р = Ро(1 + их )-1, Р = Р0(Р . Р-1У; а, р0 ах
Ро = /0 Ро0; р = /о р 0. (1)
где к - коэффициент взаимодействия; р0, р - приведенные плотности пористой среды (нижний
индекс "0" соответствует невозмущенному состоянию), определяемые через пористость среды
0
/0 и истинную плотность р .
Из (1) получим
2
э 2и
э 2и
Э,2 (1 + их )7+1 Эх
- ки.
где ао - скорость звука в двухфазной среде [2], или сокращенно
Начальные условия:
- а'и-)^ -ки-
их = и, = 0 на х = -аЛ.
(3)
(4)
Граничные условия имеют вид
и, = У0 при х = Уд, (на поршне), их - неизвестна.
Поршень находится на оси 0, . В области 1 и 2 задача решается численно методом характеристик; на плоскости хґ имеются два характеристических направления:
Жх ,тт ,
— = ±а(и ), С, v
(5)
которые имеют один и тот же вид для областей 2 и 1.
Условия на характеристиках: для области 2 -
сіи, = ±а(их)СЮх - ки,Ж; (6)
для области 1 -
Жи, =±а(их )сЮх. (7)
Обозначим границу двух сред 1 и 2 через 1(ґ). Начало координат является особой точкой, а
именно, точкой пересечения характеристик и точкой разрыва величин их и и,. Значения по-
следних будут равны нулю, если приближаться к началу координат вдоль оси 0х , и равны соответственно
и =■ 2а"
(7 -1)
1-7
1 - (1 + их) 2
и, = Уп, их =
1 - V (7 -1)
2а0
2
1-7
-1,
(8)
если приближаться к началу по оси 0/. Действительно, из условия на характеристике (7) для газовой фазы получим
| а(их )Сих
2а 0
- Ї-1
1 - (их +1) 2
0 (г -1)
Из последнего выражения определим деформацию их. Система характеристик, исходящих из одной точки, называется центрированными волнами Римана. В области 2 имеем центрированные волны для газовой фазы, границу которых определим таким образом, что на граничном луче при подходе к точке 0 со стороны оси 0, выполняется условие (8). Тангенс угла наклона граничного луча с осью , определится так:
агр = а0
1-
V (7 -1)
2а0
С =
7+1
7 -1
(9)
Задаем равномерный шаг А Ь в виде АЬ =
а0 - агр
N
где N - число "шагов". Шаг выбираем
так, чтобы А Ь составлял примерно 0,5% от а0 . Определяем соответствующие этому значению деформации и скорости (при подходе к началу координат по лучам):
их =
1 -
АЬп
7+1
2а0
1-7
1 - (1 + их) 2
(10)
(7 -1)
V у
где п = 1,2,3 ... N.
После определения этих величин в точке 0, решается смешанная задача методом характеристик. Г раницей области, откуда начинаем решение, является прямая х = —а0,. Этот результат имеет простой физический смысл, а именно: возникающие на поршне возмущения достигают через время , = -х/а0 частицу с координатой х(х<0). Иначе говоря, возмущение распространяется со скоростью звука а0 в двухфазной среде.
Таким образом, начальные условия переносятся на характеристику х = —а0і; их=и,= 0 и
2
дополняются граничными условиями и, = У0, их - неизвестно при X = 0.
Изложим основные расчетные соотношения.
Значение параметров на слоях (т = 1, 2, 3 ... ) в точках п = 0 известны (начальные значения на х = -а01). Задаем шаг Д, на оси времени и проводим прямую, параллельную оси 0х до пересечения с х = -а0,. Получаем узел (1,0), из которого проводим положительную характеристику с угловым коэффициентом, взятым в узле (1,0), до пересечения с лучом ДЬ из начала координат 0. Определяем узел (1, 1) и т.д. до узла (1, Ы).
Последовательность расчетов представлена на рис. 2 приведенными ниже формулами:
*1п =[ а (Цх1, п-1) • *1, п-1 - х1, п-1 ]Х[ а (ихп ) + а (их1, п-1) ] ;
х1п = а(иха1 ) [ х1, п-1 - а(их1,п-1) • ,1, п -1 ]х[ а(ихп ) + а(их1,п-1) ] (п = 1, 2 3 . );
У+1 у+1
а(иХа1) = а0 (1 + ихп ) 2 ; а(Ц^„-1) = а0 (1 + ^ );
В = и,п (1 - кгы) + а Цп К; А = иПп п-1 [1 - к (^ - ^ п-1)]-а (и^
Цхп = . В - (--------г; иЧп = [ а (Ц) А + а(их1,п-1)в]х[а(Цхоп) + а^-О ]-1;
аК) + а(х1,п-1)
и* = 0,5(и,п + и, 1п-1); и1 = 0,5(иЖ1,п-1 -иХы;
^1п = ^1,п - ^1,п-1; х1п = х1,п - х1,п-1;
N
и1, п+1 = -х1, п +1, и1 = ^и1п .
п=1
Начиная со слоев т = 2, 3, имеем:
хт, п хт ,п—1 а (ихт ,п-1 )(,т, п ,т, п-1);
хт,п хт-1,п а(ихт-1п )(,т,п ,т-1,п );
^п - ^„ = а(ихт^ ) [ихт_ - ихт^ ]-
- ки, (, -, .);
,т ,п-1 ' тп т, п-1/’
Ли^-У0 = а{и^+2){и^+1-и^+2).
Смещение частицы в лагранжевых коорди- Р и с. 2. Расчетная схема
натах определяем по формулам:
йи, = а (Ц) (№х; хт и = хтп - хт,п-1; = 0,5^ - и,т^); и^ = 0,5^ - ихт,-,);
N
итп и 1тп ^ тп + ихтп хтп ; Цт Цтп *
п=1
Из приведенных соотношений определяем параметры Цх, и,, х, , узла (тп), затем плотность и давление в этом узле по формулам:
р1„ = р0°(1 + иХп )-1; Рпт = р
р0 V
г тп
р0
На границе раздела пористой среды и "чистого" газа (х = 1(,)) имеют место следующие соотношения:
1) рХ = р 0 /у2;
2) рХ (У1-У2) = /0Р2 - Р + Р', Р = (1 - ЛХР - Р);
3) Р10х7
—2_+.
Р Р
Г1 г2
р10 Р 2
= /Р2X2 - ру,,
(11)
2 2 у -1
которые дополняются уравнением состояния газов в пористой области, интегралом Римана вдоль центрированной волны и условием для скорости частиц газа на поршне:
2а0
4) Р2 = р
р2 у
чРо,
; 5) Х2 =
у -1
1 -
у -1 Р[ ^ чР 00.
(12)
Из уравнений (11) и (12) получим следующее выражение:
У,2-Г,2 у
——^ + '
у -1
Рг - Р - УУЫ
Р0(2 - /0) р0 (2 - /0)
или
У02у
+
У0у у2
V;
Рг
2 (2 - /°)(у -1) (2 - /°)(у -1) 2 /°(2 - /°)(у - 1)р°Т; I р°
(У -1) Р0
^41 = 0.(13)
Р° )
Заменяем отношение плотностей (р°/ р0’)у 1 через скорость Х2, которая задается выраже-
нием
V; (У - 1)
Р2_ Р 0
У -1
2а0
Тогда для скорости У2 имеем кубическое уравнение
а¥23 + ЬГ22 + с¥2 + й = 0, где а, Ь, с, й - коэффициенты, зависящие от постоянных Ро, /°, у, У°, р°д, а° :
(14)
0,5 +
/4(7 -1)2
4а0
0,5 V;
^2У
(2 - / )(у -1)’
/з(у-1)2 , /4(7-1)
/2 —
Ку
(2 - /0)(у -1)
/з(У-1) , _У^РІУ-!)"
/0 Ро(2 - /0)
4а
’ 4 (У -1)Р0
Из решения кубического уравнения установим связь между скоростями У1 и У2, затем между плотностями и давлениями на границе раздела двухфазной пористой среды и "чистого" газа.
Для расчета параметров на границе имеем формулы (см. рис. 2):
и1, п+1 — х1,п+1? •Х1
,п+1 1,п +2
'а(их\,п+2 ) [*1,п+1 *1,п+2 ] ? х1,п+1 Х\,Ы — а(их1,N ) ё *1,п+1 *1^ ] ?
и (2) _ и
п+1 N
г(1)
а^N ) Х (их1,п+1 иъ,N ) ^иі1,N (*1,п+1 *1,N )
/0иг1;)п+1 - V; — а (и^п+2)(ия., п+1 - их., п+2 ) 1
и
ия
' — /слои 2слоїГ~^^\
Для области (1) условия на характеристиках не содержат члена Шtdt и имеют вид: dUt — а(их)dUx. Из этих формул определяем искомые величины деформации их и скоростей ии На рис. 3 приведен расчет положения границы х — / (*).
05
ал
Р и с. 3. Расчетная величина границы х — / (*)
БИБЛИОГРАФИЧЕСКИМ СПИСОК
1. Рахматулин Х.А. Газовая и волновая динамика. М.: МГУ, 1983. 280 с.
2. Кубанова А.К. О распространении одномерной плоской волны в двухфазной среде. // Проблемы динамики взаимодействия деформируемых сред. Ереван, ЕрГУ. 1987.
Поступила 10.05.2003,
а
х