Полак__Применение_вычислительной_математики_в_химической_и_физической_кинетике (972294), страница 40
Текст из файла (страница 40)
При условии замкнутости система может переходить из одного состояния в другое Все эти работы выполнены совместно лабораторией радиационной химии ИНХС АН СССР и лабораторией вычислительных методов ВНИИЯГГ на основе алгоритмов, рааработанных во ВНИИЯГГ. 183 (б) К определенной выше системе М частиц теперь можно применить формальный аппарат метода Монте-Карло. Пусть а)»,.'— полное сечение упругого рассеяния частиц ~ и ) в системе центра масс. Оно в общем случае зависит от ~п« вЂ” п,~ и определяется конкретным видом потенциала взаимодействия между частицами.
В приближении твердых сфер о(«у и 4 х где д« вЂ” диаметр»-й частицы. В лабораторной системе координат сечение упругого рассеяния ~-й частицы на у-й будет иметь вид (у) и 1э Теперь можно опре ~елнть полное сечение взаимодействия (-й частицы со всеми остальными и '=,"ч 'Е' (8) 3=1 Очевидно, что величина а; соответствует некоторому феноменологическому сечению взаимодействия ю-й частицы со средой. Однако только посредством упругих столкновений частиц. Поскольку мы ие рассматриваем конфигурационного пространства, временнбе поведение системы не является детерминированным, последовательность переходов системы из одного состояния в другое является случайным процессом, а сами эти состояния образуют марковскую цепь. Вероятности переходов между различными состояниями не зависят от времени и полностью определяются набором скоростей всех частиц.
Чтобы получить возможность описания макроскопических систем, нужно было бы положить )«' равным примерно числу Авогадро. Ввиду ограниченных возможностей современных ЭВЫ воспользуемся несколько модифицированным методом периодических граничных условий. При описании системы набором скоростей всех частиц он сводится к разбиению бесконечной системы частиц на Л' групп таким образом, что скорости всех частиц в каждой группе близки по величине и направлению друг к другу. В каждой группе выделяется «типичная» частица и считается, что остальные частицы в группе ведут себя аналогично атой частице.
Таким образом, если п — физическая концентрация частиц, величина п~Х будет соответствовать концентрации каждой из Х «типнчных» частиц. Отметим, что частицы системы могут быть разного сорта — а, р, ..., у, но при этом должно обязательно выполняться соотношение зм) з(Ю у ео у<Р) Р(«,пы(!) =. 1 — Р(Е~щ) !) = 1 — (1Р(!«) !). (10) 1=1 Учитывая выражение (9), получаем для плотности распределения минимальных времен Х Ф дР (1 .,„(1) р(! 1„) = ~'," =- 'Я с«!«1;/ехр ~ — 1,~'„о«!«11!). (11) 1=1 1=1 В этом выражении все взаимодействия частиц учтены дважды, поэтому выражение (8) надо заменить, например, на следующее —,Я~ ои. ,=1+1 Таким образом, согласно формуле (11), время можно разыгрывать по формуле !а (! — Ц! (12) г« где $ — случайное число, равномерно распределенное в интервале (01), а ь' М М э си Р = ~~'~ с« ( «11 ( = —. ~~'~ ~~'~ с«1 ( «1, — «1, !.
(13) Л (8а) ~ С»1. стр. !ао . в такой модели система не подразделяется на «среду» и «апсамбль пробных частиц». Каждую из Л' частиц можно рассматривать как пробную, но их совокупность является в то же время и средой. Выше уже говорилось, что система Х-частиц будет находиться в состоянии («1„«1„..., «1к) вплоть до момента первого столкновения. Следовательно, для построения марковского случайного процесса, заключающегося в переходах системы через последовательность состояний, нужно определить вероятности этих переходов для каждого состояния к = Рь (т = 1, 2 .....) или, что то же самое, вероятности столкновений различных пар частиц 1, ! = 1, ..., Л'.
Необходимо также ввести временную шкалу случайного процесса. Вероятность того, что время свободного пробега 1-й частицы больше 1, равна Р (11 ) т) = ехр ( — о«/и1(!). (9) Согласно формуле (9), можно в принципе разыграть времена 11 для всех 1 =. 1, ..., Х и выбрать минимальное из этих времен *. Ясно, что ! ~„является случайной величиной, закон распределения которой будет получен ниже. Пусть Р (!мсо (!) — вероятность того, что г,о„меньше !. Тогда Розыгрыш 1 по формуле (12) определяет время, через которое в системе происходит первое столкновение е. Техг самьыг определяется физическая временная шкала случайного процесса.
Пз выражения (13) следует, что величина п — бп / ог — ез ! Л Рп =- <> является вероятностью столкновения пары частиц (1,7) в заданном состоянии системы (ид, и„..., нм). Набором этих вероятностей ..., Лг; 7 = 1, ..., Лг и вероятностей изменения скоростей сталкивающихся частиц Р (тгг, тгт — ь и,', тг,:) определяется набор вероятностей перехода системы из данного состояния в другие, отличающиеся скоростями двух частиц (г, 1). Отметим, что в работе (17) монтекарловская схема, до некоторой степени аналогичная описанной выше, применялась для изучения релаксации пространственно-однородной системы молекул к максвелловскому распределению. Однако «номера» сталкивающихся молекул разыгрывались не по вероятностям типа выражений (14), а равномерным обрааом.
По этой причине схема, описанная в работе (17), неприемлема для изучения временного поведения системы. Перейдем теперь к рассмотрению практической реализации построенного случайного процесса на ЭВМ. Пусть в начальный момент времени задано начальное состояние системы (т.
е. функция распределения частиц по скоростям). Нас будет интересовать эволюция этого распределения во времени, поскольку знание функции распределения в любой момент времени дает возможность определения всех величин, представляющих физический интерес. В основных чертах блок-схема реализации случайного процесса перехода системы из одного состояния в другое сводится к следующему. 1) Блок 1 присваивает всем частицам г = 1, ..., Дг начальные скорости тгз, согласно заданному начальному распределению.
Обычно таким распределением является распределение Максвелла / (и) =- Аиз ехр ( — тих /2ЙТ), (15) где А — нормировочная константа; и — масса частиц; Т вЂ” температура; (г — постоянная Больцмана. Существует точный способ розыгрыша скоростей максвеллов- ского распределения частиц (18!. Согласно этому способу, ком- * Следует отметить, что прп описанном видоиамененни метода периодических граничных условий учитываются не только мелкомасштабные флуктуации плотности, как в работах [11, 12), но и крупномасштабные.
Этот учет происходит автоматически при розыгрыше времени, согласно формуле (12). понеыты скорости и одной частицы равны и» = У вЂ” 2о» !и 5, сое (2я$») г„= г' — 2з» ! п $, з'и (2я$») уг =-. У вЂ” 25» !и ~» соз (2я~)) (16) где ое ==-.
кТ/т; $„з„$», ь, — «лучайные числа, равномерно распределенные в интервале (0,1). По формулам (16) разыгрываются скорости для частиц разного сорта. 2) 1)лок 11 вычисляет начальные величины Р;;, определенные выражением (14), для всех пар (1, 1) и величину П, согласно вьтражепню (!3).
Блоки 1 и 11 задают начальное состояние системы, вычисления в рамках этих блоков проводятся однократно. 3) В блоке 11! производится розыгрыш физического времени, согласно формуле (12). Величина г заносится в специальный счетчик Хо После этого осуществляется розыгрыш номеров (1, 1) столкнувшихся частиц. Для этого последовательно вычисляются суммы вида яз а,„= ~~~~ Р„„ »=1 где каждый номер к соответствует одной из пар частиц (1, у), а Р— вероятности, определяемые, согласно формуле (14), и хранящиеся в памяти машины. Считается, что при а ., ~( $ ( а сталкивается пара (», 1), соответствующая номерут. 4) В блоке 1Ч разыгрывается элементарный акт столкновения частиц (1, /) с начальными скоро«тимин;, и;.
При розыгрыше надо учесть зависимость сечения от угла рассеяния. Если моделями частиц служат твердые невзаимодействующие сферы, то в системе центра масс рассеяние иаотропно. В результате столкновения частицам (», 1) присваиваются новые скорости и,', и,'.
5) В результате столкновения частиц (1, 1) система перешла в новое состояние. В соответствии с этим в блоке Ч дычисляются новые величины Рпп Р;» (к = 1, ..., Х) и вносятся соответствующие изменения в величину й. На этом один шаг в марковской цепи заканчивается, а управление передается блоку 1П. Знание мгновенной функции распределения частиц по скоростям позволяет вычислять средние значения энергии, скорости в любой момент времени.
Кроме того, кри реализации случайного процесса могут быть получены такие величины, как частота столкновений и средняя передача энергии за столкновение. Выводя все эти результаты на печать через определенный временной шаг, можно получить детальное описание кинетики процесса. Специального рассмотрения требует вопрос о точности получаемых результатов и затратах машинного времени. Точность определяется числом рассматриваемых частиц Аг е. Если бы даже и существовала возможность неограниченного увеличения числа частиц, это вряд ли было бы целесообразным. Дело в том, что точность результатов растет, как ~г У, а затраты машинного времени в такой схеме — как Ага. Поэтому в большинстве расчетов, проведенных нами по описанному методу, система состояла нз 108 частиц.
Поскольку при таком значении Аг результаты оказывались недостаточно точными (ошибка около +-10е~о), реализовывалось несколько статистически независимых цепей Маркова с одинаковым начальным состоянием системы. Это достигалось путем проведения серии аналогичных расчетов с разными начальными случайными числами. Число цепей составляло от 20 до 60. При атом ошибка усредненных результатов колебалась от -1-1 до ~-3%. Аналогичный прием использовался ранее в ряде работ (см., например, )8, 10)).