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

Инженерно-технические науки Engineering and technical sciences

УДК 519.688

МОДЕЛИРОВАНИЕ ПОВЕДЕНИЯ ГАЗА С ИСПОЛЬЗОВАНИЕМ РЕШЕТОЧНЫХ МОДЕЛЕЙ

С.П. Бобков, В.Л. Соколов

Ивановский государственный химико-технологический университет

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

Ключевые слова: решеточный газ, клеточный автомат, имитационное моделирование потоков , НРР-модель, столкновения частиц, число Кнудсена, длина свободного пробега.

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

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

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

1. В каждый момент дискретного модельного времени частицы могут находиться только в узлах решетки.

2. Все частицы имеют одинаковую массу, равную единице.

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

4. В один и тот же момент времени в каждом узле решетки не может находиться более одной частицы с одинаковыми векторами скорости.

Вид пространственной решетки, а также правила перемещения и взаимодействия частиц на ней определяются конкретной принятой моделью.

В настоящее время существует несколько разновидностей моделей решеточных газов. В данной статье мы рассмотрим наиболее простую из них - модель НРР, которая названа по первым буквам фамилий ее авторов [1].

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

Рис.1. Схема решетки и векторы скорости модели

Векторы скорости частиц модели НРР можно определить следующим образом:

С 1,3 = (± 1, 0), С2,4 = (0, ± 1). (1)

Развитие процесса поведения газа на каждом шаге по времени происходит в два этапа. Первый этап - перемещение частицы в соседние узлы (этап сдвига). Второй этап - соударение частиц в узлах (этап столкновения).

Для НРР модели принимаются следующие правила столкновений частиц в узлах:

1. Если в какой-либо узел решетки одновременно попадают две частица с разных сторон, то они сталкиваются.

До столкновения

2. При столкновении должно сохраняться количество частиц и их полный импульс.

Нетрудно показать, что для выполнения этих правил столкновения частиц должны происходить «лоб в лоб». При этом после столкновения скорости частиц должны развернуться на 90О. В остальных случаях можно считать, что столкновения не происходит, т.е. частицы пролетают мимо друга.

Направления векторов скорости частиц до и после столкновений показаны на рис. 2.

Рис. 2. Столкновения частиц в модели НРР

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

этапе столкновения частиц. При попадании частицы в узел-стенку происходит либо ее отражение в исходный узел (рисунок 3а), либо отскок частицы в соседний узел пространственной решетки (рисунок 3 б).

Рис.3. Столкновения частицы со стенкой

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

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

Рассмотрим представление НРР модели в терминах теории клеточных автоматов.

Классический клеточный автомат представляет собой множество клеток (ячеек, узлов), заполняющих дискретное О мерное пространство. Каждая клетка является конечным автоматом, входами которого являются выходы граничащих соседних клеток. В соответствии с теорией конечных автоматов они определяются кортежем А=(Х, 8, ф), у которого: Х={х],х2,...,хт} - множество входных сигналов; 8={51,52,...,5/} - множество состояний; ф - функция переходов, которая некоторым парам «состояние - входной сигнал» ставит в соответствие новые состояния автомата. Переход всех клеток (элементарных автоматов) в новое состояние происходит одновременно, переводя клеточную систему в новое глобальное состояние.

В рассматриваемом случае газ представлен клеточным массивом, в котором состояние 8 каждой клетки (узла) представлено булевым вектором длиной 4, поскольку каждая клетка (узел) имеет 4-х соседей: 8 = {51з . . ., 54}. Каждый разряд вектора определяет наличие (5г- = 1) или отсутствие (5г- = 0) в клетке частицы, которая движется со скоростью с в сторону /-го соседа.

Функция переходов ф (X, 8) для автоматов, моделирующих внутренние клети, строится исходя из следующих соображений.

На этапе сдвига происходит изменение одного из разрядов вектора состояний, то есть перемещение частицы. На этапе столкновения происходит изменение нескольких разрядов вектора состояний, чтобы направление движения частицы изменилось в соответствии с правилами функционирования модели. При этом происходит замена одного состояния клетки 8 = {51, . . ., 54} на другое 8 = {5/, . . ., 5/}. Режим работы системы - синхронный.

Для стенок и источников можно также соответствующим образом построить функции переходов. Чтобы эти требования выполнялись нужно, чтобы, во-первых, булевы векторы состояний 8 и 8 имели равное число единиц и, во-вторых, суммы векторов скоростей частиц, соответствующих единицам в 8 и 8' , были равны. То есть:

4

4

Е * = Е * X с, = х

С

7=1

7=1

7=1

7=1

Следует иметь в виду, что в узлах-стенках происходит нарушение закона сохранения импульса из-за взаимодействия частицы со стенкой, а наличие узлов-источников нарушает закон сохранения количества частиц (массы).

Процесс моделирования поведения газа можно представить как эволюцию клеточного поля пространственной решетки в дискретном времени

г = {10, Н, ...}:

П(/0), П(Ю, П(б), . . ., Щеп*) . (3)

При этом клеточное поле определяется векторами состояний клеток и их координатами г :

П = г, 4). (4)

Функция переходов р состоит из композиции функций столкновения р1 и сдвига р2 , то есть:

ф) = р2 ор1(Б). (5)

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

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

кальных конфигурации узла, при которых осуществляется столкновение:

8т = {1, 0, 1, 0} и 8п = {0, 1, 0, 1}.

При столкновении конфигурации изменяются следующим образом:

8т = {1, 0, 1, 0} ^ 8т ' = {0, 1, 0, 1}

или (6)

& = {0, 1, 0, 1} ^ $,' = {1, 0, 1, 0}.

На этом такт заканчивается, и модельное время увеличивается на величину шага.

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

Рис. 4. Моделирование простой диффузии

Рассмотренный пример иллюстри- виях таких примеров не достаточно. Моде-

рует некоторые возможности НРР- ли решеточных газов рассматривают про-

моделей, однако, для доказательства их цессы, протекающие с гипотетическими

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

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

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

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

Адекватность описанного подхода была проверена серией модельных экспериментов. Во всех случаях исследовался клеточный массив размерами 800 х 800.

Массив заполнялся частицами случайным образом. При этом в каждой серии опытов число частиц в объекте (средняя концентрация частиц) было разным. Тем самым моделировалось различное давление модельного газа.

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

80

О1-1-1-1-1-1-1-1-1-1-

О 0.02 0 04 0 06 0.08 0.1 0.12 0.14 0.16

Концентрация частиц, %

Рис. 5. Зависимость рассчитанной величины давления газа от модельной концентрации частиц

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

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

Kn = I / L, (7)

где l - средняя длина свободного пробега частиц в газе, L - характерный размер объекта (например, диаметр трубопровода, длина обтекаемого тела и пр.).

Известно, что если Кп > 1 движение газа можно считать свободномолекуляр-ным. В этом случае можно не рассматривать столкновений частиц (молекул) между собой, а учитывая лишь их соударения с ограничивающими поверхностями. Если Кп < 10-3, то есть длина свободного пробега частиц весьма мала, значительную роль начинают играть столкновения частиц между собой. При этом на длине пробега параметры газа изменяются мало и, благодаря столкновениям частиц, в окрестности каждой точки течения устанавливается локальное состояние, которое уже можно характеризовать такими макроскопическими параметрами, как плотность, температура и т.п. В этих условиях газ можно рассматривать, как сплошную среду и при расчёте течения можно использовать уравнения Эйлера или Навье - Стокса. В области зна-

о

чений числа Кнудсена 10- < Кп < 1 реализуются промежуточные между свободно-молекулярным и континуальным режимы течения разреженного газа.

Имитационный алгоритм поведения газа, использованный при проведении модельных экспериментов, позволял посчитывать не только общее число столкновений частиц (как между собой, так и со стенками объекта), но и число столкновений частиц внутри объема только друг с другом. Расчет показал, что в рассматриваемой серии опытов число Кнудсена для моделируемого газа лежит в диапазоне [0,014 ^ 0,14], то есть имеет место промежуточный режим течения газа. На рисунке 6 представлены данные, показывающие связь между числом Кнудсена для модельного газа и долей взаимного столкновения частиц от общего числа столкновений.

Анализ результатов показывает, что при увеличении числа Кнудсена доля столкновений частиц между собой существенно падает. То есть возрастает число столкновений частиц со стенками моделируемого объекта. Такие результаты полностью согласуются с теоретическими под-

ходами к описанию поведения реальных разреженных газов.

3- 0.9

0.2-1-1-1-1-1-1-1-

О 0.02 0.04 0.06 0.08 0.1 0.12 0.14 0.16

Число Кнудсена Рис. 6. Зависимость доли столкновения частиц между собой от числа Кнудсена

=е 0|-■-■-■-■-■-■-

0 10 20 30 40 60 60

Давление, Па

Рис.7. Зависимость кинематической вязкости газа от давления: сплошная линия - реальные данные; пунктир - рассчитанные по модели

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

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

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

В заключение можно отметить следующее.

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

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

ЛИТЕРАТУРА

1. Hardy, J. Molecular dynamics of a classical lattice gas: transport properties and time correlation func-tions/J. Hardy J., O. de Pazzis, Y. Pomeau // Physical Review.-1976.- V. 13.- N 5.- P. 1949-1961.

2. Бандман, О.Л.. Клеточно-автоматные модели пространственной динамики / О.Л. Бандман // Системная информатика: Сб. науч. тр.- Новосибирск: Изд-во СО РАН, 2006. - Вып. 10. - С. 59-111.

3. Бобков С.П. Анализ возможностей применения решеточных моделей для исследования процессов в газах при пониженном давлении / С.П.Бобков,

B.Л.Соколов // Вестник ИГЭУ. - 2015. - Вып. 4. -

C. 58-63.

4. Frish U. Lattice Gas hydrodynamics in two and three dimensions / U. Frish , J.P Crutchfild, B.Hasslacher // Complex Systems. -Vol.1 (1987) -P. 649-707.

Рукопись поступила в редакцию 16.10.2015.

THE GAS BEHAVIOR SIMULATION BY MEANS OF LATTICE MODELS

S. Bobkov, V. Sokolov

In article possibility of application of The two-dimensional lattice gas model with an orthogonal lattice for the analysis of behaviour of gases is considered. It is shown that the used model is simple and effective at computer realisation. It is possible to recommend it for research of processes in real gases at the lowered piessuie.

Key words: lattice gas, cellular automaton, simulation flows , HPP model, particle collisions, the Knudsen number, mean free path.