Диссертация (1149881), страница 11

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 11)

. . 40◦ . В отличие от полуэмпирической модели, приведённой выше,здесь контакт частицы с поверхностью считается точечным. Также считается, что частица не изменяет своего положения за время удара.Для рассмотрения несферических частиц необходимо ввести связанную с частицей систему координат. Она задаётся координатами её начала(xp , yp , zp ) и тремя углами Эйлера (ϕ, ψ, ϑ) (рис.

3.4).Рис. 3.4: Переход к связанной системе координат.Для перехода из базовой системы координат в связанную вводится ортогональная матрица преобразования A, или матрица направляющих косинусов.cosψcosϑsinϑ−sinψcosϑ sinϕsinψ − cosϕcosψsinϑ cosϕcosϑ sinϕcosψ + cosϕsinψsinϑA= cosϕsinψ + sinϕcosψsinϑ −sinϕcosϑ cosϕcosψ − sinϕsinψsinϑ.(3.2)Рассмотрим произвольный вектор b.

Переход из одной системы коор-72динат в другую осуществляется с помощью соотношений: bbbbξ ξ  x  x bη  = A  by  ,  by  = AT  bη  . bζbzbzbζОпишем схему вычислений параметров частиц после отскока от поверхности. Пусть Vc− и Vc+ − скорости точки контакта на поверхности частицыдо и после удара, rc − радиус-вектор из центра масс частицы в точку контакта с поверхностью, тогда:∆Vc = Vc+ − Vc− = ∆Vp + ∆! p × rc ,(3.3)1||Jp ||∆! p = rc × ∆Vc − rc × [∆! p × rc ].(3.4)mpКомпоненты вектора скорости точки контакта ∆Vc вычисляются следующим образом: ∆Vcξ∆uc∆uc  ∆Vcη  = A  ∆vc  ,  ∆vc ∆Vcζ∆wc∆wcVcτ−Vcn−Vcz− = u−c−vcwc−  , u−c−vcwc− = u−p−vpwp− = Vcτ− (aτ c − 1)Vcn− (−anc − 1)Vcz− (aτ c − 1),ωr  T  pξ1   cx  + A  ωpη1  ×  rcy  ,  ωpζ1rczгде anc и aτ c = azc − коэффициенты восстановления компонент скороститочки контакта Vcτ , Vcn , Vcz (см.

рис. 3.3).Если рассматривать уравнение (3.4) в связанной системе координат, тотензор инерции ||Jp || будет содержать только диагональные компоненты− главные моменты инерции Jpξ , Jpη , Jpζ . Тогда оно приводится к системелинейных алгебраических уравнений относительно ∆ωpξ , ∆ωpη , ∆ωpζ :J pξ + ηc2 + ζc2−ξc ηc−ηc ξcJ pη + ζc2 + ξc2−ζc ξc−ζc ηc−ξc ζc ηc ∆Vcζ − ζc ∆Vcη ∆ωpξ −ηc ζc  ∆ωpη  =  ζc ∆Vcξ − ξc ∆Vcζ  .∆ωpζξc ∆Vcη − ηc ∆VcξJ pζ + ξc2 + ηc273Конечные компоненты скорости поступательного и вращательного движения частицы определяются следующим образом:∆ωpx∆ωpξ ∆ωpy  = AT  ∆ωpη  ,∆ωpz∆ωpζ−u+p = up + ∆uc − ∆ωpy rcz + ∆ωpz rcy ,vp+ = vp− + ∆vc − ∆ωpz rcx + ∆ωpx rcz ,wp+ = wp− + ∆wc − ∆ωpx rcy + ∆ωpy rcx .3.4Моделирование разброса частиц примеси по размерамПри расчёте полидисперсной примеси считалось, что размеры частиц внабегающем потоке распределены по логарифмически-нормальному закону.

