Метод Монте-Карло на графических процессорах

Министерство образования и науки Российской Федерации
Уральский федеральный университет
имени первого Президента России Б.Н. Ельцина

К.А. Некрасов, С.И. Поташников, 
А.С. Боярченков, А.Я. Купряжкин

Метод Монте-Карло 
на графических процессорах

Учебное пособие

Рекомендовано методическим советом УрФУ
для студентов, обучающихся по направлениям подготовки
14.04.02 — Ядерная физика и технологии;
09.04.02 — Информационные системы и технологии;
14.04.01 — Ядерные реакторы и материалы

Москва
Издательство «ФЛИНТА»
Издательство Уральского университета
2017

2-е издание, стереотипное

УДК 519.245:004(075.8)
ББК 22я73+32.97я73
          М54
Авторы:
Некрасов К.А., Поташников С.И., Боярченков А.С., Купряжкин А.Я.

Рецензенты:
Институт теплофизики УрО РАН (д‑р физ.‑мат. наук, проф.В.Г. Байдаков); гл. науч. сотр. лаборатории математического моделирования Института промышленной экологии УрО РАН д‑р физ.‑мат. наук, проф. 
А.Н. Вадаксин

М54     Метод Монте-Карло на графических процессорах [Электронный 
ресурс]: учебное пособие / К. А. Некрасов, С. И. Поташников, 
А. С. Боярченков, А. Я. Купряжкин. — 2-е изд., стер. — М. : ФЛИНТА : 
Изд-во Урал. ун-та, 2017. — 58 с.

ISBN 978-5-9765-3183-3 (ФЛИНТА)
ISBN 978‑5‑7996‑1723‑3 (Изд-во Урал. ун-та)

В учебном пособии изложены основные идеи метода Монте‑Карло, принципы 
распараллеливания расчетов, организации высокоскоростных параллельных вычислений методом Монте‑Карло на графических процессорах NVIDIA архитектуры CUDA. Подробно проанализирован пример моделирования методом Монте‑Карло диффузии нейтронов через пластину.
Пособие предназначено для проведения практических занятий по программированию графических процессоров для магистрантов по направлениям подготовки 09.04.02, 14.04.02, специалистов по направлению подготовки 141401.

Библиогр.: 15 назв. Рис. 15. Прил. 1.

УДК 519.245:004(075.8)
ББК 22я73+32.97я73

© Уральский федеральный
     университет, 2016
ISBN 978-5-9765-3183-3 (ФЛИНТА)
ISBN 978‑5‑7996‑1723‑3 (Изд-во Урал. ун-та)

Введение

Одним из широко используемых подходов к численному моделированию физических систем является метод Монте‑Карло [1]–[2], для 
применения которого необходим анализ большого количества случайных состояний и вариантов поведения исследуемой системы, что связано с большими объемами вычислений.
Варианты случайного поведения системы, рассматриваемые методом Монте‑Карло, обычно генерируются и обрабатываются по одному и тому же общему алгоритму, будучи при этом независимыми друг 
от друга. Это позволяет эффективно реализовывать метод Монте‑Карло на поточно‑параллельных вычислительных системах.
В последние годы появилась практическая возможность проведения поточно‑параллельного физико‑математического моделирования 
на графических процессорах (GPU) — вычислительных устройствах 
параллельной архитектуры, изначально предназначенных для отображения графики на персональных компьютерах и рабочих станциях [3].
В настоящее время ведущие производители графических процессоров сами позиционируют новые GPU как универсальные системы, 
предназначенные не только для графических вычислений, но и для расчетов общего назначения. В частности, компанией NVIDIA выпускаются графические процессоры архитектуры CUDA [4] и соответствующее программное обеспечение [5], которые позволяют сравнительно 
легко программировать GPU на языке, являющемся расширением 
языка C, сохраняющем синтаксис последнего и включающем все его 
основные возможности.
Современные графические процессоры имеют в своем составе 
до нескольких сотен параллельных процессоров, каждый из которых 
по вычислительной производительности сравним с центральным процессором персонального компьютера [6]– [7].

ВВедение

