close

Вход

Забыли?

вход по аккаунту

?

Два статистических алгоритма для u(1) модели элементарных частиц на решетке.

код для вставкиСкачать
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
ДВА СТАТИСТИЧЕСКИХ АЛГОРИТМА ДЛЯ U(1) МОДЕЛИ
ЭЛЕМЕНТАРНЫХ ЧАСТИЦ НА РЕШЕТКЕ
Н.В. ЗВЕРЕВ, доц. каф. физики МГУЛ, канд. физ.-мат. наук
О
дной из важных задач математического моделирования и численных методов является
создание и реализация эффективных численных
методов и алгоритмов для исследования моделей
элементарных частиц. Для такого исследования
К. Вильсоном [2] предложен эффективный математический подход – метод решетки. В этом методе непрерывное пространство-время аппроксимируют дискретной совокупностью точек – узлов
решетки, а полевые и корреляционные функции
модели зависят от координат этих узлов.
Практический интерес представляет четырехмерная решеточная модель элементарных
частиц, построенная на калибровочной группе
матриц U(1) [2, 3]. При исследовании модели
частиц прежде всего необходимо найти ее многочисленные корреляционные функции (функции
Грина), которые затем используют для расчетов
физических характеристик частиц.
Исследования моделей частиц на решетке прямыми численными методами являются
нереальными из-за гигантского времени вычислений. К настоящему времени предложены два
статистических алгоритма для моделей частиц на
решетке с существенным сокращением времени
вычислений: метод гибридного Монте-Карло [4,
5] и двухшаговый мультибозонный алгоритм [6].
Данные алгоритмы отличаются многими входящими в их состав методами, процедурами и подходами.
Точность и продолжительность расчетов
и даже сама возможность вычислений с помощью любого алгоритма зависят от обоснованного
выбора его исходных параметров. А соотношения
для параметров указанных алгоритмов применительно к U(1) модели известны не были. Также
не были выполнены сравнительные исследования U(1) модели этими алгоритмами и не была
подтверждена расчетами пригодность каждого
алгоритма для данной модели. Рассмотрение и
решение этих математических проблем является
весьма актуальным.
Целью данной работы являются анализ
указанных двух статистических алгоритмов применительно к U(1) модели частиц на четырехмерной решетке, установление соотношений для
ЛЕСНОЙ ВЕСТНИК 2/2007
параметров и характеристик этих алгоритмов,
вычисление с помощью каждого алгоритма значений корреляционных функций в двух фазах U(1)
модели, выяснение пригодности рассматриваемых алгоритмов для такой модели по результатам
вычислений, а также сравнение продолжительности вычислений каждым алгоритмом.
1. U(1) модель частиц на решетке
1.1. Действие и поля
Рассматриваемая модель элементарных
частиц построена на матричной группе U(1), которая при определенных преобразованиях полей
обеспечивает инвариантность фундаментальной физической величины – действия модели.
Вследствие этой инвариантности выполняются
необходимые физические законы сохранения при
взаимодействиях и взаимопревращениях элементарных частиц. Действие S[U, ψ, ψ ] U(1) модели
на четырехмерной решетке пространства-времени
является суммой действия SG[U] калибровочного
поля, описывающего частицы – переносчики взаимодействия, и действия SF[U,ψ, ψ ] фермионных
полей, описывающих фермионные частицы [2]:
S[U,ψ, ψ ] = SG[U] + SF[U,ψ, ψ ].
Здесь
SG [U ] = β ∑ Re(1 − U x ,µν ) ,
x ,µ , ν
µ<ν
Nf
(1)
S F [U , ψ, ψ ] = ∑∑ ψ M [U ]x , y ψ ;
f =1 x , y
f
x
f
y
где Ux,µν – решеточный компактный тензор напряженностей калибровочного поля:
U x ,µν = U x ,µU x +µ^ ,νU x∗+ν^ ,µU x∗,ν ;
Ux,µ = exp(iAx,µ) – решеточное компактное калибровочное поле;
Ax,µ – вещественный потенциал калибровочного поля в интервале (–π, π];
β = 1/e02 – обратный квадрат затравочного заряда;
f
ψ x и ψ fx – фермионные поля, являющиеся антикоммутирующими переменными;
M [U] – комплексная фермионная матрица
M [U ]xy = δ xy −
{
}
−κ∑ (1 − γ µ ) U x ,µ δ x +µ€μ, y + (1 + γ µ )U y∗,µ δ y +µμ€, x ; (2)
µ
^
^
109
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
где κ = 1/(8+2m0) – хоппинг-параметр;
m0 – затравочная фермионная масса;
a – шаг решетки;
γµ – эрмитовы матрицы Дирака размером
4×4;
f – индекс поколения фермионов;
Nf – число поколений фермионов;
δxy – символ Кронекера;
x, y – узлы четырехмерной конечной решетки с целочисленными координатами
nµ = 0,1,…, Nµ −1;
µ,ν = 1,2,3,4 – направления решетки;
μ^ – единичный вектор в положительном направлении µ;
Nµ – число узлов решетки вдоль направления µ.
Полное число узлов, называемое объемом
решетки, равно V = N1N2N3N4. Здесь принята система естественных единиц измерения, в которой
квантовая постоянная  , скорость света c и шаг
решетки a выбраны равными  = c = a = 1.
Обычно рассматривают решетку с числами узлов N1 = N2 = N3 = Ns < N4. Объем такой решетки равен V = Ns3×N4.
Для полей на конечной решетке принимают одно из следующих граничных условий (со
знаком + или –)
f x + N νν€ = ± fx.
Для поля Ux,µ выбирают периодические
граничные условия (знак +) по всем направлениям
ν. Для полей ψfx и ψ xf принимают периодические
граничные условия (знак +) по пространственным направлениям ν = 1,2,3 и либо периодические, либо антипериодические (знак –) условия по
направлению времени ν = 4.
Матрица (2) удовлетворяет свойству
γ5-эрмитовости
M †[U] = γ5 M [U]γ5 ,
(3)
где γ5 = γ1γ2γ3γ4 – киральная матрица Дирака; эрмитово сопряжение † выполнено в пространстве матриц Дирака и в пространстве узлов решетки x.
ν^
1.2. Корреляционные функции и фазовая
структура модели
Для расчетов физических характеристик
частиц используют многочисленные корреляционные функции модели частиц (функции Грина).
Каждую корреляционную функцию 〈O〉 рассчитывают, главным образом, методом континуального
(функционального) интеграла. Этот метод заклю-
110
чается в усреднении соответствующей функции
O[U] по всем полям U, ψ и ψ , распределенным
с весом exp(−S[U,ψ, ψ ]). При этом после интегрирования только по фермионным полям ψ и ψ
получается следующая формула для корреляционной функции 〈O〉 модели на решетке [2]
1
N
O = ∫ O[U ]exp(− SG [U ]) det f M [U ] [dU ], (4)
Z
где Z – нормировочная постоянная.
Здесь в силу свойства (3) детерминант
фермионной матрицы (фермионный детерминант)
detM[U] является вещественной величиной.
Для упрощения расчетов значений 〈O〉 по
формуле (4) используют приближение статических фермионов, когда в диаграммах взаимодей­
ствий частиц пренебрегают вкладом фермионных
N
петель, т.е. полагают det M[U] = const. Для более
точного расчета значений 〈O〉 используют подход
динамических фермионов, когда в формуле (4)
N
учитывают зависимость det M[U] от U [7].
При исследовании моделей на решетке
обычно вычисляют следующие корреляционные
функции 〈O〉: среднюю «энергию» калибровочного поля 〈EG〉, скалярный конденсат 〈 ψ ψ〉 и «пионную» норму 〈Π〉. Эти средние величины вычисляют по формуле (4), в которой зависимость O[U]
для каждой из перечисленных функций определяется одной из следующих формул [1]
1
EG =
(5)
∑ Re(1 − U x ,µν ) ,
6V x ,µ ,ν
f
f
µ<ν
ψψ =
1
1
V
TrM -1[U ],
Tr γ 5 M -1[U ]γ 5 M -1[U ],
(6)
(7)
V
где Tr обозначает след в пространстве матриц Дирака и в пространстве узлов решетки x. В силу
свойства (3) зависимости (6) и (7) являются вещественными величинами, как и (5).
Частицы по U(1) модели на четырехмерной решетке имеют различные свойства в разных
областях значений обратного квадрата затравочного заряда β и хоппинг-параметра κ [3]. Такие
области называют фазами, из которых практический интерес представляют Кулоновская фаза и
фаза конфайнмента.
Для правильного описания свойств частиц
параметры β и κ следует выбирать вблизи линий
раздела фаз [1, 7]. На этих линиях решеточные
корреляционные функции имеют определенные
сингулярные свойства.
Π=
ЛЕСНОЙ ВЕСТНИК 2/2007
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
2. Метод гибридного Монте-Карло
2.1. Формула для корреляционных функций
Метод гибридного Монте-Карло предназначен для приближенных, но относительно быстрых вычислений корреляционных функций моделей на решетке путем статистического выбора
конечного числа конфигураций калибровочного
поля U и вспомогательных полей с весом в подходе динамических фермионов. При этом число
поколений фермионов Nf должно быть четным.
В данном методе вводят вспомогательные поля
P и χ таким образом, что формула (4) с учетом
свойства (3) для фермионной матрицы M[U] переходит в следующее соотношение [4, 5]
1
O = ∫ O[U ]exp( − H [U , P, χ]) [dU dP d χ ], (8)
Z
где H[U,P,χ] – гамильтониан, определенный формулой
1
H [U , P, χ] = SG [U ] + ∑ Px2,µ +
2 x ,µ
Nf /2
(
)
+ ∑ χ†f M † [U ]M [U ] χ f .
f =1
(9)
где Px,µ – вещественные компоненты вспомогательного поля P, называемого сопряженным импульсом;
χf – вспомогательные поля с 4V комплексными составляющими, называемые псевдофермионами.
Поля U, P и χ необходимо генерировать
случайным образом с весом exp(–H[U,P,χ]). Тогда формула (8) переходит в следующее соотношение, по которому выполняют численные расчеты
1 N
O = lim ∑ O[U ( α ) ],
(10)
N →∞ N
α=1
где O[U(α)] – значение функции O[U] при конфигурации U = U(α), состоящей из 4V комплексных составляющих поля U(α)x,µ ;
α – номер конфигурации поля: α = =1,2,…, N;
N – число генерируемых конфигураций поля
U(α).
При расчетах число N является конечным.
С увеличением N точность вычислений значений
〈O〉 возрастает.
2.2. Генерация полей
Каждую конфигурацию с номером α полей
U, P и χ с весом распределения, указанным к формуле (10), получают по следующей схеме [4, 5].
Сначала выбирают конфигурацию α сопряженного импульса P, генерируя случайным
образом компоненты Px,µ с весом exp(–P2x,µ /2),
ЛЕСНОЙ ВЕСТНИК 2/2007
а также выбирают конфигурацию α псевдофермионов χ по формуле χf = M†[U]ηf,
где ηf – вектор с 4V комплексными составляющими, генерируемый случайным образом с весом
exp(–η†fηf) (f = 1, …, Nf /2). При этом конфигурация α поля U известна.
Такую генерацию Px,µ и χf осуществляют
методом стохастического вектора. В данном методе генерируют случайные комплексные числа z
с весом exp(–C|z|2) по формуле
z = C −1 ln ξ −1 exp(iϕ) ,
где ξ и ϕ – вещественные случайные числа, равномерно распределенные в интервалах
[0,1] и (–π, π].
Каждую компоненту Px,µ сопряженного
импульса P полагают равной вещественной части
Re z при C = 1/2, а каждая из 4V комплексных составляющих вектора ηf равна z при C = 1.
Затем выполняют процедуру молекулярной динамики. Эта процедура состоит в получении промежуточной U′ и вспомогательной P '
конфигураций при известной конфигурации α
полей U, P и α по следующим формулам
∆τ
2)
Px(1/
= Px ,µ +
Fx ,µ [U ];
,µ
2
 U x( ,kµ) = U x( ,kµ−1) exp (i∆τ Px(,kµ−1/ 2) ),
 ( k +1/ 2)
= Px(,kµ−1/ 2) + ∆τ Fx ,µ U ( k )  ;
 Px ,µ
U x′,µ = U x( ,Nµ −1) exp (i∆τ Px(,µN −1/ 2) ),


∆τ
( N −1/ 2)
+ Fx ,µ [U ′].
 Px′,µ = Px ,µ

2
τ
τ
τ
Здесь ∆τ – шаг «времени», Nτ – число шагов «времени»; k = 1,… Nτ–1;
∂H
U(0)x,µ = Ux,µ; Fx ,µ [U ] = −
.
∂Ax ,µ
Рассматриваемая схема заканчивается
шагом Метрополиса для гамильтониана H. Этот
шаг заключается в выборе поля U для конфигурации α+1 из поля U′, найденного в процедуре
молекулярной динамики, и поля U конфигурации
α. Поле U′ принимают за новое (α+1) с вероятностью, равной
wacc[U ′,U] = min{1,exp(H[U,P,χ]−H[U ′,P ′,χ])}.
Для этого генерируют случайное число
z, равномерно распределенное в интервале [0,1].
Если z ≤ wacc[U′,U], то полученное поле U′ принимают за U в конфигурации α +1. В противном
случае за U в конфигурации α +1 принимают U
прежней (α) конфигурации.
111
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
В результате выполнения этой схемы для
всех N конфигураций U, P и χ эти поля оказываются распределенными с весом exp(–H[U,P,χ]),
при котором становится справедливым соотношение (10).
2.3. «Шахматное» разложение решетки
В процедуре молекулярной динамики и
в шаге Метрополиса для гамильтониана необходимо умножать псевдофермионы χf на обратные
матрицы M−1 и [M†]−1, где с учетом свойства (3)
матрицу [M†]−1 выражают через M−1. Это умножение выполняют путем решения систем уравнений вида Mξ = η с заданным η и неизвестным
ξ векторами, каждый из которых имеет 4V комплексные составляющие. Для ускорения решения
таких систем используют, согласно Т. де Гранду
[8], «шахматное» разложение конечной решетки
на четные и нечетные узлы с соответствующим
разложением величин M, ξ и η.
Для этого разложения на узлах решетки x
вводят функцию c(x), принимающую значения
+1 или –1 и удовлетворяющую соотношению
c(x± μ^ ) = −c(x) для всех узлов x и всех направлений µ. Такая функция с точностью до знака равна
∑ xν
c ( x ) = ( −1) ν .
Узел решетки x называют четным, если
сумма
∑ xν
ν
является четным числом и c(x) = 1, и нечетным,
если сумма
∑ xν
ν
нечетна и c(x) = −1.
При четном числе узлов решетки Nµ хотя
бы для одного направления µ общее количество
четных узлов равно количеству нечетных узлов.
Это позволяет разложить фермионную матрицу M
по (2) размера 4V × 4V на подматрицы одинакового размера 2V × 2V в подпространствах с четными
и нечетными узлами следующим образом [8]
 1e M eo 
M =
(11)
.
 M oe 1e 
Здесь 1e – единичная матрица в подпространствах с четными или нечетными узлами; Meo
и Moe – квадратные матрицы, первая из которых
отображает подпространство с нечетными узлами решетки, обозначаемое индексом o, в подпространство с четными узлами, обозначаемое как e,
а вторая отображает подпространство e в подпро-
112
странство o. Далее, разлагают векторы η и ξ с 4V
компонентами на векторы с 2V компонентами в
подпространствах с четными и нечетными узлами
 ηe 
 ξe 
η =   , ξ =  .
 ηo 
 ξo 
Используя разложение этих векторов, а также разложение (11) фермионной матрицы M, решают следующую систему из 4V уравнений для неизвестных векторов ξe и ξo, эквивалентную Mξ = η
Qξo = ηo − M oe ηe ,
(12)

 ξe = ηe − M eo ξo .
Здесь матрица Q размера 2V×2V, действующая в подпространстве с нечетными узлами, определена по формуле
Q = 1e − MoeMeo.
(13)
В силу «шахматного» разложения (11)
матрицы M справедливо тождество
detQ = detM.
(14)
Матрица Q обладает также свойством γ5эрмитовости, которое следует из (3)
Q† = γ5Qγ5.
(15)
Из системы (12) с учетом (13) видно, что
для нахождения вектора ξ достаточно решить линейную систему из 2V уравнений, верхнюю в (12)
с неизвестным вектором ξo, в подпространстве с
нечетными узлами. Вектор ξe элементарно выражается через ξo по системе из 2V уравнений, нижней в (12).
2.4. Метод сопряженных градиентов
Нахождение вектора ξo в верхней системе
уравнений (12) осуществляют методом сопряженных градиентов [7, 9]. Для использования этого
метода верхнюю часть системы (12) путем ее умножения на эрмитово сопряженную матрицу Q†
приводят к следующей эквивалентной системе
Q†Qξo = Q†(ηo − Moeηe).
Данная система имеет структуру в виде
Ax = y,
где A – эрмитова положительно-определенная матрица (в данном случае размера 2V×2V);
y – заданный вектор (с 2V компонентами);
x – неизвестный вектор (с 2V компонентами).
Метод сопряженных градиентов представляет собой следующую итеративную схему
построения векторов xn, сходящихся к x
2
rn
xn+1 = xn + αngn, rn = y – Axn, α n = †
,
g n Ag n
2
r
rn+1 = rn – αnAgn, β n = n +1 2 , gn+1 = rn+1 + βngn. (16)
rn
ЛЕСНОЙ ВЕСТНИК 2/2007
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
Здесь обозначено
Здесь n = 0,1,…; a = a † a – норма
1
A 12 = ∫ A[U ]exp( − SG [U ] −
вектора a; rn – векторы, называемые остатками;
(20)
Z
gn – векторы, называемые градиентами. Началь−1
†
− S B [Φ, U ]) det P2 Q [U ]Q[U ] [dU d Φ ];
ные векторы принимают равными g0 = r0 = y − Ax0,
где SG[U] – действие калибровочного поля U по (1);
где вектор x0 генерируют методом стохастическоQ – матрица размером 2V×2V, определенная
го вектора с весом exp(−||x0||2). Векторы-остатки rn
согласно (13);
являются попарно ортогональными.
SB[Φ,U] – мультибозонное действие
Выполнение схемы (16) прекращают при
n
^
^
достижении условия
S B [Φ, U ] = ∑ Φ †j M †j [U ] M j [U ]Φ j ;
j =1
||rn+1|| < δ,
(17)
Φ
–
вспомогательные
поля с 4V комплексныj
где δ – заданная малая величина: δ << 1.
ми
составляющими,
называемые мультиЭта величина для остановки метода собозонами;
пряженных градиентов в процедуре молекуляр-n
^
^
†
†
матрицы размером
S B [Φ, U ] это= ∑ Φ j M j [U ] M j [[U]
U ]Φ– jкомплексные
;
ной динамики равна δ = δmd, а для остановки
4V×4V, определяемые по формуле
j =1
го метода в шаге Метрополиса для гамильтониана
n
^
^
M eo 
 1e
†
δ = δacc. Согласно Р. Гупте и др. [10], необходимо
S B [Φ, U ] = ∑ Φ j M †j [U ] M
M€ jj [U
= ]Φ j ;
;
выполнение условия δacc << δmd.
j =1
 γ 5 M oe ( γ 5 − ρ j )1e 
Для векторов-остатков rn справедливо слематрицы 1e, Meo и Moe определены в (11); ρj – компдующее неравенство
лексные корни вещественного многочлена P1(x2):
n
ζ −1
λ max
2ρ
P1(ρj2) = 0, удовлетворяющие при j = 1,…,n1/2 усrn ≤
r
,
ρ
=
,
ζ
=
,
0
ловиям в виде ρ j + n /2 = −ρ∗j и Imρj > 0; P1(x), P2(x) и
1 + ρ2 n
λ min
ζ +1
где λmax и λmin – максимальное и минимальное собсP4(x) – вещественные положительные многочлетвенные значения положительно-опредены степеней n1, n2 и n4, которые используют для
ленной матрицы A (в данном случае Q†Q);
грубой, промежуточной и точной аппроксимаций
−N /2
ζ – число обусловленности матрицы A.
функции x
на интервале x∈[ε, λ]
−N /2
−N /2
Из указанного неравенства для ||r n|| в
P1 ( x ) ≈ x
, P1 ( x ) P2 ( x ) ~
−− Nx / 2 ,
случае ζ >> 1 получим следующую оценку для
P1 ( x) P2 ( x) P4 ( x) = x
(21)
числа итераций N(CG) схемы метода сопряженИнтервал [ε, λ] содержит минимальное
ных градиентов (16) до достижения условия ос〈λmin〉 и максимальное 〈λmax〉 средние собствентановки (17)
ные значения матрицы Q†Q при их усреднении
по всем конфигурациям калибровочного поля U с
ζ 2
N
N ( CG ) ~
ln .
(18)
−
весом exp(−SG[U]) det M[U].
2
δ
Поля U и Φ необходимо генерировать случайным
образом с весом
3. Двухшаговый мультибозонный алгоритм
exp(−SG[U]−SB[Φ,U])det−1P2(Q†[U]Q[U]).
3.1. Формула для корреляционных функций
Тогда формула (20) переходит в следуюДвухшаговый мультибозонный алгоритм
щее соотношение, по которому выполняются чисявляется статистическим методом приближенных
ленные расчеты
и относительно быстрых вычислений корреляци1 N
A
=
lim
(22)
∑ A[U ( α ) ],
онных функций моделей на решетке в подходе
12
N →∞ N
α=1
динамических фермионов при произвольном знагде U(α), α и N – те же, что и в формуле (10).
чении числа поколений фермионов Nf. В данный
3.2. Генерация полей
алгоритм введены вспомогательные поля – мультибозоны, мультибозонное действие и аппроксиКаждую конфигурацию с номером α помирующие многочлены. При этом формула (4)
лей U и Φ с весом распределения, указанным к
преобразуется с учетом свойств (14) и (15) к слеформуле (22), получают путем выполнения двух
дующему виду [6]
шагов.
(
)
1
1
1
1
f
f
f
f
f
O =
(
)
Nf
O[U ]det −1 P Q † [U ]Q[U ] signdet Q[U ]
(
)
Nf
det −1 P Q † [U ]Q[U ] signdet Q[U ]
ЛЕСНОЙ ВЕСТНИК 2/2007
12
12
. (19)
Для этого вводят матрицы Π+ и Π− – проекторы на четное и нечетное подпространства узлов решетки при ее «шахматном» разложении
113
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
1
Π+ =  e
 0e
0e 
 0e
, Π− = 

0e 
 0e
0e 
,
1e 
где 0e – нулевая матрица в подпространствах с
четными или нечетными узлами.
Эти проекторы имеют следующие элементы в пространстве узлов решетки
∑x
1 ± ( −1)
δ xy .
(Π ± )xy =
2
С помощью данных проекторов матрицы
n
^
^
†
†
принимают
вид
= ∑ Φ j M j [U ] M nj [U ]Φ ^j ;
^
j =1
†
†
= ]Φ
Πj+;M + Π−γ5M −ρjΠ− .
S B [Φ, U ] = ∑ Φ j M j [U ] M j [U
j =1
Затем действие SB[Φ,U] представляют в
следующем виде:
2
V j , x ,r
S B [Φ , U ] = A j , x , r Φ j , x , r +
+ SB0 ,
Aj , x , r
где SB0 – слагаемое, не зависящее от величины
Φj,x,r;
µ
µ
1
1
V j , x , r = κ∑ {U ∗
^
x −µ , µ
µ,s
([ρ j , x γ 5 (1 + γ µ ) − 2]rs Φ
+[ γ 5 (1 + γ µ )]rs χ
−γ µ ) − 2]rs Φ
+
1 − ( −1)
χ j , x ,µ , r
^
j , x −µ ,µ , s
^
j , x +µ ,µ , s
Aj , x , r = 1 + 4 × 4 κ 2 +
∑x
)} ,
µ
{| ρ
µ
j
ν,s
ν≠µ
+[ γ 5 (1 + γ ν )]rs U ∗
^
x −ν , ν
1 + ( −1)
}
|2 − 2[ γ 5 ]rr Re ρ j ф ,
2
= κ∑ {[ γ (1 − γ )] U Φ
5
ν rs
x,ν
ρ j,x =
+
) + U x , µ ([ρ j , x γ 5 (1 −
+ [ γ 5 (1 − γ µ )]rs χ
^
j , x +µ , s
^
j , x −µ , s
∑x
µ
µ
ρj +
Φ
^
^
j , x+ ,s
+
ν
},
j , x− ,s
ν
1 − ( −1)
∑x
µ
µ
ρ∗j .
2
2
Здесь j = 1,…,n1; x – узел решетки; r,s
– номера строк и столбцов матриц Дирака:
r,s = 1,2,3,4; матрица γ5 диагональна по r и s, и все
ее ненулевые элементы [γ5]rr = ±1.
Далее сумму действий SG[U] + +SB[Φ,U]
представляют в виде
SG[U] + SB[Φ,U] = −Re(Ux,µFx,µ) + SGB0 ,
где SGB0 – слагаемое, не зависящее от Ux,µ;
Fx ,µ = β∑ (U
ν
ν≠µ
^
U∗
x +µ ,ν
^
U x∗,ν + U ∗
x + ν ,µ
^
^
U∗
x − + ,ν
ν µ
^
U
x − ,µ
ν
^
x − ,ν
)−
ν
∗
−2 κ ∑ {[ γ 5 (1 − γ µ )]rs × (ρ j , x Φ j , x , r Φ j , x + ^ , s +
∗
+Φ j , x ,r χ
114
j ,r ,s
^
j , x + ,µ , s
µ
+ χ∗j , x ,µ , r Φ
µ
^
j , x+ ,s
µ
∗
) − 2δ rs Φ j , x , r Φ
^
j , x+ ,s
µ
}.
В первом шаге каждую новую (α+1) конфигурацию поля Φ и промежуточную конфигурацию U′ получают из предыдущих (α) конфигураций U и Φ путем выполнения следующих четырех
последовательных процедур методами «тепловой
бани» (HB) и верхней релаксации (OR):
1. Метод HB для мультибозонов с получением 4Vn1 составляющих поля Φ′ из Φ для каждой совокупности индексов (j,x,r) по формулам
V j , x ,r
Φ′j , x , r = −
+ ξ j , x ,r ,
Aj , x , r
где Vj,x,r = V[Φ,U];
ξj,x,r – комплексное случайное число, распределенное с весом exp(−Aj,x,r|ξj,x,r|2) и генерируемое методом стохастического вектора (п. 2.2).
Найденное поле Φ′ принимают за Φ.
2. Метод OR для мультибозонов с получением 4Vn1 составляющих поля Φ′ из найденного
в процедуре 1 поля Φ для каждой совокупности
(j,x,r) по формулам
V j , x ,r
Φ′j , x , r = −Φ j , x , r − 2
.
Aj , x , r
Найденное поле Φ′ принимают за поле Φ
в α+1 конфигурации.
3. Метод HB для калибровочного поля с
получением 4V составляющих поля U′ из U для
каждой совокупности узлов и направлений решетки (x,µ) при поле Φ в α+1 конфигурации по
формулам
U′x,µ = (Fx,µ / |Fx,µ|)−1exp(iηx,µ),
где Fx,µ = F[Φ,U];
ηx,µ – вещественное случайное число в интервале (−π,π], распределенное с весом
exp(|Fx,µ|cos ηx,µ).
Генерацию таких чисел осуществляют с использованием вспомогательного веса и метода Метрополиса [7]. Найденное поле U′ принимают за U.
4. Метод OR для калибровочного поля с
получением 4V составляющих поля U′ из поля U,
найденного в процедуре 3, для каждой совокупности (x,µ) при поле Φ в α+1 конфигурации по
формулам
U′x,µ = U∗x,µ (Fx,µ / |Fx,µ|)−2.
Далее выполняют второй шаг алгоритма
– метод Метрополиса для детерминантов [6]. Этот
шаг заключается в выборе поля U для конфигурации α+1 из поля U′, найденного в процедуре 4 первого шага, и поля U конфигурации α. Поле U′ принимают за новое (α+1) с вероятностью, равной
ЛЕСНОЙ ВЕСТНИК 2/2007
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
 det −1 P2 (Q † [U ′]Q[U ′])
wacc [U ′, U ] = min  1,
,
 det −1 P2 (Q † [U ]Q[U ])  (23)


где, согласно И. Монтвею [6], отношение детерминантов определяют по формуле
det −1 P2 Q † [U ′]Q[U ′]
=
det −1 P2 Q † [U ]Q[U ]
(
(
)
)
= exp( −ξ† P2 (Q † [U ′]Q[U ′])ξ + η† η) . (24)
η
Здесь обозначено
1
B η = ∫ B[η]exp( −η† η) [d η];
(25)
Z
где ξ – вектор с 2V комплексными составляющими;
ξ = P3(Q†[U]Q[U])η;
P3(x) – вещественный положительный многочлен степени n3, аппроксимирующий с
высокой точностью функцию P2−1/2(x) на
интервале x∈[ε,λ]
P3(x) = P2−1/2(x).
(26)
При расчетах вектор η генерируют случайным образом с весом exp(−η†η) методом стохастического вектора (см. п. 3.2), затем вычисляют вектор ξ и отношение детерминантов по (24).
Далее генерируют случайное число z, равномерно распределенное в интервале [0,1]. Полученное
поле U′ принимают за U в конфигурации α, если
z ≤ wacc[U′,U]. В противном случае за U в конфигурации α+1 принимают U прежней (α) конфигурации.
В результате выполнения первого и второго шагов алгоритма для всех N конфигураций
U и Φ эти поля оказываются распределенными с
требуемым весом
exp(−SG[U]−SB[Φ,U])det−1P2(Q†[U]Q[U]).
При вычислениях величин 〈O〉 по (19) и
(22) значения det−1P4(Q†Q) находят, согласно И.
Монтвею [6], по формуле
det−1P4(Q†Q) = 〈exp(−η†P4(Q†Q)η + η†η)〉η, (27)
где усреднение выполнено по (25).
3.3. Аппроксимирующие многочлены
Апппроксимирующие многочлены P1(x) с
заданными степенями n1, где k = =1,2,3,4, вычисляют, согласно И. Монтвею [11], по следующей
схеме.
Выбирают интервал аппроксимаций [ε, λ]
так, чтобы он содержал средние минимальное
〈λmin〉 и максимальное 〈λmax〉 собственные значения
матрицы Q†Q: [〈λmin〉,〈λmax〉] ⊂ [ε,λ]. Эти значения
〈λmin〉 и 〈λmax〉 вычисляют по формулам (19) и (22),
ЛЕСНОЙ ВЕСТНИК 2/2007
где в (19) величина O заменяется на λmin или λmax.
Далее, интегральным методом наименьших квадратов находят коэффициенты аппроксимирующих многочленов P1(x), P2(x) и P4(x) при
заданных степенях n1, n2 и n4
Nf /2
N /2
1 P
min,
x x P1 ( P
x1)(−x)1−
→→
min,
f
1
Nf /2
P1
N /2
1 P
min,
x x P1 ( P
x1)(Px2)(Px2)(−x)1−
→→
min,
ф1 иксирован
P
фиксирован
P
, ,
Р фиксирован,
f
1
2
P2
1
Nf /2
N /2
1 P
min,
x x P1 ( P
x1)(Px2)(Px2)(Px4)(Px4)(−x)1−
→→
min,
PP1 , P2 Р , Р фиксированы,
P1 , фиксированы
фиксированы
, ,
(28)
2
1
2
f
P
где ||f(x)|| – интегральная норма функции f(x):
1/ 2
 1 λ

f ( x) = 
dx | f ( x) |2  .
∫
λ−ε ε

(29)
Минимизацию каждого интеграла (29) в
(28) выполняют по коэффициентам многочлена
P1(x), P2(x) и P4(x) соответственно.
Нахождение коэффициентов многочлена
P3(x) заданной степени n3 осуществляют итерационным методом, аналогичным методу Ньютона
касательных
P3( j +1) ( x ) =
1
2
(P
( j)
3
)
( x) + P3( j ) ( x) , j = 0,1,….
Здесь P3(j)(x) и P3( j ) (x) – многочлены степени n3 для всех j. Многочлен P3( j ) (x) аппроксимирует функцию (P2(x) P3(j)(x))−1. Его коэффициенты
находят интегральным методом наименьших
квадратов
( j)
P
3
P ( x ) P ( j ) ( x ) P ( j ) ( x ) − 1 →
min,
2
3
3
(j)
(30)
P2 , P Р2, Р3 фиксированы,
фиксированы
.
Минимизацию интеграла (29) для (30) выполняют по коэффициентам многочлена P3( j ) (x).
В пределе j → ∞ многочлен P3(j)(x) → P3(x).
После вычисления коэффициентов многочленов находят методом Лежандра корни ρj
многочлена P1(x2), необходимые для первого шага
данного алгоритма.
Вычисление матриц Pk(Q†Q), где k = 2,3,4,
осуществляют, согласно И. Монтвею [11], по следующим формулам, в которых x заменяется на Q†Q
( j)
3
nk
Pk ( x ) = ∑ c (j k ) O (j k ) ( x ) , k = 2, 3, 4;
j =0
(k )
j +1
O ( x ) = ( x + β(jk ) )O (j k ) ( x ) + γ (jk−)1O (j −k 1) ( x).
Здесь Oj(k)(x) – многочлены степени j,
удовлетворяющие условию ортогональности
λ
∫ dx ρ
(k )
( x) Ol( k ) ( x) O (j k ) ( x) = 0 , l ≠ j ,
ε
115
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
где
(k )
0
( 2)
(k )
−1
2
1
2
3
O ( x ) = 1, γ
Nf
= 0,
ρ ( x) = x P ( x) ,
ρ(3) ( x ) = P22 ( x) P ( x) ,
Nf
ρ( 4) ( x ) = x P12 ( x) P22 ( x).
3.4. Явно повторяемый метод Ланцоша
Для выбора границ интервала [ε, λ] необходимо вычислять для каждой конфигурации
поля U минимальное λmin и максимальное λmax
собственные значения матрицы Q†Q. Эти собственные значения находят явно повторяемым методом Ланцоша [12].
Этот метод в общем случае эрмитовой
матрицы A представляет собой итерации по следующей схеме:
1. При первой итерации n = 1 принимают
единичный вектор e0(n) произвольного направления. Задают целое число M ≥ 2.
2. Для данного номера итераций n вычисляют M единичных векторов ej(n) (j = 0,1,…,M−1)
по формуле
Ae (j n ) = β(jn ) e (j n+1) + α (jn ) e (j n ) + β(jn−)1e (j n−1) ,
где
αj(n) = (ej(n))†Aej(n),
(n)
(n)
βj = ||Aej −αj(n)ej(n)− −βj−1(n)ej−1(n) ||, e−1(n) = 0.
3. Строят эрмитову трехдиагональную матрицу T(n) размером M × M с элементами αj(n) (j = 0,1,…,M−1) на главной диагонали и
βj(n) (j = 0,…,M−2) на соседних диагоналях.
4. Вычисляют методом Ритца либо минимальное, либо максимальное собственное значение
λ(n) матрицы T(n) и соответствующий единичный
собственный вектор-столбец s(n): T(n)s(n) = λ(n)s(n).
5. Проверяют выполнение следующего
условия для номера итераций n>1:
|λ(n) −λ(n−1)| ≤ δ,
(31)
где δ – заданное малое число.
Если это условие выполнено, то вычисления останавливают и принимают полученное
λ(n) либо за λmin, либо за λmax собственное значение
матрицы A. В противном случае переходят к следующему шагу.
6. Строят новый исходный вектор e0(n+1)
для следующей (n+1) итерации:
M −1
e0( n +1) = ∑ s (j n ) e (j n ) ,
j =0
где sj(n) – компоненты вектора-столбца s(n).
7. Возвращаются к шагу 2 при n → n+1 и
e0(n+1) → e0(n).
116
В случае данного алгоритма матрица
A = Q Q.
†
3.5. Вероятность принятия новых
калибровочных полей
Средняя вероятность 〈wacc〉 принятия новых калибровочных полей U в двухшаговом мультибозонном алгоритме равна отношению числа
полей U, принятых за новые, к полному числу генерируемых полей. Применяя методы А. Борелли
и М. Кройца, с учетом выражений (23) и (24), получим при 〈(∆E)2〉 <
~ 1 следующее соотношение
1
wacc = 1 −
( ∆E ) 2 ,
(32)
2π
где
〈(∆E)2〉 ∼ V||P2(x) −1||2,
(33)
где ||f(x)|| – интегральная норма по (29);
(
)
ξ = P (Q [U ]Q[U ])η ;
∆E = ξ† P2 Q † [U ′]Q[U ′] ξ − η† η ;
†
3
1
C[U ′, Φ′, U , Φ, η] ×
Z∫
× exp − SG [U ] − S B [Φ, U ] − η† η ×
× p12 [U ′, Φ′, U , Φ ] [dU d Φ dU ′ d Φ′ d η];
C =
(
)
p12[U′,Φ′,U,Φ] – плотность вероятности перехода от полей U и Φ конфигурации α к
полям U′ и Φ′ конфигурации α+1.
4. Исходные параметры и производительность
двух алгоритмов для U(1) модели
Для исследования U(1) модели каждым из
двух описанных алгоритмов и их сравнения необходимо обоснованно выбрать значения исходных
параметров этих алгоритмов. В противном случае время вычислений существенно увеличится
и результаты расчетов могут быть неправильными. Ниже получим соотношения для значений
исходных параметров и для производительности
рассматриваемых алгоритмов при одинаковом заданном числе поколений фермионов Nf = 2.
4.1. Параметры метода гибридного
Монте-Карло
Для расчетов методом гибридного Монте-Карло должны быть заданы следующие его
исходные параметры: шаг «времени» ∆τ, число
шагов «времени» Nτ, малые параметры в условии
(17) остановки метода сопряженных градиентов
δmd в молекулярной динамике и δacc в шаге Метрополиса для гамильтониана.
Для получения необходимых соотношений рассмотрим положительно-определенную
ЛЕСНОЙ ВЕСТНИК 2/2007
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
матрицу Ω размером 4V×4V, задаваемую формулой
 ∂ 2 H  †  ∂ 2 H
Ω =  2   2
 ∂A   ∂A


 
1/ 2
,
где H – гамильтониан по (9);
A – векторный потенциал поля U; усреднение выполнено по полям U, P и χ с весом
exp(−H[U,P,χ]).
Оценим норму этой матрицы ||Ω|| при
Nf = 2. Для этого используем свойства симметрии
U(1) модели относительно дискретных сдвигов и
поворотов, формулу для средних значений псевдофермионов в виде 〈χ†χ〉 = M†M и формулу для
собственных значений λj матрицы M†M в виде λj
(M†M) = λj2[(M†M)1/2]. Кроме того, сумму по собственным значениям λj заменим интегрированием,
а также учтем, что для матриц M†M и Q†Q справедливы, согласно де Гранду [14], оценки в виде
λ max ( M † M )
(34)
Здесь ζ – число обусловленности матрицы
Q†Q. В результате получим следующую оценку:
†
Ω ~ Tr ( M M )
−1
λ
1
∑λ
~
j
~
2
j
max
λ
∫
min
dλ
λ2
~ ζ , (35)
где λj, λmin и λmax – собственные значения матрицы
(M†M)1/2.
При исследовании модели свободного
бозонного поля на решетке методом гибридного Монте-Карло А. Кеннеди и Б. Пендлетон [13]
получили следующие соотношения, обеспечивающие высокие точность и производительность
вычислений при средней вероятности принятия
новых полей 〈wacc〉 ∼ 1
V Ω
где
2
(∆τ )
4
< 1,
~
Ω=
Ω N τ ∆τ ~ 1,
(36)
∂2 H
∂A2
– матрица, не зависящая от свободного бозонного поля A;
H – гамильтониан бозонного поля;
V – объем решетки.
Мы применим к U(1) модели соотношения (36), в которых величину ||Ω|| возьмем по (35).
Тогда получим [14]
1
∆τ <
, N > V 1/ 4 ,
(37)
~ (V ζ )1/ 4 τ ~
ЛЕСНОЙ ВЕСТНИК 2/2007
N ( CG ) ~
ζ ln V .
(39)
4.2. Параметры двухшагового
мультибозонного алгоритма
Для расчетов двухшаговым мультибозонным алгоритмом должны быть заданы следующие
его исходные параметры: границы ε и λ интервала
−N /2
аппроксимаций функции x
по (21) и степени
n1, n2, n3 и n4 аппроксимирующих многочленов.
Границы ε и λ интервала аппроксимаций
выберем следующим образом [14]
ε = 0,5〈λmin〉, λ = (1.2 – 1.4)〈λmax〉,
(40)
где 〈λmin〉 и 〈λmax〉 – средние минимальное и максимальное собственные значения матрицы
Q†Q, усредненные, как в (19).
Такой выбор учитывает необходимость
нахождения крайних значений 〈λmin〉 − σ λ
и
min
〈λmax〉 + σ λ
внутри интервала [ε, λ], где σ λ и
max
min
σλ
– среднеквадратичные отклонения.
max
Согласно И. Монтвею [6], высокие точность и производительность данного алгоритма
обеспечиваются при условиях
〈wacc〉 ∼ 1, 〈det−1P4(Q†Q)〉 ∼ 1,
где усреднение вероятности wacc выполнено так
же, как в (32), а усреднение величины det−1P4(Q†Q)
– как в (19). При этих условиях из (27), (32) и (33)
получим следующие оценки для многочленов
P2(x) и P4(x)
1
1
P2 ( x ) − 1 ~
, P4 ( x ) − 1 ~ ,
(41)
V
V
где норма ||f(x)|| определена по (29).
Точные аппроксимации (21) и (26) обеспечиваются, согласно И. Монтвею [6], при условиях
в виде
1
P2 ( x ) P32 ( x ) − 1 ~ 2 ,
V
1
N /2
x P1 ( x) P2 ( x) P4 ( x) − 1 ~ 2 .
(42)
V
f
λ max (Q †Q )
~ζ=
,
λ min ( M † M )
λ min (Q †Q )
λ max ( M † M ) ~ 1.
Для значений параметров δmd и δacc принимаем, согласно Р. Гупте и др. [10], следующие
оценки
1
1
δ md ~ , δ acc ~ 2 ,
(38)
V
V
Итак, получены соотношения (37) и (38)
для исходных параметров метода гибридного
Монте-Карло применительно к U(1) модели при
Nf = 2.
Из соотношений (18) и (38) имеем следующую оценку для значения среднего числа итераций 〈N(CG)〉 в методе сопряженных градиентов
f
117
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
Используя (41) и (42), получим следующие оценки
1
N /2
x P1 ( x) − 1 ~
,
V
1
N /2
x P1 ( x) P2 ( x) − 1 ~ .
(43)
V
Для оценки значений степеней nk многочленов Pk(x) введем, удобные для анализа многочлены Pk (x) степеней nk , составленные из
многочленов Чебышева Tn +1 (x), где k = 1,2,3,4.
k
Многочлены Pk (x) при k = 1,2,4 аппроксимиру−N /2
ют функцию x
на отрезке [ε′,λ′] с одинаковым небольшим отклонением от 1 функций
N /2
x Pk (x), а многочлен P3 (x) аппроксимирует
функцию P2(x)P32(x)
−N /2
P ( x ) ≈ P ( x) ≈ x
,
f
f
f
f
f
1
1
−N /2
P2 ( x ) ≈ P1 ( x) P2 ( x) ~
x
,
−
2

P3 ( x ) ≈ P2 ( x) P3 ( x) = 1,
−N
P ( x ) ≈ P ( x) P ( x) P ( x) = x
f
4
1
Здесь
2
(44)
f
/2
4
1  Tn +1( λ′−ε′ )
Pk ( x ) =  1 − k
 ,
λ ′ +ε ′
T
−
x 
n +1 ( λ ′ −ε ′ ) 
k
(45)
(x)=cos[( nk +1)arccos x].
Согласно нашим расчетам при Nf = 2, многочлены P1 (x) и P1(x) практически одинаково ап−N /2
проксимируют функцию x
при выполнении
следующих условий
′
n1 = n1 , ε′ ~
(46)
− ελ , λ = λ .
Для многочленов Чебышева справедлива
оценка в виде
k
+1
f
Tn
k
+1
2( n +1)
 λ ′ + ε′ 
k
− ′ ′ ~ e
 λ −ε 
ε′ / λ ′
.
(47)
Из определения числа обусловленности ζ
по (34) с учетом (40) следует оценка в виде
λ
~ ζ.
(48)
ε
Из (43), первого соотношения (44) и (45)
– (48) получим
−1
−1/ 4
1
−2( n +1) ζ
N /2
 λ ′ + ε′ 
1

x P1 ( x) − 1 ~ Tn +1  −
~
e
~
,

 λ ′ − ε′ 
V
откуда имеем оценку степени n1 [14]
n1 ∼ ζ1/4lnV.
(49)
Расчеты при Nf = 2 показали, что многочлены P2 (x) и P1(x)P2(x) практически одинаково
−N /2
аппроксимируют функцию x
при условиях
′ = λ.
n2 = n1 + n2 , ε′ ~
ε
,
λ
(50)
−
Из (43), второго соотношения (44), а также (45), (47), (48) и (50) получим
f
k
f
118
Nf /2
 λ ′ + ε′ 
P2 ( x) − 1 ~ Tn +1  −

 λ ′ − ε′ 
2
−1
~ e −2( n + n +1) ζ
1
2
−1/ 2
~
1
V
,
откуда имеем оценку для степени n2 [14]
n2 ∼ ζ1/2lnV.
(51)
Результаты аналогичных расчетов значений многочленов P3 (x), P (x), P2(x)P32(x) и
P1(x)P2(x)P4(x) с выполнением соотношений (44)
и последующий анализ приводят к следующим
условиям для степеней n3 и n4 [14]
n3 > n2 , n4 > n2 .
(52)
~
~
Итак, получены соотношения (40, 49, 51,
52) для исходных параметров двухшагового мультибозонного алгоритма применительно к U(1)
модели при Nf = 2. В расчетах корреляционных
функций 〈O〉 используются только аппроксимирующие многочлены Pk(x), степени которых nk выбирают по (49, 51, 52).
4.3. Производительность алгоритмов
.
2 x − ( λ ′ +ε ′ )
где Tn
x
Производительностью алгоритма вычислений конкретной корреляционной функции 〈O〉
модели на решетке назовем величину P [14], равную
P = (〈Nop〉〈τint〉)−1,
(53)
где 〈Nop〉 – среднее число компьютерных операций, необходимых для перехода от одной
конфигурации полей к последующей;
〈τint〉 – интегральное время автокорреляций,
равное среднему расстоянию между ближайшими статистически независимыми
выборочными значениями O[U(α)] этой
корреляционной функции [7].
По способу суммирования [9] это время
вычисляют по формуле
1 1 W Γ (α )
τint = + ∑ O
,
(54)
2 W α=1 Γ O (0)
где W – некоторое обрезающее число вычислений: 1 << W << N ;
ΓO(α) – автокорреляционная функция
Γ O (α ) =
1
N −α

∑ O
N −α 
j =1
j
−
1
N −α
∑O
N −α
k =1
k
N
1


Ok  .
∑
 O j +α −
N
−
α


k =α+1
Здесь Oα = O[U(α)]; остальные обозначения те же, что и в формуле (10).
Время вычисления величины 〈O〉 определяется формулой
tN
t = tc N op N = c i ,
P
где Ni – число статистически независимых
значений Oα: Ni = N / 〈τint〉;
tc – время выполнения одной компьютерной
операции.
ЛЕСНОЙ ВЕСТНИК 2/2007
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
Лучшим алгоритмом является такой, у которого большая производительность P при одинаковых параметрах модели Nf, V, β, κ и числе N
генерируемых конфигураций полей.
В случае метода гибридного Монте-Карло
(HMC) на модели свободного бозонного поля на
решетке А. Кеннеди и Б. Пендлетон [13] получили следующую оценку времени автокорреляций:
〈τintHMC〉 ∼ (Nτ ∆τ)−2.
(55)
Мы так же, как и в п. 4.1, применим это
соотношение к U(1) модели, в которое подставим
наши оценки (37).
В случае двухшагового мультибозонного алгоритма (TSMB), согласно М. Пиардону, К.
Янсену и А. Бориси, для времени автокорреляций
справедливо соотношение




−1
τint ~  min Im ρ 2j  ,
j
(56)
где ρj2 – j-й корень аппроксимирующего многочлена P1(x): P1(ρj2) = 0.
Для рассматриваемой U(1) модели при
Nf = 2 заменим, согласно (44), многочлен P1(x) на
P1 (x) и с учетом (45) найдем корни ρj2 многочлена
P1 (x), при которых P1 (ρj2) =0
πj
2πj
ρ 2j = (λ′ + ε′) sin 2
+ i λ′ε′ sin
,
n1 + 1
n1 + 1
j = 1, 2, … , n1 .
Подставляя это выражение в (56) и учитывая (40, 46, 48 и 34), получим следующую оценку
−1
1/ 2

n1  λ′ 
2π 
TSMB
′
′
τint
~  λ ε sin
 ~ ′ ′  ~
n1 + 1 
λ ε 

1/ 2
n1  λε 
1/ 4
(57)

 ~ n1ζ .
λ ε 
Для каждого алгоритма вычислений в
U(1) модели получим следующие оценки значений среднего числа операций 〈NopHMC〉 и 〈NopTSMB〉,
необходимых для перехода от одной конфигурации полей к последующей
~
(
N opHMC = 4VN τ 1 + 2 N CG
)+
+8V N CG ~ V N CG N τ ,
(58)
N opTSMB = 2V (128n1 + n2 + n3 ) ~ Vn2 .
(59)
Из (53) с учетом (37), (39), (55) и (58) для
метода HMC, и с учетом (49), (51), (57) и (59) для
алгоритма TSMB получим следующие оценки
производительностей этих алгоритмов применительно к U(1) модели при Nf = 2
1
1
PHMC ~
, PTSMB ~
. (60)
5/ 4
ζV ln V
ζV ln 2 V
ЛЕСНОЙ ВЕСТНИК 2/2007
Отсюда для отношения производительностей этих алгоритмов вычислений в U(1) модели имеем оценку [14]
V 1/ 4
PTSMB PHMC ~
.
(61)
ln V
Видно, что при увеличении объема решетки V U(1) модели двухшаговый мультибозонный
алгоритм становится лучшим по производительности, чем метод гибридного Монте-Карло.
5. Результаты вычислений и сравнение
двух алгоритмов
5.1. Параметры модели и параметры
алгоритмов
Нами были разработаны, отлажены и
применены для расчетов программы вычислений
средних значений корреляционных функций рассматриваемой U(1) модели по Вильсону на четырехмерной решетке в подходе динамических
фермионов как методом гибридного Монте-Карло (HMC), так и двухшаговым мультибозонным
алгоритмом (TSMB). По этим программам были
вычислены указанные в разделе 2 величины 〈EG〉,
〈 ψ ψ〉 и 〈Π〉 данной U(1) модели в областях Кулоновской фазы и фазы конфайнмента. Вычисления
проведены на суперкомпьютере VPP-300.
Были заданы следующие параметры U(1)
модели, одинаковые при вычислениях каждым
алгоритмом [14]: число поколений фермионов
Nf = 2; объем решетки V = 63×12; параметры
β = 2 и κ = 0.130 для Кулоновской фазы, β = 0 и
κ = 0.238 для фазы конфайнмента; число генерируемых полей U(α) N = 10000; периодические граничные условия по пространственным направлениям
ν = 1,2,3 и антипериодические по направлению
времени ν = 4 для фермионных полей.
Выбранные значения β и κ находятся
вблизи линии раздела фаз, что необходимо для
правильного описания свойств частиц. Кроме
того, интересно узнать характеристики алгоритмов при «трудном» для них значении β = 0.
Для расчета значений исходных параметров каждого алгоритма в каждой фазе были вычислены средние собственные значения 〈λmin〉 и
〈λmax〉 матрицы Q†Q в приближении статических
фермионов путем усреднения значений λmin и λmax
по калибровочным полям U с весом exp(−SG[U]),
как в формуле (22). Для каждой конфигурации
поля U значения λmin и λmax рассчитаны явно повторяемым методом Ланцоша (п. 3.4) при значениях M = 10 и в формуле (31) δ = 10−6. Генерация
119
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
конфигураций полей U с указанным весом выполнена методом «тепловой бани» при Φ = 0 (п. 3.2).
По найденным 〈λmin〉 и 〈λmax〉 в приближении статических фермионов были вычислены для
каждой фазы модели по (34), (37) и (38) предварительные значения числа обусловленности ζ и
следующих исходных параметров метода HMC:
Nτ, ∆τ, δmd и δacc.
Затем методом HMC были сгенерированы конфигурации калибровочных полей U, расN
пределенные с весом exp(−SG[U]) det M[U]. Для
каждой из таких конфигураций были снова найдены явно повторяемым методом Ланцоша значения λmin и λmax. Эти значения были усреднены
N
по U с новым весом exp(−SG[U]) det M[U], как в
формуле (10), и в результате для них были получены следующие значения, одинаковые в каждом
алгоритме [14]: 〈λmin〉 = 0.13(1) и 〈λmax〉 = 1.63(1)
в Кулоновской фазе, и 0,0005(1) и 6,59(1) в фазе
конфайнмента.
Используя эти окончательные значения
〈λmin〉 и 〈λmax〉 в подходе динамических фермионов,
по соотношениям (34, 37, 38) для метода HMC,
и по (40, 49, 51, 52) для алгоритма TSMB были
рассчитаны значения исходных параметров каждого алгоритма. Эти значения приведены в табл.
1 [14].
В случае алгоритма TSMB интегральным
методом наименьших квадратов (п. 3.3) найдены
коэффициенты аппроксимирующих многочленов
Pk(x), k = 1,2,3,4, а также найдены корни P1(x).
Расчетные значения данных многочленов удовлетворяют оценкам (41–43).
f
f
5.2. Значения корреляционных функций и
характеристики алгоритмов
Каждым из расматриваемых алгоритмов
с использованием формул (5–7, 10, 19, 20, 22)
были рассчитаны корреляционные функции 〈EG〉,
〈 ψ ψ〉 и 〈Π〉 данной U(1) модели в двух ее фазах.
Кроме того, были рассчитаны значения величин 〈wacc〉, 〈τint〉 способом суммирования по (54),
PTSMB / PHMC и значения среднего компьютерного
времени 〈tCPU〉 перехода от одной конфигурации
полей к последующей, равного
〈tCPU〉 = tc〈Nop〉.
(62)
Значения отношения PTSMB / PHMC были рассчитаны по формуле, следующей из (53) и (62)
PTSMB PHMC =
120
HMC
tCPU
HMC
τint
TSMB
tCPU
τTSMB
int
.
Таблица 1
Значения исходных параметров двух
алгоритмов для двух фаз
Фаза
Кулоновская
конфайнмента
∆τ
0,025
0,01
ε
0,025
Кулоновская
конфайнмента 0,000225
HMC
δmd
Nτ
40
10–3
10
10–3
TSMB
λ
n1
n2
n3
2,5
6
36
48
9
50
360 450
δacc
10–7
10–7
n4
200
500
Таблица 2
Значения корреляционных функций 〈O〉 в
двух фазах U(1) модели, вычисленные двумя
алгоритмами, и характеристики алгоритмов
〈EG〉
〈OHMC〉
〈OTSMB〉
〈waccHMC〉
〈waccTSMB〉
〈τintHMC〉
〈τintTSMB〉
〈tCPUHMC〉, с
〈tCPUTSMB〉, с
PTSMB / PHMC
〈OHMC〉
〈OTSMB〉
〈waccHMC〉
〈waccTSMB〉
〈τintHMC〉
〈τintTSMB〉
〈tCPUHMC〉, с
〈tCPUTSMB〉, с
PTSMB / PHMC
〈 ψ ψ〉
Кулоновская фаза
0,1332(1)
0,9381(1)
0,1331(1)
0,9379(1)
0,94(1)
0,48(1)
3,2(3)
2,0(2)
3,0(3)
2,8(2)
15,1(2)
8,96(2)
1,8(2)
1,2(2)
фаза конфайнмента
0,939(1)
0,95(1)
0,938(1)
0,96(1)
0,72(1)
0,68(1)
65(7)
60(7)
120(20)
125(15)
76(1)
69(1)
0,6(1)
0,5(1)
〈Π〉
1,378(1)
1,376(1)
25(4)
50(8)
0,8(2)
13,9(2)
13,7(2)
35(5)
45(5)
0,9(1)
Результаты выполненных расчетов приведены в табл. 2 [14]. Значения каждой из трех
корреляционных функций в каждой из двух фаз
модели, полученные с применением двух разных
алгоритмов, оказались одинаковыми. Это свидетельствует о пригодности каждого из рассматриваемых алгоритмов, отличающихся многими
сложными процедурами и формулами, для вычислений с существенным сокращением времени в U(1) модели по Вильсону на четырехмерной
решетке. Кроме того, это подтверждает правильность установленных соотношений для исходных
параметров каждого алгоритма. Особо отметим
совпадение результатов вычислений каждым алгоритмом и пригодность этих алгоритмов для
«трудного» случая β = 0 в фазе конфайнмента.
ЛЕСНОЙ ВЕСТНИК 2/2007
ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ
Полученные средние значения 〈wacc〉
для каждого алгоритма находятся в пределах
0,48 – 0.94, что соответствует принятой оценке
〈wacc〉 ∼ 1. Значения 〈τint〉 для каждой корреляционной функции оказались одного порядка величины для разных алгоритмов, причем для величин
〈EG〉 и 〈 ψ ψ〉 в фазе конфайнмента при β = 0 они
оказались на полтора порядка больше, чем в Кулоновской фазе при β = 2, а для корреляционной
функции 〈Π〉 – одного порядка величины. Значения 〈tCPU〉, пропорциональные продолжительности
вычислений, оказались одного порядка величины
в каждой фазе при разных алгоритмах, но в фазе
конфайнмента они были в 5–8 раз больше, чем в
Кулоновской фазе.
Отношения производительностей алгоритмов для шести рассмотренных случаев (три
корреляционные функции и две фазы) оказались
в пределах PTSMB / PHMC = 0,5 –1,8. Этот диапазон
соответствует нашей оценке (61), по которой для
данного объема решетки V = 63×12 отношение
PTSMB / PHMC ≈ 0,91. При увеличении объема решетки производительность каждого алгоритма,
согласно (60), уменьшается, а время вычислений,
соответственно, увеличивается. Но отношение
PTSMB / PHMC, согласно (61), возрастает, и двухшаговый мультибозонный алгоритм становится
лучшим по производительности, т.е. с меньшим
временем вычислений, чем метод гибридного
Монте-Карло.
6. Заключение
В данной работе рассмотрены два алгоритма для U(1) модели элементарных частиц на
четырехмерной решетке пространства-времени:
метод гибридного Монте-Карло и двухшаговый
мультибозонный алгоритм. Каждый из этих алгоритмов является статистическим методом получения калибровочных полей с весом, включающим
фермионный детерминант. Такие поля необходимы для вычисления корреляционных функций
модели. В каждый алгоритм входит метод стохастического вектора, и используется «шахматное»
разложение узлов решетки на четные и нечетные
узлы.
Вместе с тем, данные алгоритмы отличаются многими входящими в их состав методами,
процедурами и подходами. В метод гибридного
Монте-Карло входят использование вспомогательных полей – псевдофермионов и сопряженного импульса, процедура молекулярной динамики,
ЛЕСНОЙ ВЕСТНИК 2/2007
метод сопряженных градиентов и шаг Метрополиса для гамильтониана. Двухшаговый мультибозонный алгоритм включает использование
вспомогательных полей – мультибозонов, методы
«тепловой бани» и верхней релаксации, применение аппроксимирующих многочленов, явно повторяемый метод Ланцоша и шаг Метрополиса для
детерминанта.
Установлены и применены в расчетах соотношения для параметров и характеристик этих
алгоритмов применительно к рассматриваемой
U(1) модели. Это зависимости для шага «времени», числа шагов «времени» и параметров остановки метода сопряженых градиентов в методе
гибридного Монте-Карло, соотношения для границ интервала аппроксимации и для степеней четырех аппроксимирующих многочленов в случае
двухшагового мультибозонного алгоритма, а также зависимости для производительности каждого алгоритма. Данные соотношения необходимы
для обеспечения достаточно точных и быстрых
вычислений.
Разработаны, отлажены и реализованы
программы вычислений корреляционных функций U(1) модели каждым из указаных алгоритмов.
С помощью каждого алгоритма в областях Кулоновской фазы при β = 2 и фазы конфайнмента при «трудном» значении β = 0 вычислены
значения следующих корреляционных функций:
средней «энергии» калибровочного поля, скалярного конденсата и «пионной» нормы. Значения
корреляционных функций в каждой фазе, рассчитанные разными алгоритмами, оказались одинаковыми. Это свидетельствует о пригодности
каждого из рассмотренных альтернативных алгоритмов для исследований данной U(1) модели при
существенном сокращении времени вычислений,
а также подтверждает правильность полученных
соотношений для параметров алгоритмов.
Установлено аналитически и подтверждено численными расчетами, что для выбранной решетки 63×12 время вычислений корреляционных
функций с помощью каждого алгоритма является
примерно одинаковым. При увеличении объема
решетки время вычислений каждым алгоритмом
возрастает, но алгоритмом с меньшим временем
вычислений становится двухшаговый мультибозонный алгоритм.
Указанные алгоритмы с параметрами по
установленным зависимостям целесообразно ис-
121
Документ
Категория
Без категории
Просмотров
5
Размер файла
564 Кб
Теги
решетки, элементарные, алгоритм, статистический, части, модель, два
1/--страниц
Пожаловаться на содержимое документа