Диссертация (1103257), страница 4
Текст из файла (страница 4)
Поэтому представляется важным метод вычисления кинетической энергии,не опирающийся на данную теорему. Прямое усреднение по формуле (29) "наивного" решеточного представления кинетической энергииm (xn+1 − xn )22a2(61)как известно расходится в непрерывном пределе a → 0. Далее будет показано, что "наивное" выражение для кинетической энергии (61) в действительности не является верным, ив правильное выражение решеточной наблюдаемой для кинетической энергии входит ещеодно слагаемое, сокращающее расходимость "наивного"среднего.Вычислим явно среднюю кинетическую энергию при помощи континуального интеграла (рассмотрим одномерный случай, обобщение тривиально).RZZDxK(x)e−S(x)/~11−Ĥτ/~hKi = R=hx|K̂e|xidx =K(x, x0 )hx0 |e−Ĥτ /~ |xidxdx0 .ZZDxe−S(x)/~(62)Последнее выражение сдержит ядро оператора кинетической энергии, которое, как можнопоказать, равноK(x, x0 ) =−~2δ(x − (x0 − ∆x)) − 2δ(x − x0 ) + δ(x − (x0 + ∆x))lim,2m ∆x→0(∆x)2(63)что при подстановке в (62) дает1 −~2hKi =Z 2mZ∂ 2 hx0 |e−Ĥτ /~ |xi ∂x2dx0 .(64)x=x0Выберем представление hx0 |e−Ĥτ /~ |xi на решетке следующим образом:Zm(x − xN )2−a m(x1 − x0 )200 −Ĥτ /~+ V (x ) + .
. . ++ V (xN ) .hx |e|xi = dx1 . . . dxN exp~2a22a2(65)Подставляя последнее выражение в (64) и беря производные, получаемZ~~aV 00 (x) V 02 (x)a2m(x − xN )21−S(x)/~0dx1 . . . dxe++ (x − xN )V (x) +−−.hKi =Z2a22a2m2m(66)Ясно, что в непрерывном пределе (a → 0) мы получаем следующее простое соотношение:hKi = h−m (xn+1 − xn )2~+ i.22a2a(67)Как уже говорилось, среднее от первого слагаемого в данной формуле расходится. Однако второе слагаемое явным образом выделяет и сокращает расходящуюся часть. Этопозволяет вычислять кинетическую энергию любых систем напрямую.3.4Периодические граничные условияПри использовании рассматриваемых методов для моделирования среды, то есть большогоколичества частиц, взаимодействующих друг с другом, при отсутствующем или относительно слабом внешнем потенциале, возикает вопрос о постановке и реализации в численных рассчетах граничных условий в пространстве.
Граничные условия типа жестких стенок радикально ухудшают сходимость к непрерывному пределу, поэтому представляетсяразумным использовать периодические условия, тем более, что они позволяют сформировать модель, более близкую к бесконечной среде. Реальное моделирование бесконечнойсреды, разумеется невозможно, в первую очередь потому, что конечная плотность в системах с отталкивающим потенциалом (такая система будет описана в главе 3) достигаетсятолько за счет граничных условий, которые не позволяют частицам разлететься на бесконечное расстояние.На (Рис.
2) изображено стремление на непрерывному пределу энергии свободной частицы, ограниченной жесткими стенками с расстоянием L = 1 между ними (основное состояние, черным) и с периодическими граничными условиями L = 2 (первое возбужденноесостояние, красным). Зеленым отмечен точный ответ π 2 /2 ≈ 4, 9 (в атомных единицах).76E543210,000,020,040,060,080,100,120,14aРис. 2: Непрерывный предел в эквивалентных задачах при различных граничных условияхВ данном параграфе будет описан метод корректной и быстрой генерации траекторий в пространстве с периодической топологией.
Надо отметить, что именно изменениетопологии пространства и, следовательно, возможных траекторий – причина некоторойнетривиальности при введении граничных условий. Рассмотрение будет ограничено одномерным случаем, обобщение тривиально.Пусть траектория формируется поточечно (особенности многоуровневого алгоритмабудут описаны несколько позже). Задача, таким образом, состоит в правильной генера-ции координаты xn , если известны xn−1 и xn+1 . Пространственные граничные условия –периодические с периодом L, то есть −L/2 < x < L/2, x+L ≡ x. Учитывая периодическиеусловия по времени, двумерное пространство время имеет топологию тора.Ясно, что, в отличие от бесконечного пространства, координат x не достаточно для полного описания конфигурации, так как все множество траекторий, соединяющих xn1 и xn2включает пути с разным топологическим числом (количеством "намоток" на пространствовремя).
Для описания этих существенно различных (не переходящих друг в друга непрерывно) траекторий вводится новый набор координат {on } – топологических чисел научастках (xn−1 , xn ), принимающих, разумеется, целые значения.Решеточное действие выглядит теперь следующим образом:Nstep −1 S(x) =Xn=0m (xn+1 − xn + Lon+1 )2+ V (xn ) a.2a2(68)Так как использование кинетической части статистического веса в качестве пробного распределения весьма перспективно, надо получить аналог пригодный для практических численных рассчетов аналог формулы (40) для действия (68).
При этом, очевидно, должныгенерироваться значения как xn , так и on и on+1 , которые теперь являются полноправными переменными. Одновременное преобразование вместе с координатой обоих соседнихтопологических чисел введена потому, что при изменении одного из них изменение второго более вероятно (когда одна из точек переходит через границу, эффект усиливаетсяпри увеличении L). Выражениеn m oT (o0n , x0n , o0n+1 ) ∼ exp −(xn+1 − x0n + Lo0n+1 )2 + (x0n − xn−1 + Lo0n )22a~(69)не дает простого способа генерации.Вводя новые переменныеx0∞ = x0n + Lo0n ,(70)o0s = o0n+1 − o0n ,(71)и подставляя их в (69), получим:T (x0∞ , o0s )n m o0200 2.(x∞ − xn+1 − Los ) + (x∞ − xn−1 )∼ exp −2a~(72)Это выражение можно разложить, вводя условную вероятность:T (x0∞ , o0s ) = T (x0∞ |o0s )T (o0s ),(73)гдеT (o0s )Z=dx0∞ T (x0∞ , o0s ).(74)Теперь алгоритм понятен.
Сначала, в соответствии сon m2[xn+1 − xn−1 + Lo0s ]T (o0s ) ∼ exp −4a~(75)генерируется o0s , так как это – дискретная переменная, для этого требуется явное вычисление вероятностей каждого значения. По этой причине при численных рассчетах топологическое число приходится ограничивать, что при достаточно большом обрезании1mL2 2o≈ 5 − 104a~ max(76)не сказывается на их результатах в силу быстрого убывания функции распределения (75).Вероятность выбора x0∞ при уже известном o0s()0 2mx+x+Lon+1n−1sT (x0∞ |o0s ) ∼ exp −x0 −a~ ∞2(77)является гауссовой и позволяет использовать алгоритм (41). Наконец, с учетом условия−L/2 < x0n < L/2, формулы (70) и (71) позволяют однозначно определить x0n , o0n и o0n+1 .Важно отметить, что в многоуровневом алгоритме генерация топологического числадолжна производиться только на первом уровне.
В дальнейшем при делении отрезка пополам суммарное количество "намоток" должно сохраняться, то есть за o0s принимаетсязначение топологического числа на данном интервале, полученное на предыдущем уровне.Этот факт на может не радовать, так как применение формулы (75) заметно замедляетсчет. В результате рассмотренный алгоритм радикально повышает скорость генерациинескоррелированных конфигураций (см. (Рис. ??)).В заключение надо сказать, что при использовании периодических граничных условийпотенциал взаимодействия двух частиц должен быть записан в видеV(|x(1)n−x(2)n |)→∞XV(|x(1)no=−∞−x(2)n+ Lo|) ≈oVX(2)V (|x(1)n − xn + Lo|).(78)o=−oVРазумеется приближение бесконечного ряда конечной суммой при произвольной функ(1)(2)ции V (|xn − xn |) невозможно, но для для любого физически осмысленного потенциалавзаимодействия это можно сделать.1,00,90,8correlation0,70,6N0,5level00,4120,330,20,10,0-0,10200400600n8001000confРис.
3: Корреляция конфигураций в зависимости от количества уровней алгоритма, оптимальное (согласно (135)) действие уровня3.5Расчётная технология: распараллеливание, CUDA, работа спамятьюCUDA- архитектура некоторых видеокарт Nvidia, позволяющая выполнять вычисления навидеокарте. В процессе своего развития видеокарты становились способными выполнятьвсё более и более сложные операции по выводу и обработке изображения. Характернойоперацией является вращение сложной трёхмерной фигуры, составленной из большогочисла простых примитивов. Отсюда становится понятным, почему в видеокартах сейчаснаходится большое количество сравнительно слабых вычислительных ядер. В ходе дальнейшего развития стало возможным выполнять произвольные (а не только графические)вычисления на видеокартах (GPGPU-General Purpose GPU).
Для научных вычисленийважным этапом стало появление поддержки операций над числами двойной точности.8ln Lcorr64200246Рис. 4: Длина автокорреляции конфигураций в зависимости от количества уровней алгоритма, в логарифмической шкалеНа данный момент существуют 3 основных системы для выполнения расчётов на GPU:Nvidia CUDA, AMD FireStream и OpenCL. Первые 2 представляют собой совокупностьпрограммных и аппаратных средств и применимы только на видеокартах соответствующих производителей. OpenCL- язык программирования плюс API, ориентированный наразработку программ для гетерогенных платформ (т.
е. имеющих разнородные вычислительные модули в своём составе, вообще говоря не обязательно CPU+GPU). Главным егодостоинством является поддержка GPU как от Nvidia, так и от AMD. В качестве недостатка можно указать меньшую производительность [22].В силу разных причин на данный момент наибольшее распространение имеют системына основе решений Nvidia (для наглядности можно посмотреть рейтинг 500 самых мощныхкомпьютерных систем [?]).Архитектура видеокарт такова, что вычислительные ядра работают не независимодруг от друга, а группируются в пучки (warps), которые выполняют одинаковые операции над разными данными. Это так называемая SIMD (Single instruction, multiple data)архитектура построения арифметико-логических устройств. Если поведение программыразлично в зависимости от номера потока в пределах одного пучка, то сначала выполнятся одинаковые действия для одной группы потоков, потом для другой и так далее.То есть такие действия будут выполняться последовательно, а не параллельно, что плохоскажется на производительности.На логическом уровне есть 2 ступени распараллеливания при выполнении операцийна CUDA-совместимой видеокарте.
При вызове функции, исполняющейся на GPU, формируется grid, состоящий из указанного числа блоков, в каждом из которых находитсяопределённое число нитей-тредов, к-ые являются минимальной вычислительной единицей. Блоки исполняются независимо друг от друга, нити внутри блока группируются впучки, внутри которых исполняются синхронно. Адресацию блоков внутри grid-а можнозадавать одно- или двумерной, адресацию тредов внутри блока- трёх-, двух-, или одномерной. (см. рис. 5)Существует несколько типов памяти на GPU, которые различаются между собой размером, скоростью работы и возможностями доступа.