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

_______УЧЕНЫЕ ЗАПИСКИ Ц А Г И

Т о м IV 197 3

№ 4

УДК 533.7

СТАТИСТИЧЕСКИЙ МЕТОД РЕШЕНИЯ ПРИБЛИЖЕННОГО КИНЕТИЧЕСКОГО УРАВНЕНИЯ

Ю. И. Хлопков

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

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

1. Описание метода. Рассмотрим приближенное, так называемое модельное кинетическое уравнение [1]

§=■^1 <1.1> ЯГ Т

здесь /о — локально равновесная функция распределения,

3/2 ,

т (? — и)2

/о=я<ійИ ехр

п, Т, и — макропараметры — плотность, температура и скорость, зависящие ог координат физического пространства; •

1

х — — среднее время между столкновениями, зависящее от сорта молекул;.

. а0 — среднее сечение столкновения молекул;

£■—относительная скорость.

Уравнение (1.1) является также нелинейным интегро-дифференциальным уравнением, так как в функцию /0 функция распределения / входит нелинейным образом. Несмотря на то что уравнение (1.1) не имеет строгого математического обоснования, оно дает качественно правильные результаты для широкого круга задач при произвольной длине пробега молекул и возмущениях любой интенсивности. При решении этого уравнения, так же как и при решении полного уравнения Больцмана, необходимо интегрирование по скоростям молекул, что при численном решении на ЭЦВМ приводит к необходимости запоминать функцию распределения. Чтобы избежать это, другими словами, сократить необходимый объем оперативной памяти ЭЦВМ, воспользуемся методом Монте-Карло в приложении к уравнению (1.1).

Решение нелинейных задач методом Монте-Карло обычно сводится к решению последовательности линейных задач. В нашем случае уравнение (1.1), записанное в интегральном виде, представляет собой нелинейное уравнение

/(П=/ К(У, Уи /т)/(Уг)аУ1+/т(У), (1.2)

где У — фазовое пространство (х, е); /да — функция распределения на границе. Это уравнение решается методом последовательных приближений по схеме

/к(У) = /К(У, у,, Гк-1(У^))/к(У1)аУ1+и(У), к = 1,2,... (1.3)

С уравнением (1.3) связывают однородную цепь Маркова, заданную плотностью начального распределения /т(У) и переходной плотностью К(У, У^-Задача считается решенной, если вычислен функционал

У = //(У)Ф(У)йУ, (1.4)

соответствующий макропараметрам в поле течения, где Ф1 = 1, 5, ... ный признак.

Известно, что математическое ожидание случайной величины С,

-молекуляр-

: ^ ^ равно 1=1

функционалу (1.4). Здесь I—число рассеяний частицы вдоль траектории с переходной вероятностью, равной ядру уравнения (1.3). Принципиальное отличие модельного уравнения (1.1) от полного уравнения Больцмана ^заключается в предположении, используемом при составлении (1.1), о том, что функция распределения молекул после столкновения является наиболее вероятной при заданных числах сталкивающихся частиц, их импульсе и энергии. Уже в результате одного столкновения система с произвольной функцией распределения переходит в состояние, близкое к равновесному. Таким образом, мы избавляемся от необходимости в каждом приближении запоминать распределение молекул по скоростям, потому что после столкновения плотность вероятности скорости пробной частицы В точке столкновения есть /о-

2. Течение Куэтта. Метод был апробирован на одной из простейших задач

о течении и теплопередаче между двумя параллельными бесконечными пластинами, движущимися друг относительно друга (фиг. 1). На этом сравнительно простом течении опробованы почти все известные методы решения уравнения Больцмана.

Рассмотрим течение Куэтта при малых относительных скоростях пластинок т± и малых отношениях температур пластинок Если ш>± и Тт мало от-

личаются друг от друга, то функция распределения / мало отличается от некоторого максвелловского распределения /од:

Фиг. 1

здесь /„о = и0

т

3/2

(2.1)

, 2тс£7'0

і (х, і) — малая добавка.

ехр —

ті2 \

' '2кТ0] ;

Ю9

Плотность частиц и температура могут быть представлены в виде п(х) = пй[\ +->(*)]; Г(х) = Го!1 + 0И],

