ВЕСТНИК ТОМСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА
2008 Управление, вычислительная техника и информатика № 4(5)
УДК 539.2, 536.1
В.П. Маркова КЛЕТОЧНО-АВТОМАТНОЕ МОДЕЛИРОВАНИЕ РАСПРОСТРАНЕНИЯ ВОЛНЫ В НЕОДНОРОДНЫХ СРЕДАХ1
LGСA-модели (Lattice Gas Cellular Automata) основаны на микроскопической модели физического процесса и могут рассматриваться как дополнение к традиционному моделированию пространственной динамики. В работе экспериментально исследована зависимость скорости распространения звуковой волны от параметров простых LGA-клеточных моделей на двумерной решетке с четырьмя соседями. Приведен способ определения параметров клеточно-автоматной модели по физической скорости звуковой волны в некотором веществе. Промоделирован волновой процесс в неоднородной среде, состоящей их двух веществ (алюминий, эбонит). Граница раздела среды не требует дополнительных правил переходов.
Ключевые слова: клеточный автомат, волновой процесс, модельная скорость.
LGCA-модели [1, 2] описывают природные явления множеством частиц, которые движутся в дискретном пространстве. Пространство разбито на множество одинаковых клеток, организованных в регулярную структуру. Каждой клетке ставится в соответствие конечный автомат, состояние которого определяется находящимися в нем частицами. Таблица переходов задается набором правил взаимодействия частиц. Набор правил выбирается таким образом, чтобы при каждом столкновении суммарные масса и импульс всех частиц пространства оставались неизменными. Каждый такт работы автомата разбит на два этапа: столкновение и сдвиг. Все клетки LGCA вычисляют новое состояние синхронно и параллельно, в результате чего происходит изменение глобального состояния LG-автомата. Итеративная смена глобальных состояний LG-автомата описывает динамику моделируемого явления.
В работе исследуется возможность клеточно-автоматного моделирования распространения звуковой волны в неоднородных средах. В качестве моделей используются LG-модели на двумерной решетке с четырьмя соседями: HPP [1] и HPPrp [2]. В отличие от HPP-модели, клетки HPPrp-модели, кроме движущихся частиц, могут содержать частицы покоя (rest particles (rp)), которые позволяют моделировать распространение волны с разной скоростью. В данной работе экспериментально исследована зависимость скорости распространения звуковой волны от параметров HPPrp-модели (количества частиц покоя, вероятностей их присутствия в начальном распределении частиц, вероятности создания (разрушения) частиц покоя и плотности движущихся частиц). Приведен способ определения HPPrp-модели и ее параметров по физической скорости звуковой волны в твердом веществе. Промоделировано распространение звуковой волны в среде алюминий - эбонит и эбонит - алюминий. Граница раздела среды не требует дополнительных правил переходов.
1 Работа поддержана программой Рособразования «Развитие научного потенциала ВШ», проект РНП.2.2.1.1.3653.
1. ЬСЛ-клеточные модели
В НРР-модели состояние каждой клетки (НРР-клетки) определяется вектором скорости движущихся частиц V = (у4 , у3 , у2 , v1 ). Значение (нуль (единица)) 1-го компонента, 1 = 1, 2, 3, 4, вектора скорости означает отсутствие (наличие) движущейся частицы по направлению к 1-му соседу. Скорость движущихся частиц одинакова и равна по модулю 1. На рис. 1, а показана НРР-клетка с вектором у = (0,1,1,0).
2 к
і к 1
3 г
> м
-> V = (0, 1,1, 0)
2
і к 1
3 Г
> г 4
О 8
04 V = (0, 1 1 °)
^ 4 т = (0, 1, 0)
О 2
б
а
Рис. 1. Примеры клеток: HPP-клетка (а) и HPPrp-клетки (б)
Частицы в НРР-клетках сталкиваются по следующему правилу.
Правило (лобовое столкновение). Две частицы, прилетевшие в клетку с противоположными направлениями скоростей, покидают ее, изменяя направления движения на 90° (рис. 2).
±
-к
а б
Рис. 2. Правила столкновения частиц в HPP-клетке
В НРРгр-модели каждая клетка (НРРгр-клетка) может содержать движущиеся частицы и частицы покоя. Состояние НРРгр-клетки определяется вектором масс т и вектором скорости у (рис. 1, б). Частицы покоя имеют одинаковую скорость, равную нулю, и разную массу, равную 2, 4, 8 и т.д. Длина вектора масс т = (тг,тгтх) не зависит от структуры пространства и равна количеству
частиц покоя. Значение компонента тк, к = 1, 2,...,г, в векторе т определяет наличие (отсутствие) частицы покоя с массой 2к в клетке. Далее НРРгр-модель с одной частицей покоя будем обозначать через НРР1гр, модель с двумя частицами покоя будем обозначать через НРР2гр и т.д. Правила столкновения условно можно разбить на три группы.
Группа 1 (лобовое столкновение). Независимо от наличия или отсутствия частиц покоя, движущиеся частицы сталкиваются по правилам НРР-клетки (рис.2).
Группа 2 (создание частицы покоя). Если в клетке для двух (четырех) движущихся частиц справедливы правила столкновения и частица покоя с массой 2 (4) разрушена, тогда в клетке создается частица покоя с массой 2 (4) и исчезают движущиеся частицы (рис. 3, а).
Группа 3 (разрушение частицы покоя). Если в клетке существует частица покоя с массой 2 (4) и отсутствуют 2 (4) движущихся частиц, для которых справедливо одно из правил столкновения, то частица покоя разрушается и появляются 2 (4) движущихся частицы (рис. 3, б).
О 8 О 4 О 2
О 8 О 4 О 2
О 8
О 4 О 2
б
Рис. 3. Правила создания (а) и разрушения (б) частиц покоя
О 8
О 4 О 2
Частицы покоя создаются (разрушаются) с определенной вероятностью Рк, к = 1,2,...,г, причем, при создании (разрушении) частиц покоя должны выполняться следующие условия:
Р к+1^ Р к
и I Р < 1.
к=1
(1)
2. Свойства ИГГ^-моделей
В терминах клеточно-автоматного моделирования процесс распространения звуковой волны в среде представляется клеточным массивом W (рис. 4). Несколько столбцов массива, начиная с 1-го столбца, образуют источник для моделирования внешнего возмущения. Клетки источника - НРР-клетки. Клетки остальной части массива - НРРгр-клетки. Каждой клетке массива Ж ставится в соответствие конечный автомат. Таблица переходов автомата имеет 24+г, г = 0,1,..,4, наборов состояний. Правила переходов образуют правила сдвига и столкновения, для них должно соблюдаться условие (1).
0 1 2
^ 0
г 1
М-1
N-1
Ж
Рис. 4. Клеточный массив для моделирования волнового процесса
а
г
Клетки-источники с некоторой вероятностью в течение одного такта генерируют движущиеся частицы. Исходные состояния НРРгр-клеток генерируются в соответствии с набором вероятностей присутствия частиц покоя в начальном распределении частиц. Ниже набор вероятностей будем обозначать через Ргр ~ (Рг, Рг-1 >—> Р\), где рк , к = 1,2,., г, - вероятность присутствия к-й частицы
покоя. Краевые условия следующие. Левая граница массива (столбец Wo) - стенка с одним правилом столкновения: частица, столкнувшаяся со стенкой, отражается от нее обратно. Правая граница массива (столбец WN-1) открыта, и WN_l = по-
сле каждой итерации. Краевые условия для верхней и нижней границ массива -периодические. Ниже массив из НРРгр-клеток, будем называть HPPгp-средой. Клетку среды с координатой (г, у) будем задавать плотностью движущихся
частиц и плотностью всех частиц d j . Пусть Уу = (у4 (г/), у3 (г/), у2 (/), V (г/)) и ту = (тг (г/), тг-1 (г/),..., т1 (г/)) - вектор скорости и вектор масс клетки с координатой (г,у) соответственно. Тогда dij = XV* (г/), с1у = X1=1 V* (г/) + X*=1(*/). Волновой процесс в НРРгр-среде будем характеризовать осредненной плотностью движущихся частиц (р0), осредненной плотностью среды (р1) и осредненной скоростью распространения звуковой волны (у5). Далее (р1) и (у5) будем называть модельной плотностью и модельной скоростью соответственно. Очевидно, что величина плотности (р1) определяется набором вероятностей присутствия частиц покоя Ргр = {рг,рг_х,...,р1). В месте действия внешнего возмущения
(р) = (р1) + (р2), где (р1) - плотность начального распределения частиц в среде, (р2) - плотность возмущения. Для заданного радиуса осреднения г, осредненные значения модельной плотности для у-го, у = 0, 1,., N-(^+1), столбца массива Ж
ЕГ ^М-1 ^
/=-г 2 1=0 У т> /г 1
вычисляются следующей суммой ф у; =---------------------—. В табл. 1 приведены зна-
' (2 г +1 )М
чения плотности (р1) для всех сред при условии, что вероятности присутствия всех частиц покоя в начальном распределении равны 0,5. Очевидно, что для Нрр-среды (р1) = (р0).
Таблица 1
Модель Hpp Hpp1rp Hpp2rp Hpp3rp Hpp4rp
<р‘> 2 3 5 9 17
В наших экспериментах процесс распространения звуковой волны моделируется в массиве размером 200 х 800 клеток. 20 столбцов массива образуют источник внешнего возмущения, каждая клетка источника с вероятностью 0,95 генерирует 3 движущихся частицы. Определение модельной скорости звука в среде продемонстрируем на примере. На рис. 5 показано изменение функции плотности частиц в НРР2гр-среде при воздействии возмущения через каждые 100 итераций. В данной среде частицы покоя с массой 2 и 4 в начальном распределении присутствуют с одинаковой вероятностью, а именно, p = p2 = 0,5, величина осредненной плотности движущихся частиц (р°) = 2. Количество клеток массива, на которое передвинулось максимальное значение плотности частиц в среде за единицу времени, будем считать модельной скоростью звука в среде с данными характеристиками.
На рис. 6 приведены зависимости скорости звуковой волны от плотности движущихся частиц, 0 < (р°) < 4, для всех НРРгр-сред. Кривые скорости, независимо от вероятности присутствия частиц покоя, симметричны относительно (р°) = 2. Во всех средах звуковые волны достигают минимальной скорости и максимально отличаются друг от друга по величине для (р°) = 2. Модельная скорость звука может быть увеличена (уменьшена) следующими способами.
Длина массива Рис. 5. Изменение функции плотности в Hpp2rp-среде
0,8
н
о
Л
§
о
0,6 -
§ 0,4
ч
а>
I
0,2 -
0 0,2
0,6
3,4
1 1,4 1,8 2,2 2,6 3
Плотность движущихся частиц
Рис. 6. Зависимость модельной скорости от плотности движущихся частиц
3,8
3. Изменение начального распределения частиц покоя
При фиксированной плотности движущихся частиц скорость звуковой волны может быть увеличена (уменьшена) за счет уменьшения (увеличения) вероятности присутствия частиц покоя в начальном распределении. Для заданной НРРгр-сре-ды, звук распространяется с наибольшей модельной скоростью в той НРРгр-среде, у которой частицы покоя с максимальной массой имеют минимальную вероятность в начальном распределении частиц (рис.7).
8
0 а § О
3
я
4
1
Плотность движущихся частиц
Рис. 7. Зависимость скорости звуковой волны в HPP2rp-среде от плотности при различных вероятностях присутствия частиц покоя с массой 4
4. Изменение правил столкновения
Модельная скорость волны может быть увеличена (уменьшена) за счет уменьшения (увеличения) вероятности создания частицы покоя из двух (четырех) движущихся частиц. Очевидно, что наибольшую модельную скорость звуковая волна достигает в НРРгр среде с таким набором правил столкновения, в которых вероятность создания частиц покоя уменьшается с увеличением их массы. Так при уменьшении вероятности создания частицы покоя с массой 8 в НРРЗгр среде с 0,5 до 0,25, скорость звуковой волны в эксперименте достигла (рис. 8) величины модельной скорости ум звуковой волны, полученной по формуле (2), приведенной в работе [2],
га (а -1)
2
V = -
т
(2)
4^-1) + 2 ткркак(1 - с1к) к=1
где d - средняя плотность движущихся частиц в клетке, т.е., d = (р0), тк - масса к-й частицы покоя, рк - вероятность наличия к-й частицы покоя в начальном распределении частиц, dk = dтк /тк + (1 - d)тк - плотность ^й частицы покоя в клетке среды НРРгр.
Плотность движущихся частиц
Рис. 8. Зависимость скорости звуковой волны в НРРЗф среде от плотности при различных вероятностях создания частицы покоя с массой 8
Г
5. Соответствие физических и модельных скоростей звуковой волны
Для клеточно-автоматного моделирования звуковой волны в неоднородных средах, состоящих из нескольких веществ, каждому веществу необходимо поставить в соответствие НРРгр-среду по физической скорости Уф распространения звуковой волны в заданном веществе и наоборот. Для этого нужно определить масштабный коэффициент пересчета скоростей (физической скорости в модельную скорость и наоборот). Пусть звуковая волна распространяется в твердых веществах со скоростью 2000 - 6320 м/с. Максимальную физическую скорость звуковая волна (м/с) достигает в алюминии, максимальную модельную скорость (клетка/итерация) звук достигает в НРР-среде. Следовательно, алюминию можно поставить в соответствие НРР-модель. В общем случае в качестве базовой модели можно взять НРРгр-модель при условии, что отношение физических скоростей звука в веществах не превосходит отношения скоростей в моделях, соответствующих этим веществам. Ниже в качестве базовой модели будем использовать НРР-модель, тогда
масштабный коэффициент определяется следующим соотношением:
Ц(ал) = — = 8,94-103. v
’ Л/Г
м
Например, звуковая волна в никеле распространяется со скоростью 5400 м/с; согласно масштабному коэффициенту, модельная скорость звука составляет величину, равную 0,634. Такую скорость звуковая волна может достигать во всех НРРгр-моделях с частицами покоя, но при разной плотности движущихся частиц и заданной вероятности присутствия частиц покоя в начальном распределении (рис.6). Выбор модели зависит от задачи. Чаще всего выбирается та модель, которая при полученной скорости имеет плотность, близкую к равновесному состоянию среды. Если такую модель выбрать нельзя, то следует изменить вероятность присутствия частиц покоя и (или) вероятность создания (разрушения) частиц покоя. Так, для нашего примера, никелю соответствует НРР1гр-среда с вероятностью 0,25 присутствия частицы покоя с массой 2 .
В экспериментах в качестве неоднородных сред рассматриваются два твердых вещества: алюминий (легкая среда) и эбонит (тяжелая среда). Эбониту (уф = 2500 м/с), согласно коэффициенту пересчета р(Ш), соответствует среда с модельной скоростью, равной 0,3. Модельную скорость ум = 0,3 звук достигает в НРР2гр-среде, которая характеризуется следующими параметрами: плотность (р°) = 2, вероятности присутствия частиц покоя равны р = 0,5, р2 = 0,44, правила столкновения равновероятны. Алюминий, как было определено выше, соответствует НРР-среде. Процесс распространения звуковой волны в неоднородной среде представляется клеточным массивом (размер массива 800 х 200 клеток), состоящим из источника и двух массивов, один массив образован НРР-клетками, второй массив - HPP2rp-клетками. Граница раздела среды не требует никаких дополнительных граничных условий и правил столкновения.
Пример 1. (Моделирование волнового процесса в среде алюминий - эбонит.) В эксперименте граница раздела среды проходит по 250-му столбцу клеточного массива. На 300-й итерации звуковая волна наталкивается на границу раздела сред (рис. 9) и происходит преломление и отражение волны.
6. Результаты КА-моделирования процесса распространения звуковой волны в неоднородных средах
8
0
100
200
-------500
-------600
0 50 100 150 200 250 300 350 400
Длина массива
Рис. 9. Распространение волны в среде алюминий - эбонит
Эксперимент показал, что звуковая волна распространяется в тяжелой среде (эбоните) со скоростью, равной модельной скорости звуковой волны в среде, соответствующей эбониту (рис. 10). Скорость отраженной волны равна скорости падающей волны.
£ 0,8 О
0 Л §0,6 3
1 0,4
<и £
Е 0,2
74 142 210 262 291 320 350 380 409
Длина массива
Рис. 10. Изменение модельной скорости волны в среде алюминий - эбонит
Пример 2. (Моделирование волнового процесса в среде эбонит - алюминий.) В эксперименте граница раздела сред проходит по 180-му столбцу. На 500-й итерации волна доходит до границы раздела эбонита и алюминия (рис. 11). Далее волна плавно переходит в легкую среду. Моделирование показало, что скорость звуковой волны после границы раздела сред, т.е. в алюминии (рис. 12), практически совпадает с модельной скоростью звуковой волны в НРР-среде, которая в нашем случае соответствует алюминию.
Длина массива
Рис. 11. Распространение волны в среде эбонит - алюминий
Длина массива
Рис. 12. Изменение модельной скорости волны в среде эбонит - алюминий
Заключение
В работе экспериментально исследована зависимость скорости распространения звуковой волны от параметров HPPrp-моделей (плотности движущихся частиц, количества частиц покоя, вероятностей присутствия частиц покоя в начальном распределении частиц и вероятностей создания частиц покоя). Результаты подтверждаются расчетными данными из работы [2]. Приведен способ определения HPPrp-модели и ее параметров по физической скорости звуковой волны в некотором веществе. Промоделирован волновой процесс в неоднородной среде, представляющей два разных вещества (эбонит, алюминий), где переходный слой не требует дополнительных правил переходов. Моделирование показало, что наличие границы раздела в неоднородной среде не изменило скорость звуковой волны в среде за границей раздела (эбоните и алюминии) по сравнению с полученной модельной скоростью звука в однородной среде, соответствующей эбониту и алюминию.
ЛИТЕРАТУРА
1. Hardy J., Pomean Y., de Pazzis O. Time evolution of a two-dimensional model system // J. Math. Phisics. 1973. V. 14. P. 1746.
2. Zhang M., Cule D., Shafai L., Bridges G. and Simons N. Computing Electromagnetic Fields in Inhomogeneous Media Using Lattice Gas Automata // Proc. of 1998 Symposium on Antenna Technology and Applied Electromagnetics, Aug.14 - 16. Ottawa, Canada, 1998. P. 1.
Статья представлена кафедрой программирования факультета прикладной математики и кибернетики Томского государственного университета и оргкомитетом 7-й Российской конференции с международным участием «Новые информационные технологии в исследовании сложных структур», поступила в научную редакцию 3 октября 2008 г.