Модификация метода крупных частиц для исследования течений газовзвесей

Отсутствие в природе чистых веществ требует активного развития математических моделей многокомпонентных сред, достоверно описывающих физические процессы, применяемые в различных отраслях науки и техники, с одной стороны. С другой стороны, развитие вычислительной техники позволяет получать решения для новых [1], все более сложных математических моделей многокомпонентных сред. Более того, есть такие проблемы, когда математическое моделирование является единственным средством предварительного изучения явлений (например, [2]). Адекватность математических моделей многокомпонентных сред физическим процессам предъявляет достаточно жесткие требования к математическим моделям: с одной стороны, уравнения сохранения должны быть инвариантны относительно преобразования Галилея [3], с другой стороны, должны выполняться законы сохранения для смеси [4]. В работах [1, 5] было показано, каким образом можно выполнить оба эти условия.

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

1.    Математическая модель газовзвеси

Рассмотрим одномерный плоский случай математической модели течения газа с твердыми частицами (аэровзвесь), которая описывается системой уравнений сохранения [5]. Данная

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ система уравнений двухфазной аэровзвеси [5] без химических превращений имеет следующий вид

dp 1 , dp 1 v 1

∂t ∂x d1v1

^{ρ 1} dt

∂ρ 2    ∂ρ 2 v 2

, dt ⁺ dx ’

-

∂p         d 2 v 2

∂x         dt

∂n ∂nv

+      = ⁰ ’

∂t ∂x

• -a 2 yp + nf, ∂x

(1 - 1)

(1 - 2)

d 1 e 1 ^{ρ 1} dt

Р^а 1 d 1 Pl

⁺ nf ^{( v 1} —^{v 2)}

d 2 e 2 pa 2 d 2 p 2

• ^{P 2} IT =--

• _-

nq,

(1 - 3)

p 2 dt

+ nq.

(1 - 4)

p = p 1(P1 ,T1) = p2(P2,T2), e 1 = e 1(P1 ,T1),e2 = e2(P2,T2), v2

(1 - 5)

P 1 = P° 1 a 1 ,p 2 = p° 2 a 2 , a 1 + a 2 = 1 ,E i = e i + у ( i = 1 , 2) , f = nd ² p° C d ( v 1 - v 2 ) |v 1 - v 2 1/ 8 , q = ndA 1 Nu ( T 1 - T 2 ) -

(1 - 6)

Система уравнений (1.1) - (1.6) замыкается уравнениями состояния газовой фазы и частиц

e 1 = c v 1 ⁽ T 1 ^- T 0 ) + C 0 , e 1 =