В этом случае функция массового распределения частиц по размерамимеет вид:" 2 #1lnr−lnrgo√pg∞exp −.(rp ) = √2πrp ln s2 ln s(3.5)Здесь rg = rpg exp(ln2 s), ln s − среднеквадратичное отклонение логарифма радиуса частиц, rpg − наиболее вероятный размер частиц. На рис.3.5 представлен график функции массового распределения, использовавшейся в расчётах.Алгоритмическая реализация задания размеров частиц выполняетсяследующим образом.Этап 1. Получение нормально распределённой величины с помощьюметода полярных координат Бокса−Мюллера−Марсальи [23]:шаг 1.

Получение случайных равномерно распределённых в интервале(0,1) величин u1 и u2 .шаг 2. Вычисляются величины v1 = 2u1 − 1, v2 = 2u2 − 1.74Рис. 3.5: Массовая функция распределения частиц по размерам в невозмущённом потоке.шаг 3. Вычисляется величина S = v12 + v22 .шаг 4. Вычисляется нормально распределённая величинаpRN = v1 −2 ln(S)/S.Этап 2. Вычисляется случайный радиус произвольной частицы (rp∗ ) изоблака с помощью соотношения:ln(r∗p ) = RN · ln(S) + ln(rg ).При расчёте радиуса частиц накладывалось ограничение, отсекающеечасть кривой на рис.

3.5 ниже уровня g∞ (rp ) = 10−4 . Если вычисленныйслучайный радиус частицы меньше чем 10−4 , то происходит его повторноевычисление.3.53.5.1Численный методЧисленный метод расчёта движения одной частицыДля интегрирования уравнений движения частиц используется методтипа предиктор−корректор второго порядка точности по времени:vp∗ = vpn +∆tp n(fD + fMn ),mp!∗p = !np +75∆tp nL ,Jp pr∗p = rnp + ∆t vpn ,vpn+1 = vpn +∆tp n(fD + fMn + fD∗ + fM∗ ),2mp!n+1= ! np +p∆tp n(Lp + L∗p ),2Jp∆t n(vp + vp∗ ).2Это метод используется как в случае расчёта бесстолкновительной приrn+1= rnp +pмеси, так и в случае учёта столкновений между частицами.3.5.2Метод Монте-Карло (DSMC) для расчёта столкновительной примесиВ общем случае количество моделирующих частиц может отличатьсяот количества частиц в реальном течении.

Тогда размер моделирующихpчастиц рассчитывается по правилу: rm = rr · Nr /Nm , где rr и Nr − радиуси счётная концентрация моделируемых частиц.Использованная "техника" прямого статистического моделированиястолкновений между частицами основана на методе мажорантной частоты[18], усовершенствованном в работе [137]. Рассмотрим ячейку, в которойнаходится N частиц.

Если в ячейке находится только две частицы − i и j,то частота их столкновения друг с другом равна ϑij = π(ri + rj )2 |vj − vi |.Введём мажоранту ϑmax таким образом, что ϑmax = const > ϑij для любыхпар частиц i и j. В каждой ячейке осуществляется следующий вычислительный цикл:(1) Рассчитываются все возможные комбинации ϑij , наибольшее значение из них и будет ϑmax . Далее рассчитывается мажоранта столкновениядля всех N частиц в ячейке.νmax =N (N − 1) ϑmax,2Ωгде Ω-объём ячейки.(2) С помощью генератора случайных чисел получается случайная величина β в интервале [0,1] и вычисляется шаг по времени τk = −lnβ/νmax ,являющийся промежутком времени до следующего столкновения.76(3) Проверяется условие прерывания цикла, которое приведено ниже,если оно выполняется, то переходим к пункту (5), иначе выполняем пункт(4).(4) Разыгрывается возможное столкновение.

