CC BY
66
ТЕСТИРОВАНИЕ ГЕНЕРАТОРОВ СЛУЧАЙНЫХ ЧИСЕЛ С ПОМОЩЬЮ ЧИСЛЕННОГО МОДЕЛИРОВАНИЯ ТОЧНО РЕШАЕМОЙ ЗАДАЧИ
B. П. Мемнонов
C.-Петербургский государственный университет, канд. физ.-мат. наук, ст. научн. сотрудник, pokusa@star.math.spbu.ru
1. Введение
Статистические методы Монте-Карло широко применяются для численного моделирования различных задач в науке и технике. Существенной частью этих методов является интенсивное использование псевдослучайных последовательностей для моделирования различных статистических распределений и исходов событий. Обычно эти последовательности вырабатываются специальными компьютерными программами, генераторами случайных чисел (ГСЧ). Поэтому они являются только некоторыми приближениями к последовательностям настоящих случайных чисел. Степень этого приближения проверяется с помощью различных локальных и глобальных [1] тестов, которые хороший генератор должен выдерживать с необходимостью. В то же время как достаточные в некотором практическом смысле можно рассматривать специальные проверки путем решения типичных задач, допускающих независимую оценку результатов с помощью аналитических или численных методов и относящихся к области дальнейшего применения этих ГСЧ [2-3]. Причем для параллельных программ, запускаемых на большом числе процессоров, нужно столь же много по возможности независимых последовательностей. В работе [4], а затем и в [5] доказывается, что в генерируемых псевдослучайных последовательностях имеются корреляции дальнего действия, которые могут влиять отрицательно на их качество при использовании отдельных частей псевдослучайной последовательности одного и того же генератора для работы на разных процессорах в параллельных программах. В работе [6] показано, что для оценки степени независимости между псевдослучайными последовательностями, получаемыми с помощью '^сЬшапп-НШ генераторов, может быть использован спектральный тест [7]. В заметке [8] указаны способы получения воспроизводимых псевдослучайных последовательностей для использования в многопроцессорных компьютерах. Некоторые предложения по снижению влияния корреляций для генератора специального вида высказаны в статье [9]. Дополнением к этим исследованиям непосредственно генерируемых последовательностей может быть численное моделирование точно решаемой задачи с использованием испытуемых генераторов, разных на каждом из процессоров работающей параллельной программы. Сравнение результатов с точным решением позволит в этом случае оценить степень влияния качеств генераторов для практического применения в некотором круге проблем, так или иначе связанных с этой задачей. С этой целью в работе методом прямого статистического моделирования Монте-Карло (ПСМ) решается задача о релаксации вращательной энергии молекул из первоначального неравновесного состояния посредством межмолекулярных столкновений, для расчета которых применяются псевдослучайные последовательности испытуемых генера-
© В. П. Мемнонов, 2009
торов. В конечном, равновесном состоянии должно установиться распределение Больцмана с известной и одинаковой для всех степеней свободы температурой Т. Именно по величине отклонений от этого теоретического значения с учетом ошибок самого метода ПСМ, в наибольшей степени связанных со статистическим рассеянием из-за неизбежной ограниченности используемых статистических выборок, и оценивалось качество тестировавшихся мультипликативных генераторов. Также по выборочным коэффициентам корреляции оценивались истинные коэффициенты корреляции между конечными температурами, полученными на разных процессорах и в случае, когда на каждом процессоре производились дополнительные выборки в конечном состоянии для увеличения общей выборки. Кстати наличие подобных корреляций и их возможное влияние на точность получаемых результатов не учитывалось в работе [10], в которой изучалась зависимость статистических ошибок метода ПСМ от параметров в некоторых прикладных задачах.
2. Постановка задачи
Рассматриваются молекулы, равномерно распределенные в объеме с зеркально отражающими стенками и имеющие только поступательные и вращательные степени свободы. Причем в начальный момент создается неравновесное состояние системы, когда вся энергия сосредоточена только на поступательных степенях свободы. В конечном состоянии молекулы должны иметь распределение Больцмана для всех степеней свободы. Отклонения поступательной Т4г и вращательной Тг температур от известного теоретического значения в равновесии Т могут существовать, но при идеальных генераторах это возможно только вследствие ограниченного размера статистической выборки N. Связанную с этим ошибку гN метода Монте-Карло [11], например, для вращательной температуры 6ТГ можно записать в виде
гм = х(В(Тг )1/2, (1)
где В(ТГ) —дисперсия вращательной температуры, N — величина выборки, а параметр х принимался равным 3, соответствуя коэффициенту доверия 0.997. Поэтому, если в результате релаксации со статистической выборкой N получаются отклонения от Т4,
5ТГ = | Тг — Тг |, (2)
больше ГN, то это значит, что испытуемый генератор с вероятностью 0.997 плохой.
2.1. Параметры задачи и генератора, модель взаимодействия молекул
В настоящей работе последовательная программа Берда [2] для расчета вращательной релаксации в однородном газе была распараллелена с помощью технологии МР1 и в дальнейшем использовалась только эта параллельная версия. Для расчета вращательной релаксации использовалась модель Ларсена—Боргнакке (ЛБ) [12], как и у Берда. После выполнения распараллеливания и добавления в программу еще специально разработанных процедур для оценки корреляций получился удобный инструмент испытания генераторов для параллельных вычислений. Тестировались линейные мультипликативные генераторы вида
Zi+1 = А * Zi (шо(ШоМ), М = 231 — 1, (3)
с различными множителями А у. Начальное значение Zо полагалось равным 1, а множители Ау брались из работы [1], где представлены 207 множителей, прошедших всесторонние теоретические тесты. Для моделирования с помощью методов Монте-Карло
последовательность Zi нормируется и превращается в последовательность U,
U =[Ui = Zi/M; i = 1, 2,...], (4)
элементы которой представляют собой как бы результат моделирования из равномерного распределения на интервале [0,1).
Для использования такого генератора с разными множителями Aj в параллельной программе нужно в рабочей подпрограмме иметь массив с этими множителями и каждый процесс, используя свой номер rank в коммуникаторе, выбирает себе по соответствующей позиции в этом массиве уникальный множитель. Заметим, что с помощью введения в генератор счетчика вызовов его программой и условия перехода его на другой множитель при приближении к периоду генератора M — 1 можно получать последовательности гораздо большей длины. Действительно, в рамках интерфейса MPI процессам известно общее число запущенных процессов Nproc и, применяя это число и номер процесса rank, можно для каждого процесса последовательно выбирать уникальные позиции Pk,
Pk = rank + 1 + (NProc — 1) * k, (5)
в массиве до полного его исчерпания. В принципе таким же путем можно получать и использовать массив гораздо большего размера, если в каждой его ячейке сохранять не только множитель Aj, но и начальные значения более коротких подпоследовательностей, на которые заранее будет разбиваться последовательность псевдослучайных чисел, связанная с каждым множителем Aj.
Относительно взаимодействия частиц предполагалось, что молекулы сталкиваются как твердые сферы, а в неупругих вероятностях вращательных переходов при этом используется значение 5 для параметра числа столкновений, что соответствует вероятности неупругого перехода 0,2. Расчетное поле разбивается на одномерные ячейки, в которых производятся процедуры столкновений метода ПСМ [2]. Временной шаг т, так же как ив [2], был равным половине среднего времени между столкновениями для одной молекулы, а вся расчетная область содержала Nm молекул, причем во всех расчетах Nm равнялось 100001. После достижения равновесного состояния через промежутки времени At,
At = Ntf х т, (6)
выполнялись еще Kts дополнительных выборок параметров системы, причем минимальное число Ntf временных шагов т между выборками необходимо как-то оценивать. В практических расчетах часто нужны подобные добавочные осреднения для увеличения статистической выборки. Поэтому, имея в данной работе возможность оценки как возникающих при этом кореляций, так и ошибок численного решения путем сравнения с точным, в дальнейшем проведем исследование условий, когда это действительно можно делать. В итоге суммарный объем статистической выборки N получается тогда равным:
N = Nm х Kts х Nproc, (7)
где Nproc —число используемых процессоров параллельного кластера.
В начальном состоянии, как уже говорилось выше, вся энергия была сосредоточена только на поступательных степенях свободы молекул. Поступательная температура Ttr полагалась равной 500K, вращательная же Tr была равна нулю. Предполагалось, что молекулы обладают только двумя вращательными степенями свободы. В модели ЛБ строго выполняются законы сохранения энергии и импульса и условие детального
баланса [12]. Поэтому в процессе вращательной релаксации должно устанавливаться конечное равновесное состояние с одинаковыми значениями температуры Т, равными 300К, как это следует из закона сохранения полной энергии для данного случая. Несмотря на то, что модель ЛБ только приближенно описывает во времени вращательную релаксацию, конечное равновесное состояние получается точным. И сравнение результатов численного моделирования проводится именно по этому конечному состоянию.
3. Оценка свойств независимости с помощью перехода к нормальным случайным величинам
В практических расчетах из-за ограниченности числа доступных процессоров на параллельном кластере часто приходится в конечном состоянии на каждом процессоре использовать еще эти К.^а дополнительные выборки параметров задачи, поэтому очень важно установить величину промежутков между ними или Ntf в единицах т так, чтобы они были, по возможности, независимы. Также очень существенно оценить, насколько независимы реализации, получаемые с помощью ГСЧ с разными множителями А у. Для этого можно использовать оценки коэффициентов корреляции р с помощью выборочных коэффициентов корреляции К, получаемых из результатов численного моделирования нашей задачи. При этом более надежные результаты получаются, если рапределение рассматриваемых величин достаточно близко к нормальному [13]. Дело в том, что в отличие от других наиболее употребительных распределений, только в случае двумерного нормального распределения обращение в нуль коэффициента корреляции эквивалентно независимости рассматриваемой пары случайных величин [14]. Поэтому ниже производится переход к нормированным суммам вращательных температур, а затем с помощью неравенства Бери—Эссеена [15] оценивается близость их распределения к нормальному закону.
В распределении Больцмана, например, для вращательной энергии Ег удобно перейти к переменной Тг = Ег/х, где х — постоянная Больцмана, и получить экспоненциальную плотность функции распределения ] (Тг):
/ (Тг) = ехр(-Тг/Тг)/Тг, (8)
где Т — получаемое из закона сохранения полной энергии значение конечной температуры, равное 300 К, в рассматриваемом случае молекул с двумя вращательными степенями свободы. При этом
Е (Тг) = Ти (9)
В(ТГ) = Т2, (10)
е(тЗ) = 6Т3. (11)
Каждая из Nm молекул обладает одинаковой плотностью распределения (8) и всеми
моментами. Пусть Рнт обозначает функцию распределения нормированной суммы SNm
для Nm случайных величин — вращательных температур Тгу, записанную с учетом (9) и (10) в виде
у1*™ (Тгу — Т)
= Д^ (12)
т Тгл/Йш
Обозначим также через Ф стандартную нормальную функцию распределения:
1 Г X
Ф(ж) = __ е~г /2сЫ. (13)
Согласно центральной предельной теореме последовательность функций распределения F^m таких нормированных сумм S^m при Nm ^ ж равномерно сходится к стандартной нормальной функции распределения Ф(ж):
e(FNm , Ф) = SUP lFNm (x) - Ф(х)| ^ 0. (14)
X
При конечных значениях Nm с помощью неравенства Бери—Эссеена [15] для ошибки нормальной аппроксимации f3(FNm, Ф) в случае конечного третьего момента (11) и с учетом (10) можно получить оценку [16]
/?(^т,Ф)<4=, (15)
в которой константа С\1 по-видимому, близка к значению 1/а/27г [16]. Тогда для Nm = 100001 получаем
l3(FNm, ф) < 0, 757 * 10-2. (16)
Такая точность аппроксимации достаточна для дальнейших оценок корреляции. Заметим, что на нее была потрачена весьма существенная часть общей выборки (7).
3.1. Коэффициенты корреляции
На каждом процессоре производится осреднение конечных температур поступательного Ttr и вращательного Tr движений всех Nm молекул и, тем самым, получаются приближенно нормально рапределенные случайные величины, например, для вращательных температур. Обозначим их Tr. Их общее число NK равно Kts х Nproc. С помощью операции MPI — Reduce они собираются последовательными частями от каждого процессора размером Kts в специальном массиве, будем называть его массивом пк. В дальнейшем разные их половины величиной П2 = NK/2 можно рассматривать как выборки для двух нормально распределенных случайных величин X и Y. Затем с помощью их выборочного коэффициента корреляции (В.К.К.) rxY оценим истинный коэффициент корреляции pxy, а вместе с ним зависимость или независимость этих последовательностей и порождающих их генераторов. Итак, определим два типа В.К.К. по позиции j для Tr в одномерном массиве пк. Сравнивая попарно четные и нечетные позиции в нем, получим rodev:
Г- = (17,
Sx * Sy
Здесь выборочная ковариация COV (X,Y) дается формулой
E”21 (Tr (2j — 1) — Xav) * (Tr (2j) — Yav)
COV(X, Y) = ^0=lK y J K K J------—, (18)
п2
представленной через выборочные дисперсии SX2 и s2 ,
C2 _ £”2l(Tr (2j — 1) — Xav )2 C2 _ j (Tr (2j) — Yav )2
п2 — 1 п2 — 1
и выборочные средние Xav и Yav,
„ E”2i % (2j — 1) E”2i % (2j)
п2 п2
где j меняется от 1 до П2.
(20)
Когда на каждом процессоре используются дополнительные выборки и К^ > 1, причем Къа четное число, выборочный коэффициент корреляции То^еу оценит для ближайших соседних из них, разделенных промежутком (см. (6)), возможные попарные статистические связи. Если же К^ = 1, т. е. выполняется только одна выборка на процессоре, то можно оценить статистические связи между случайными величинами Тт, порождаемые попарно соседними Г.С.Ч. с множителями Л у, расположенными в некотором специальном порядке, например, как в статье [1].
В другом типе В.К.К. тъе, рассматриваемом также в дальнейшем, используются Тт из верхней и нижней половин одномерного массива пк попарно:
= ^!Ы). ,21,
ЭХ * Эу
В этом случае выборочная ковариация СОУ(Х,У) дается формулой
ЕП21 (Тг(з) - хау) * (Тт(мк - з + 1) - Уау)
СОУ(Х, У) = ^3=1[ Шк к-----------------——-------------------------------—, (22)
П2
представленной через выборочные дисперсии и
о2_ТГМ^)-х^? с2 _Т^ит^к-з + 1)-га^
п2 - 1 п2 - 1
и выборочные средние Хау и Уау,
ь.аЁЕ±», (м,
п2 п2
где ] также меняется от 1 до п2. Отличное от 1 число Ка, как уже говорилось выше,
полагалось четным. Четным же полагалось и число использумых процессоров Мртос.
В таком случае случайные величины X и У для В.К.К. тье принадлежат разным процессорам и, соответственно, разным множителям Лу.
Истинный коэффициент корреляции рху между случайными величинами X и У можно оценить с помощью В.К.К. тоаеу или тье. Для больших выборок среднее квадратичное отклонение аг коэффициента, например, то^еу от рху, согласно [13] может быть представлено в виде
1 Т2
огк,8г =-----2^, (25)
л/П2
при этом тоаеу приближенно следует нормальному закону с параметрами (рху ,&т). И для уровня значимости q процентов, имея неравенство
Тойеу * &т < рХУ < ^ойеу + tq * &т ? (26)
где для q =1 и П2 > 100, по [13] получаем tq = 2, 58. Таким образом, при выполнении неравенства
1 Тойеу |> ~tq * ®т (27)
коэффициент корреляции рху отличен от нуля и между рассматриваемыми случайными величинами есть значимая статистическая связь.
4. Результаты численного моделирования
Численное моделирование вращательной релаксации выполнялось на 40 процессорном кластере Санкт-Петербургского университета. Однако практически можно было использовать одновременно не больше 16 процессоров. Поэтому в приводимых ниже данных часто сначала накапливалась необходимая выборка шагами по 16, а затем она обрабатывалась специальной программой. В частности, именно таким путем были получены результаты для 204 процессоров, представленные на рис. 1 при использовании только одной выборки Kts = 1 на каждом из них. Кривые 1 и 2 там изображают зависимость от числа «истраченных» случайных чисел генератора в виде доли Pp его периода M — 1 (в рассматриваемом случае M = 2147483647) для отношений RI абсолютных величин коэффициентов корреляции Todev и The к своим граничным значениям rbodev = tq *<Jr и rbbe. Последняя получается аналогично предыдущей с заменой в формуле (25) для ar коэффициента rodev на Tbe. Причем приближенно нормальные случайные величины Tr располагались так же, как и соответствующие множители Aj в статье [1], т. е. в порядке убывания для них значений величины Si,k. Эта величина характеризует максимальные расстояния между гиперплоскостями в единичном гиперкубе размерности к. И множитель Aj, который соответствует минимальному значению этого максимальные расстояния, обеспечивая меньшую величину пустых областей, более полезен для использования в генераторе. Это и послужило в [1] основным критерием для выбора приведенных множителей. Одновременно там отмечается, что этот критерий соответствует спектральному тесту. Как видно из рис. 1 обе кривые 1 и 2 располагаются ниже горизонтальной прямой при RI = 1, которая служит границей появления значимости для истинного коэффициента корреляции р (см. раздел 3.1). Поэтому оцениваемых таким образом корреляций между соседними или, наоборот, последовательно крайними в этом порядке, не существует. Кривые 3 и 4 на рис. 1 представляют отношения абсолютных величин отклонений поступательных 5Ttr = |Ttr — Tt | и вращательных 5Tr = |Tr — Tt| температур от теоретического равновесного значения Tt к ошибке tn метода Монте-Карло (1). Эти кривые также располагаются гораздо ниже 1 при всех Pp, обеспечивая надежность одновременного применения 204 генераторов с такими множителями.
На рис. 2 представлены результаты, полученые для 16 процессоров с дополнительными 14 выборками на каждом из них. В этих группах последовательно в порядке
п к
Рис. 1. Отношения Ш в зависимости от по- Рис. 2. Коэффициенты |го^еу | —3, \гье | —5 и
траченной части периода Рр генератора. гЬо^е« —4, гЬье —6; отклонение \5Т^,г|тах — 1 и
ошибка гN —2 для разных Б\,к•
согласно статьи [1] используются по 16 множителей Лу. Для каждой группы таких множителей вычисляется среднее значение величин и именно оно, 6\,к, откладывается по оси абсцисс. Сплошная кривая 1 изображает масимальное отклонение поступательных или вращательных температур ^Т^^ |тах от теоретического равновесного значения Т(, а прямая 2 показывает статистическую ошибку тN = 0.19 (1) метода Монте-Карло при постоянной для каждой группы суммарной выборке N = 224 * 10001. Начиная с середины, точнее с 65 номера в последовательности [1], на кривой 1 появляются участки, заметно превосходящие эту границу, и соответственно эти группы множителей не следует использовать в методе ПСМ с дополнительными выборками на каждом процессоре. Заметим, что они же при К^ = 1 не портили результаты на рис. 1. Кривая 3 на рис. 2 изображает абсолютную величину выборочного коэффициента корреляции \rodev | и повсюду превосходит кривую 4, его граничное значение тЬо^еу. Это означает, что попарные «четные — нечетные» корреляции, т. е. корреляции между соседними дополнительными выборками на каждом из 16 процессоров в данном случае имеют место, несмотря на сравнительно большой применявшийся между 14 выборками в конечном состоянии интервал N1/, равный 24 средним временам между соударениями одной молекулы 2т. Действительно, за это время молекула при заданной вероятности неупругого соударения 0,2 в среднем 5 раз изменяет свою вращательную энергию. Кривая 5 изображает абсолютную величину выборочного коэффицента |тье | и почти нигде не превосходит свою границу — кривую 6. Опять имеем попарную независимость между случайными величинами, получаемыми с помощью 16 генераторов с разными множителями Л у .
На рис. 3 более подробно расматриваются результаты для двух первых множителей Л1 = 742938285 и Л2 = 950706376. В этом случае использовались 112 дополнительных выборок в конечном состоянии и исследовался вопрос о влиянии на результат величины интервала N1/ между выборками. Во всех рассмотренных случаях здесь поддерживалось одно и тоже граничное значение достижения равновесного состояния, равное 1920 итерациям или 960 средним временам между соударениями одной молекулы.
Здесь кривая 1 представляет отношение Ш абсолютной величины выборочного коэффициента |rodeV | к его граничному значению тЬ^еи и, как видно из рис. 3, вплоть до значений N1/ = 12 превосходит 1. Так что попарные корреляции между величинами, полученными с помощью дополнительных выборок, на этом интервале существуют, но после этого значения они исчезают, в отличие от соответствующих результатов на рис. 2, где они есть и при N1/ = 24. Может быть влияют остальные, несколько худшие 14 множителей, которые участвуют в эксперименте помимо первых двух. Причина такого расхождения пока неизвестна. Таким образом, только для первой трети множителей статьи [1], несмотря на наличие корреляций, с некоторой осторожностью можно производить дополнительные выборки.
Кривая 2 представляет отношение абсолютной величины выборочного коэффициента |гье | к его граничному значению тЬье и все время меньше 1, еще раз подтверждая попарную независимость случайных величин, получаемых с помощью генераторов с разными множителями. Кривые 3 и 4 представляют также, как и на рис. 1, отношения абсолютных величин отклонений поступательных 5Т1Г = — Т^ и вращательных
5ТГ = |ТГ — Т^ температур от теоретического равновесного значения Т к ошибке TN метода Монте-Карло. За исключением начального участка с N1/ = 1 обе кривые меньше
1, обеспечивая хорошее качество этих множителей.
На рис. 4 расматриваются результаты для двух множителей Л97 = 964028288 и Л98 = 1493834601, которые входят в «плохую» группу на рис. 2.
О 5 10 15 20 25
№
Рис. 3. Отношения Ш1 коэффициентов | |, |гЬе| к гЬойеу
и гЬье —кривые 1 и 2, и отклонений 1$^г | и 15ТГ | к ошибке — кривые 3 и 4 в зависимости от интервала Ntf между выборками.
Рис. 4. Отношения В1 коэффициентов | Го^е-о |, |гье| к гЬойещ и гЬье — кривые 1 и 2, и отклонений |^Т;Г| и 15ТГ| к ошибке г^ — кривые 3 и 4 в зависимости от интервала Ntf между выборками.
Как видно на всем интервале значений N4/ кривая коэффициента |г0^е« | располагается выше 1, свидетельствуя о попарной статистической зависимости в дополнительных выборках. При этом Монте-Карловские ошибки \5Т^Г | и \5ТГ | также больше допустимого уровня. Поэтому для этих множителей нельзя пытаться увеличивать статистическую выборку с помощью дополнительных выборок в конечном состоянии. Возможно, в данном случае это связано с какой-то спецификой структуры гиперплоскостей в многомерном пространстве, на которые попадают все псевдослучайные числа в таких генераторах.
5. Заключение
Проведенное численное моделирование точно решаемой задачи позволило оценить качество рекомендованных в работе [1] множителей в линейном мультипликативном генераторе на всем протяжении его периода в отличие от ограниченного интервала для выполненных тестов. Отмечены некоторые группы множителей, использование которых в методе ПСМ может привести к ошибкам. Выполненные численные эксперименты устанавили попарную независимость последовательностей случайных величин, порождаемых мультипликативными генераторами с разными множителями, приведенными в [1], и, наоборот, выявили существование попарной статистической связи при использовании дополнительных выборок на каждом процессоре, которая может сопровождаться большой ошибкой и поэтому ее не следует делать для множителей [1], начиная с 65.
Литература
1. Fishman G. S., Moore L. R. An Exhaustive Analysis of Multiplicative Congruential Random Number Generators with Modulus 231 — 1 // SIAM J. Sci. Stat. Comput. 1986. Vol. 7, N 1. P. 24-45.
2. Bird G. A. Molecular Gas Dynamics and the Direct Simulation of Gas Flows. Oxford, Clarendon Press, 1994. 458 p.
3. Ермаков С. М., Михайлов Г. А. Статистическое моделирование. М.: Наука, 1982. 296 с.
4. De Matteis A., Pagnutti S. Parallelization of Random Number Generators and Long-Range Correlations // Numer. Math. 1988. Vol. 53. P. 595-608.
5. Eichenauer-Herrman J., Grothe H. A Remark on Long-Range Correlations in Multiplicative Congruential Pseudorandom Number Generators // Numer. Math. 1989. Vol. 56. P. 609-611.
6. Maclaren N. M. Generation of Multiple Independent Sequences of Pseudorandom Numbers // Appl. Statist. 1989. Vol. 38, N2. P. 351-359.
7. Кнут Д. Искусство прогарммирования для ЭВМ / D.E. Knuth. The art of computor programming Vol. 2 / Пер. с англ. Под ред. К. И. Бабенко. М.: Мир, 1977. 724 с.
8. Koniges A. E., Leith C. E. Parallel Processing of Random Number Generation for Monte Carlo Turbulence Simulation // J. Comput. Phys. 1989. Vol. 81. P. 230-235.
9. Brody T. A. Random Number Generation for Parallel Processors // Comput. Phys. Commun. 1989. Vol. 56. P. 147-153.
10. Chen G., Boyd I. D. Statistical Error Analysis for the Direct Simulation Monte Carlo Technique // J. Comput. Phys. 1996. Vol. 126. P. 434-448.
11. Соболь И. М. Численные методы Монте-Карло. М.: Наука, 1973. 311 с.
12. Borgnakke C., Larsen P. S. Statistical collision model for Monte Carlo simulation of polyatomic gas mixture // J. Comput. Phys. 1975. Vol. 18. P. 405-420.
13. Смирнов Н. В., Дунин-Барковский И. В. Курс теории вероятностей и матеатической статистики. М.: Наука, 1969. 511 с.
14. Кендал М., Стьюарт А. Статистические выводы и связи / M. G. Kendal, A. Stuart. The Advance Theory of Statistics. Vol. 2. Inference and Relationship, sec. ed. London, 1969 / Пер. с англ.; Под ред. А. Н. Колмогорова, М.: Наука, 1973. 899.
15. Феллер В. Введение в теорию вероятностей и ее приложения / An Introduction to Probability Theory and its Applications / Пер. с англ.; Под ред. Ю. В. Прохорова. М.: Мир, 1984. Т. 2. 751 с.
16. Королев В. Ю., Шевцова И. Г. О точности нормальной аппроксимации // Теория веро-ятн. и ее применения. 2005. Т. 50. Вып. 2. С. 353-366.
Статья поступила в редакцию 18 сентября 2008 г.