( Г- Л р.e 2 = c 2 ( ^{Т 2 -} Т ⁰⁾ ,

(1 - 7)

Здесь индексы 1, 2 относятся соответственно к газу и частицам; p°,a i ( i = 1 , 2) - истинные плотности и объемные содержания фаз; p i ,v i ,T i ,e i , E i - парциальная плотность, скорость, температура, внутренняя и полная энергия i-й фазы; p - давление, n - число частиц в единице объема смеси; c v 1 и c 2 - теплоемкости фаз: C о - постоянная для нормирования внутренней энергии газовой фазы: А 1 - теплопроводноеть газовой фазы; R 1 - универсальная газовая постоянная; C d и Nu - коэффициент трения и число Нуссельта, определяемые числами Рейнольдса ( Re) и Прандтля (Pr) относительного движения фаз соответственно: k - показатель адиабаты Пуассона; d - диаметр частиц.

Уравнения (1.1) - уравнения неразрывности газа и частиц и уравнение сохранения числа частиц в единице объема смеси; (1.2) - уравнения импульса газа и частиц; (1.3) и (1.4) -уравнения сохранения внутренней энергии газа и частиц соответственно; (1.6) - уравнения, определяющие члены теплового ( q ) и силового ( f ) взаимодействия между фазами: (1.7) -уравнения состояния фаз.

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

Умножая уравнение сохранения импульса газовой фазы на v 1, а уравнение сохранения импульса конденсированной фазы на v ₂, получим уравнения сохранения кинетической энергии газа и частиц соответственно

∂ρ 1 v 1

^{v 1[} ~вГ ⁺

^дР 1 ^{v 2} ] = ∂x

∂p

-v 1 α 1     - nfv 1 ,

∂x

∂ρ2v2   ∂ρ2v2          ∂p v2[ Щ +      ] = -v2a2   - nfv2,

∂t ∂x          ∂x которые после простых преобразований принимают следующий вид

2 ∂ρ 1 ^v ₂ ¹

∂t

др 1 v 1 v 2 1 ∂x

∂p

-α 1 v 1    - nfv 1 ,

∂x

(1 - 8)

Д.С. Грищенко, Ю.М. Ковалев, Е.А. Ковалева

∂ρ 2 ^v 2 ^{2 2}

дt ⁺

∂ρ 2 v 2 ^v ₂ ²²

∂x

-

∂p α 2 v 2    - nfv 2 .

∂x

(1 - 9)

Преобразуем левые части уравнений сохранения внутренней энергии газа (1.3) и частиц (1.4) к дивергентному виду. С учетом равенств (1.1) они могут быть представлены в виде др 1 e 1 , др 1 e 1 v 1

∂t ∂x

pa 1 d 1 р °

р1

∂t

+ nf ( v 1 - v 2 ) - nq,

(1 - 10)

др 2 e 2 ₊ др 2 e 2 v 2

∂t ∂x

pa 2 d 2 р₂

= ~ 2 + ⁺ nq-р2 ^dt

(1 - 11)

Из уравнений неразрывности газовой и конденсированной фаз (1.1) легко получить следующие равенства

α 1

α 2

d 1 р 1 dt d 2 р 2 dt

-

-

∂α ∂α v р^ ST + da 2 da 2 v 2

^{р 2(} ГТ ⁺ ГГ

) ,

) -

Подставляя данные

выражения в уравнения (1.10) и (1.11) соответственно, получим

∂ρ 1 e 1 ∂t

+ др 1 e 1 v 1 ∂x

∂α 1 ∂α 1 v 1

^{-p (} S ⁺ ~^{) + nf ( v 1} - v ²⁾

-

nq,

(1 - 12)

др 2 e 2    др 2 e 2 v 2

∂t ∂x

∂α 2    ∂α 2 v 2

^{-p (} S ⁺ —^{) + nq}-

(1 - 13)

В случае несжимаемости конденсированной фазы уравнения сохранения внутренней энергии газовой (1.3) и конденсированной (1.4) фаз, легко преобразуются к виду

др 1 e 1 , др 1 e 1 v 1

∂t ∂x

∂α 1 v 1 ∂α 2 v 2

^{-p (} ST ⁺ ~^{) + nf ( v 1} - ^{v 2)} - ^nq-

(1 - 14)

∂ρ 2 e 2 ∂ρ 2 e 2 v 2

— ⁺ -gS = ^nq-

(1 - 15)

Для получения уравнения сохранения полной энергии смеси просуммируем левые и правые части уравнений (1.8), (1.9), (1.14), (1.15). В результате получим уравнение сохранения полной энергии смеси в виде

д ( р i E 1 + р 2 E 2 )

∂t

∂

+ тН р 1 v 1 E 1 + р 2 v 2 E 2 + ( a 1 v 1 + a 2 v 2 ) p ] = 0 -∂x

(1 - 16)

Система уравнений (1.1), (1.2),(1.5) - (1.7), (1.14) - (1.16) представляет собой замкнутую систему уравнений для описания течений газовзвесей, инвариантную относителвно преобразования Галилея.

2.    Модификация метода крупных частиц для расчета течений газовзвеси идеологией метода крупных частиц [6] систему законов сохранения , (1.5) - (1.7), (1.14) - (1.16) на эйлеровом этапе можно представить

В соответствии с газовзвеси (1.1), (1.2) следующим образом

д е=о =о %=о ,

(2 - 1)

∂v 1

^{р 1} д^ =

-

∂p

α 1

∂x

-

∂v2       ∂p f2 at = -a2 дХ + f

(2 - 2)

МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ

∂e 1         ∂α 1 v 1 ∂α 2 v 2 P 1 dt = p ( dx 1 ax 1 + nf ( v 1 v 2) nq’	(2 . 3)
∂e 2 P 2 -di = nq,	(2 . 4)
(              ) + я [( а 1 v 1 + а 2 v 2 ) p ] = 0 . ∂t         ∂x	(2 . 5)

Учитывая несжимаемость конденсированной фазы р2 = const, запишем уравнения (2.1), (2.3), (2.5) в более удобном для представления на эйлеровом этапе виде дР1   п ◦ да2 _ n dn _ о ◦ да 1 _ о

(2 . 6)

^{а 1} dt    ^{0,P 2} dt ⁰’ dt ^{0,P 1} dt ⁰

∂e 1

ρ 1 _∂t

-p ( а 1 ^ + а 2^) + nf ( v 1 - v 2 ) - nq, ∂x ∂x

(2 . 7)

P1   + P2   + a 1 ЛИ + а 2 ^ =0 ■(2

∂t ∂t ∂x∂x

Подставляя уравнение состояния газовой фазы (1,7) в уравнение (2.7) получим следующее базовое соотношение для определения давления на эйлеровом этапе

(2 . 9)

at = -     ^{p ( а 1} -X ^{+ а 2} dX )⁺ Х О Т^{( nf ( v 1} —^{v 2) -} nq ) .

Используя явные разностные представления для равенства (2.9), легко получить выражения для определения предварительных значений давления на новом m + 1 временном слое на границах i — 1 / 2 и i + 1 / 2 для ячеек i — 1, ini + 1

m+1 __ +1+1 + pi (k — 1) m m m m m m A t pi+1 / 2 =      5     (1 — Cm     ( ° 1,1+1 / 2( V1,1+1 — V1,1) + ° 2,1+1 / 2( V 2,1+1 — V 2,1 ^АД) —

^{2           а} 1 ,1 +1 / 2                                                    ^{A X}

• — ( ( Л ( n ”1 / 2 9 1 +1 / 2 )A t +        ( n ”1 / 2 f +1 / 2 ( v 1 m +1 / 2 - v m +1 / 2 )A t ) .        (2 . 10)

^а 1 ,1 +1 / 2                         ^а 1 ,1 +1 / 2

Здесь A t - шаг по времени, A x - шаг по пространству. Полученные значения давления используются для определения промежуточных величин скоростей на эйлеровом этапе:

m vm+=v m1 - m (pm++/2 - p m+12)^ - -m fmA t,             (2.11)

^P 1 ,1                    ^{A x   P} 1 ,1

m

• v m+ = v m - m ( p m+ ⁺ 1 / 2 - p m_ + 12)^ - -m f m ^a t.             (2 . 12)

^P 2 ,1                    ^{A x   P} 2 ,1

Для получения промежуточных значений скоростей газовой и конденсированной фаз можно использовать еще одну модификацию эйлерова этапа метода крупных частиц, связанную с частично неявной аппроксимацией силы межфазного взаимодействия. В этом случае равенства (2.11) и (2.12) можно представить в виде ^m

= v m - m ( p + / 2 - p m + / 2 )    - -m ( nd ² p 1 C d lv 1 - v 2 1/ 8) m ( Щ+ ¹ - v m )A t

^P 1 ,1                     ^{A x    P} 1 ,1

m v^ = vm -    (p 12 - p“У/2)^ + -m(nd2p1 Civ 1 - v2\/8)”(vm - v^^At.

^P 2 ,1                     ^{A x P} 2 ,1

Д.С. Грищенко, Ю.М. Ковалев, Е.А. Ковалева

Из полученных уравнений промежуточные значения скоростей легко определяются в явном виде mm+1    m  1_ (mm +1   ~m+1 ^t v 1 ,i = (v 1 ,i   p1 m (Pi+1 /2 Pi-1 /2) Ax +

(2 . 11) '

(2 . 12) ‘

^m                                                 ^m

+n^ (nd2 p1 civ 1 - v 21/8) mvmi at)/ (i +    (nd2 P1 cd\v 1 - v 21/8) mat), p 1 ,i                                                      p 1 ,i

~^{m +1} m 1 (mm ⁺¹ _ ~m ^{+1 A} t .

v 2 ,i =⁽ v 2 ,i p ^{( p} i +1 / 2 ^p i - 1 / 2 )A x ⁺

^m                                             '                ^m

+ n m ( nd 2 p 1 c_d\v 1 - v 2 1/ 8) m v mi a t ) / (1+( nd ² P 1 c_d\v 1 - v 2 1/ 8) m a t ) .

P 2 ,i                                                       P 2 ,i

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

v 1 ,i +1 / 2 = ⁽ v 1 ,i ⁺ v 1 ,i +1 ) / ² ’v 2 ,i +1 / 2 = ⁽ v 2 ,i ⁺ v 2 ,i +1 ) / ² •                  ⁽² • ¹³⁾

Теперв можно определить промежуточные значения внутренней энергии конденсированной фазы tm+ = emi +   ^qm At                        (2 • 14)

P 2 ,i и полной энергии смеси

„m ^ m +1    m fr m +1 _ m mm   m mm    n     m +1   m +1

^P 1 ,i ^E 1 ,i ^{+ P} 2 ,i ^E 2,i = P 1 ,i ^E 1 ,i ⁺ P 2 ,i ^E 2,i ~ ^{( a} 1 ,i +1 / 2 ^v 1 ,i +1 / 2 ^Pi +1 /2

• ⁿ ~m + ¹ m + ^{1 A} t n ~m + ¹ m +¹ _ n m + ¹ m +¹ ^      ⁽² • ¹⁵⁾

• ^a 1 ,i - 1 / 2 v 1 ,i - 1 / 2 ^P i - 1 / 2 ) A x ^{( a} 2 ,i +1 / 2 v 2 ,i +1 / 2 ^P i +1 / 2 ^a 2 ,i - 1 / 2 v 2 ,i - 1 / 2 ^P i - 1 / 2 ) A x'

На этапе Лагранжа и заключительном этапе метода крупных частиц для каждой фазы были использованы формулы, приведенные в монографии О.М. Белоцерковского и Ю.М. Давыдова [6].

Заключение

• 1.    Тестирование предложенной модификации метода крупных частиц проводилось на решении задач о распространении ударных волн в «замороженной» газовзвеси [8] и в облаке газовзвеси [9].

• 2.    Было показано, что применение на этапе Эйлера уравнений (2.10) - (2.15) более эффективно, чем применение метода крупных частиц [6] и модификации метода [10], при решении задач о распространении ударных волн в «замороженной» газовзвеси [8] и в облаке газовзвеси [9].

• 3.    Применение на этапе Эйлера уравнений (2.11)’ и (2.12) позволяет проводить расчеты задач [8, 9] при больших значениях числа Куранта.

Авторы выражают свою благодарность профессору В.Ф. Куропатенко за полезные обсуждения и интерес к работе.

Работа выполнена при поддержке РФФИ грант № 13 - 01 - 00072.