Для этого равновероятно выбираем две частицы i и j из всего множества частиц в ячейкеи вычисляем величину ϑij = π(ri + rj )2 |vj − vi |. Проверяем условиеβ < ϑij /ϑmax . Если оно выполняется, то находим компоненты вектораpnij = cosχi + sinχcosεj + sinχsinεk, где cosχ = 2β − 1, sinχ = 1 − cos2 χ,ε = 2πβ. Если выполняется условие gij · nij <= 0, то рассчитываем параметры частиц после столкновения друг с другом по формулам (3.6). Затемпереходим к пункту (1) для моделирования следующего возможного столкновения в ячейке.(5) Увеличиваем счётчик времени на сумму интервалов τ n+1 = τ n +Pτkи переходим к расчёту в следующей ячейке.(6) После расчёта столкновений на заданном интервале времени во всехячейках, производится расчёт поступательного и вращательного движениячастиц.Условие прерывания циклаПусть движение частиц в ячейке было рассчитано до момента времениtn .

Теперь надо рассчитать их движение на интервале [tn , tn+1 ]. На этоминтервале происходит К столкновений между частицами с различнымиинтервалами времени τk между ними. Прерывание цикла выполняется привыполнении одного из двух условий:nt +KXn+1τk < tn≤t +k=1nt +K−1XK+1Xτk ,k=1n+1τk < tk=1n≤t +KXτk .k=1В первом случае мы прерываем цикл, если следующее столкновениепроизойдёт после момента времени tn+1 . Во втором случае рассчитываемэто столкновение и затем выходим из цикла. Либо первый, либо второй77случай выбираются равновероятно. Такой подход позволяет улучшить сходимость метода.

Впервые он был предложен в работе [134]. После вычисления поля течения несущего газа, будет происходить расчёт столкновенийчастиц от момента времени последнего столкновения до текущего моментавремени по газу. Иными словами на каждом временном шаге расчёта несущей среды последнее рассчитанное столкновение может как отставать, таки опережать по моменту времени рассчитанный временной шаг для газа.3.5.3Модель соударения двух частицЭта модель была предложена для сталкивающихся частиц различногоразмера в предположении, что они не деформируются при парных столкновениях, однако соударение считается неполностью упругим. В даннойработе она используется для частиц одинакового размера. Параметры частиц до и после столкновения обозначены верхними индексами "−" и "+"соответственно.

Нижний индекс означает номер частицы − "1" или "2". Запишем уравнения сохранения импульса и момента импульса для системыиз двух частиц [9]:−−++mp vpi+ mp vpj= mp vpi+ mp vpj,−+−Ip (! +pk − ! pk ) = mp rp ek × (vpk − vpk ), k = {i, j}, ei = nij , ej = −nij .Данная система является незамкнутой, поскольку содержит четыре++неизвестных вектора: vpi, vpj,!+pi, !+pj .Введём относительную скоростьповерхностей частиц в точке контакта:U = g + rp n × (! 1 + ! 2 ).Тогда относительную скорость поверхностей частиц в точке контактапосле столкновения можно вычислить следующим образом:−−U+ij = −apn Uij(n) + apt Uij(τ ) ,Uij(n) = (Uij · nij ), Uij(t) = (Uij − Uij(n) ).78В расчётах, коэффициенты восстановления нормальной и касательнойсоставляющих относительной скорости контакта частиц U были приняты равными apn = 0.5 и apt = 0.9.

Они характеризуют потери относительной скорости вследствие неупругого характера столкновения и шероховатости поверхности частиц. Их значения могут находиться в пределах0 ≤ {apt , apn } ≤ 1.Запишем конечные соотношения для вычисления скоростей частиц после соударения:v1+ = v1− + J/mp , v2+ = v2− − J/mp ,!+k = !−k +mpJ=2(3.6)rp×J, k = {1, 2},Ip2(1 + apn )(u− · n)n + (1 − apt )(u− − (u− · n)n) ,7где J − импульс сил, действующих в точке контакта поверхностей частиц.3.6Выводы по третьей главеОписаны модели течения и алгоритмы расчёта параметров дисперснойфазы, основанные на лагранжевом описании бесстолкновительной примесии на кинетическом подходе при учёте столкновений между частицами.

Характеристики

Список файлов диссертации