Опыт показывает, что при решении хорошо распараллеливаемых 
задач реальная производительность GPU оказывается не менее чем 
на два порядка выше, чем у современных центральных процессоров 
ПК. Отметим, что затраты на приобретение и эксплуатацию одного 
вычислительного комплекса с GPU примерно совпадают с затратами на обычный персональный компьютер, тогда как кластер эквивалентной производительности занимал бы намного больше места, будучи при этом намного сложнее в эксплуатации и на порядки дороже.
В работах [8]–[10] представлены результаты применения графических процессоров для моделирования кристаллов диоксида урана 
методом молекулярной динамики. Использование графических процессоров позволило ускорить решение данной задачи (по сравнению 
с центральным процессором) более чем в 600 раз. Настоящее пособие 
посвящено оценке возможностей графических процессоров архитектуры CUDA (Compute Unified Device Architecture [3]–[5]) на примере метода Монте‑Карло, в частности при моделировании одномерной диффузии нейтронов через пластину.
Архитектуры графических процессоров, выпускаемых ведущими 
производителями — компаниями NVIDIA и AMD/ATI, несколько 
различны, что особенно проявляется при неграфических вычислениях. Наиболее удобной для таких вычислений мы считаем архитектуру 
CUDA [3]–[5], предлагаемую компанией NVIDIA.
В настоящем пособии будут рассмотрены особенности параллельной реализации метода Монте‑Карло для высокопроизводительных 
вычислений на графических процессорах. При этом мы ограничимся 
программированием для CUDA, так как в настоящее время эта платформа представляется наиболее эффективной, а реализованный в ней 
подход — наиболее перспективным.

1. Применение метода Монте-Карло  
для моделирования физических систем

1.1. Содержание метода Монте-Карло

1.1.1. Метод Монте-Карло на примере вычисления площадей

Простейшей иллюстрацией метода Монте‑Карло [1] является вычисление площади плоской фигуры (S) по принципу, показанному 
на рис. 1.1 [2].
Фигура заключена в квадрат, по которому 
случайным образом разбрасывается большое 
количество точек. Если число точек достаточно велико, то отношение площади S к площади квадрата L 2 будет с заданной точностью равняться отношению количества точек 
ў
N , 
попавших внутрь фигуры, к полному количеству точек N

S
L N

N
»
ў
2
.

Точность оценки площади возрастает вместе с количеством точек N. Погрешность вычислений обычно пропорциональна N ‑1/2 .
Площадь фигуры можно интерпретировать как ее объем в двухмерном пространстве. Очевидно, что метод, показанный на рис. 1.1, 
может быть применен и для вычисления объема трехмерной фигуры, 
заключенной внутрь куба. Аналогично можно рассчитывать объемы 
многомерных фигур в пространствах бόльшей размерности. Этот подход может быть без модификации применен и для многомерного численного интегрирования, поскольку определенные интегралы могут 
быть представлены как объемы многомерных тел.

L

S

Рис. 1.1. Вычисление 
площади плоской фигуры с криволинейной 
границей методом Монте‑Карло [2]

1. Применение метода монте-Карло для моделироВания физичесКих систем 

1.1.2. Распределения случайных величин

Одним из основных этапов любой вариации метода Монте‑Карло является генерирование случайных величин, распределенных 
по необходимому закону. При вычислении площади фигуры в параграфе 1.1.1 этими случайными величинами были координаты точек, 
разбрасываемых по площади квадрата. Каждой точке соответствовали две координаты (x, y), каждая из которых должна была быть случайной величиной, равномерно распределенной на интервале (0; L).
Пусть некоторая случайная величина x определена на интервале 
(a; b). Плотностью распределения x называют функцию w (x), такую, что

W x
w x dx
P
x

a

x
( ) =
( )
=
<
(
)
т
x
 —

вероятность того, что x окажется меньше, чем x. Сам интеграл W(x) называется интегральной (кумулятивной) функцией распределения величины x, значения W(x) всегда принадлежат интервалу (0; 1). 
Поскольку P
b
g Ј
(
)  = 1, для любой плотности распределения выполняется условие нормировки

 
w x dx

a

b
( )
=
т
1. 
 (1.1)

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

w x
b
a
x
a b

x
a b
( ) =
О(
)

П(
)

м

нп

оп

1
, при
;

0,
при
;

;

.

Такое распределение имеет прямоугольную форму, как это показано для интервала (0; 1) на рис. 1.2.

Кроме равномерно распределенных, 
широко используются и часто рассматриваются нормально распределенные (распределенные по Гауссу) случайные веw(x)

1

x
0
1
0

Рис. 1.2. Плотность вероятности для равномерного распределения на интервале (0; 1)