где v {х) = _L Г/00 cprfS; 6 (*) = J_ 1 Г 2j7oo<P<ft - v (*).

«о J 3 Ц /0 J 2

Скорость потока по определению равна '

«у = — С /оо £у n0 J

Можно показать, что в этом случае задача может быть разбита на две: на задачу о чистом сдвиге при Tw+ = Tw_, для которой ч = 0 = 0, и задачу о теплопередаче при = 0, для которой иу — 0. В обоих случаях был использован способ введения статистических весов [2].

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

= • /2 2)

Г0 Р0 РоП"” ,

здесь и ^ — вероятности вылета пробной молекулы с границы с абсолютно равновесной и с возмущенной функциями распределения соответственно; р0, р— вероятности пролета расстояния х и столкновения в элементе Лх; <50, (?— вероятности обладания скоростью ^ пробной частицы после столкновения.

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

Макроскопические величины Ф1, находятся в виде

1 " {У

<“>

где УУ—число траекторий; а — номер траектории; р (а) —число пересечений траекторией а границы данной физической ячейки.

Число Кнудсена вводится как отношение длины пробега X к расстоянию Между пластинами: Кп = \/й.

На фиг. 2, а представлены профили скорости, отнесенной к скорости стенки между пластинами в зависимости от числа Кнудсена, на фиг. 2, б —- профили температуры, отнесенной к температуре стенки. Точками на фиг. 2,а обозначены результаты, полученные в данной работе; для сравнения сплошными кривыми представлены результаты точного численного решения модельного уравнения, взятые из работы [1].

3. Решение задач внешнего обтекания. Рассмотрим обтекание сферы диаметром £> сверхзвуковым потоком газа. На бесконечности* практически на расстоянии нескольких калибров (границы течения определялись экспериментально), задается невозмущенная функция распределения с параметрами на бесконечности

(

m

/о-«оо \2-KkTU ехр

3/2

гп (с - UV

2kT.

(3.1)

ПО

Поскольку мы имеем дело с осесимметричным течением, то для удобства выбирается цилиндрическая система координат с началом в центре сферы. Ось х направлена вниз по потоку, а угол 9 отсчитывается от оси у в плоскости, перпендикулярной оси х (фиг. 3). Поле течения разбивается на ячейки по осям х и у. Таким образом, расчет течения происходит внутри цилиндра радиусом у+ и высотой х~ -\-х+. На поверхности цилиндра задается распределение (3.1), а на поверхности сферы условие непротекания, т. е. поток Задающих на поверх* ность сферы молекул I = (I — нормаль к поверхности сферы) равняется

потбку молекул, отраженных от поверхности /та = У /та(£/)^:

У>0

/=/„,. (3.2)

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

■ (т. \з/2 / с2 \

Л»-Я«(ак*7’в1) ехр(~2И^); (3'3)

здесь находится из условия непротекания (3.2).

Поле течения в силу осевой симметрии можно разделить на равные части с углом раствора <р, которые не будут отличаться друг от друга, а на плоскостях сечения задать закон зеркального отражения молекул. В данном случае рассматривался сектор с углом раствора Дер, равным одной ячейке <р (см. фиг. 3).

Этот случай удобен тем, что в сечении, перпендикулярном набегающему потоку, сектор с углом раствора Д<р можно заменить треугольником и тогда поле течения будет представлять собой клин высотой у+, длиной х~ -\-х+ и углом раствора Д'-?^:Ьх1у+. Сечение клина разбивается на ячейки по х и у с ребрами ячеек Ах и Лу. Влет молекул происходит с тор- у ца И основания клина, а от его образующих молекула отражается зеркально.

Траектория заканчивается, когда молекула вылетает из клина. Число 5 вводится как отношение скорости потока и к

тепловой скорости молекул: 5 => и |/ ■ а число Кп = где Х^— дли-

на свободного пробега молекулы в невозмущенном течении.

Были заданы следующие параметры: £^=1,0; Кп = 1 и температура поверхности сферы, равная температуре на бесконечности, т. е. Тт = . Сетка

в поле течения бралась постоянной с размерами ячеек Ддг = Ду = Л/6. Численным экспериментом было установлено, что поле течения можно ограничить следующими параметрами: х~ — расстоянием вверх по потоку, равным длине трех калибров; х~^ — расстоянием вниз по потоку, равным длине четырех калибров: у+ — расстоянием в направлении, перпендикулярном потоку, равным длине трех калибров. В первом приближении в поле течения задается, как и в работе (3], свободномолекулярное распределение параметров или постоянное распределение с параметрами, равными параметрам на бесконечности. Блуждание частиц происходит на заданном или полученном из предыдущего приближения поле течения до тех пор, пока распределение параметров я, 5 и Т в поле течения не установится. Практически для установления течения достаточно 30 тыс. траекторий. Затем переходят к следующему приближению. Макропараметры в каждой ячейке вычисляются по формуле (1.4). На фиг. 4 показано распределение плотности в поле течения в окрестности холодной сферы. Сплошными линиями нанесены линии равной плотности, отнесенной к плотности набегающего потока, полученные в данной работе. Для сравнения пунктиром нанесены изолинии, полученные в работе [3]. Коэффициент сопротивления сферы сх оказался равным 3,5. Следует отметить, что, несмотря на то что в работе [3] молекулы рассматривались как твердые сферы, а в данном случае как псевдомаксвелловские, результаты получились довольно близкими друг к другу.

Расчет проводился на ЭЦВМ БЭСМ-6. Для установления решения необходимы три-четыре итерации. Время одной итерации примерно 45 мин.

При тех же параметрах проводился расчет обтекания цилиндра. Поле плотности, отнесенной к параметрам на бесконечности, представлены соответственно на фиг. 5. В случае обтекания цилиндра техническая схема расчета существенно упрощается, поскольку поле течения из осесимметричного становится плоским и расчет осуществляется в прямоугольнике длиной х~ + х+ и высотой у+. Влет молекул происходит с границ прямоугольника х = х~, х — х+ и у = , а на линии у = 0 происходит зеркальное отражение

молекул. Верхняя половина цилиндра помещается в начале координат. Количество итерации остается таким же, как и для сферы, а время на одну итерацию примерно 35 мин. Коэффициент сопротивления цилиндра ^ = 3,8.

Статистическая ошибка метода равняется примерно 3%. Проверка была произведена по известному точному решению—из области течения удалялось тело и вычислялось равновесное распределение параметров.

Расчет был произведен также для цилиндра с адиабатической поверхностью при 5^ = 5. На фиг. 6, а представлены линии равной плотности я/ях в поле течения над цилиндром при Кп=1. Для сравнения пунктиром нанесены изолинии.

Фиг. 4 Фиг. 5

Фиг. 6

взятые из работы [4]. Сравнительно небольшую разницу в параметрах течения’ по-видимому, можно объяснить различием исходных уравнений. (В работе [4] методом Монте-Карло решается полное уравнение Больцмана).

На фиг. 6, б представлены линии равной плотности при 5^ = 5 и Кп=0,33. Коэффициент сопротивления оказался равным 1,9.

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

Автор выражает признательность М. Н. Когану за руководство и большую помощь, оказанную в процессе работы.

ЛИТЕРАТУРА

1. Коган М. Н. Динамика разреженного газа. М. „Наука”, 1967.

2. Власов В. И., Горелов С. Л., Коган М. Н. Математический эксперимент для вычисления коэффициентов переноса. Докл. АН СССР, 1968, т. 176, № 6.

3. Л а р и н а И. Н. Обтекание сферы разреженным газом. ПММ, т. 33, вып. 5, 1969.

4 Vogenitz F. W., Bird Q. А., В г о a d w е 11 I. Е., R u п g а 1-dier Н. Theoretical and experemental study of rarfied supersonic flows about several simple shapes. AIAA Juornal, No 12, 1968.

Рукопись поступила 2/Х 1972