Два статистических алгоритма для u(1) модели элементарных частиц на решетке Текст научной статьи по специальности «Математика»

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

ДВА СТАТИСТИЧЕСКИХ АЛГОРИТМА ДЛЯ U(1) МОДЕЛИ ЭЛЕМЕНТАРНЫХ ЧАСТИЦ НА РЕШЕТКЕ

Н.В. ЗВЕРЕВ, доц. каф. физики МГУЛ, канд. физ.-мат. наук

Одной из важных задач математического моделирования и численных методов является создание и реализация эффективных численных методов и алгоритмов для исследования моделей элементарных частиц. Для такого исследования К. Вильсоном [2] предложен эффективный математический подход - метод решетки. В этом методе непрерывное пространство-время аппроксимируют дискретной совокупностью точек - узлов решетки, а полевые и корреляционные функции модели зависят от координат этих узлов.

Практический интерес представляет четырехмерная решеточная модель элементарных частиц, построенная на калибровочной группе матриц U(1) [2, 3]. При исследовании модели частиц прежде всего необходимо найти ее многочисленные корреляционные функции (функции Грина), которые затем используют для расчетов физических характеристик частиц.

Исследования моделей частиц на решетке прямыми численными методами являются нереальными из-за гигантского времени вычислений. К настоящему времени предложены два статистических алгоритма для моделей частиц на решетке с существенным сокращением времени вычислений: метод гибридного Монте-Карло [4,

5] и двухшаговый мультибозонный алгоритм [6]. Данные алгоритмы отличаются многими входящими в их состав методами, процедурами и подходами.

Точность и продолжительность расчетов и даже сама возможность вычислений с помощью любого алгоритма зависят от обоснованного выбора его исходных параметров. А соотношения для параметров указанных алгоритмов применительно к U(1) модели известны не были. Также не были выполнены сравнительные исследования U(1) модели этими алгоритмами и не была подтверждена расчетами пригодность каждого алгоритма для данной модели. Рассмотрение и решение этих математических проблем является весьма актуальным.

Целью данной работы являются анализ указанных двух статистических алгоритмов применительно к U(1) модели частиц на четырехмерной решетке, установление соотношений для

параметров и характеристик этих алгоритмов, вычисление с помощью каждого алгоритма значений корреляционных функций в двух фазах U(1) модели, выяснение пригодности рассматриваемых алгоритмов для такой модели по результатам вычислений, а также сравнение продолжительности вычислений каждым алгоритмом.

1. U(1) модель частиц на решетке

1.1. Действие и поля

Рассматриваемая модель элементарных частиц построена на матричной группе U(1), которая при определенных преобразованиях полей обеспечивает инвариантность фундаменталь-

ной физической величины - действия модели. Вследствие этой инвариантности выполняются необходимые физические законы сохранения при взаимодействиях и взаимопревращениях элементарных частиц. Действие S[U, у, у ] U(1) модели на четырехмерной решетке пространства-времени является суммой действия Sq[U] калибровочного поля, описывающего частицы - переносчики взаимодействия, и действия 5ДЦ,у, у ] фермионных полей, описывающих фермионные частицы [2]:

аду, у ] = Se[U] + 5ад,у, у ].

Здесь

Sq[U] = вЕ Re(1 -U^),

x,|x,v

Д<У

_ Nf —

SP [U, ^ у] = ЕЕ уХм [U ]x,, уf;

f=1 x, y

(1)

где Ux^v - решеточный компактный тензор напряженностей калибровочного поля:

U = U U U* a U* ;

Ux^ = exp(iAx ^) - решеточное компактное калибровочное поле;

Ax^ - вещественный потенциал калибровочного поля в интервале (-п, п];

Р = 1/е02 - обратный квадрат затравочного заряда;

уf и уX - фермионные поля, являющиеся ан-

тикоммутирующими переменными;

M [U] - комплексная фермионная матрица

M [U ] = 5 -

xyxy

-кЕ{(1 -Y )U 5 . + (1 + y )U* 5 . }; (2)

x,^ x + ц, y v у,ц y+p ,x у ’ V /

Д

ЛЕСНОЙ ВЕСТНИК 2/2007

109

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

где к = 1/(8+2т0) - хоппинг-параметр;

т0 - затравочная фермионная масса;

а - шаг решетки;

Y - эрмитовы матрицы Дирака размером 4x4;

f - индекс поколения фермионов;

Nf - число поколений фермионов;

5 - символ Кронекера;

х, у - узлы четырехмерной конечной решетки с целочисленными координатами n = 0,1,...,N -1;

p,v = 1,2,3,4 - направления решетки; ц - единичный вектор в положительном направлении ц;

N - число узлов решетки вдоль направления ц.

Полное число узлов, называемое объемом решетки, равно V = N1N2N3N4. Здесь принята система естественных единиц измерения, в которой квантовая постоянная Й, скорость света c и шаг решетки а выбраны равными Й = c = a = 1.

Обычно рассматривают решетку с числами узлов N1 = N2 = N3 = Ns < N Объем такой решетки равен V = Ns3xN4.

Для полей на конечной решетке принимают одно из следующих граничных условий (со знаком + или -)

^Х+NvV U

Для поля U выбирают периодические граничные условия (знак +) по всем направлениям v. Для полей у^ и у. принимают периодические граничные условия (знак +) по пространственным направлениям v = 1,2,3 и либо периодические, либо антипериодические (знак -) условия по направлению времени v = 4.

Матрица (2) удовлетворяет свойству у5-эрмитовости

Mt[U] = Y5M[U\y5 , (3)

где y5 = Y1Y2Y3Y4 - киральная матрица Дирака; эрмитово сопряжение t выполнено в пространстве матриц Дирака и в пространстве узлов решетки x.

1.2. Корреляционные функции и фазовая структура модели

Для расчетов физических характеристик частиц используют многочисленные корреляционные функции модели частиц (функции Грина). Каждую корреляционную функцию (O) рассчитывают, главным образом, методом континуального (функционального) интеграла. Этот метод заклю-

чается в усреднении соответствующей функции O[U] по всем полям U, у и у , распределенным с весом ехр(-£[и,у, у ]). При этом после интегрирования только по фермионным полям у и у получается следующая формула для корреляционной функции (O) модели на решетке [2]

O = Z j O[U]exp(-SG [U])detNf M[U] [dU ], (4)

где Z - нормировочная постоянная.

Здесь в силу свойства (3) детерминант фермионной матрицы (фермионный детерминант) deM[U] является вещественной величиной.

Для упрощения расчетов значений (O) по формуле (4) используют приближение статических фермионов, когда в диаграммах взаимодействий частиц пренебрегают вкладом фермионных петель, т.е. полагают det f M[U] = const. Для более точного расчета значений (O) используют подход динамических фермионов, когда в формуле (4) учитывают зависимость det fM[U] от U [7].

При исследовании моделей на решетке обычно вычисляют следующие корреляционные функции (O): среднюю «энергию» калибровочного поля (Eg), скалярный конденсат (у у) и «пион-ную» норму (П). Эти средние величины вычисляют по формуле (4), в которой зависимость O[U] для каждой из перечисленных функций определяется одной из следующих формул [1]

Eg =~^ I Re(1 - Ux,v), (5)

6V X ,Д,У

Д<У

уу = —1— TrM l[U], (6)

4V

П = — Tr y 5M 4[U\y 5M 4[U ], (7)

4V

где Tr обозначает след в пространстве матриц Дирака и в пространстве узлов решетки х. В силу свойства (3) зависимости (6) и (7) являются вещественными величинами, как и (5).

Частицы по U(1) модели на четырехмерной решетке имеют различные свойства в разных областях значений обратного квадрата затравочного заряда Р и хоппинг-параметра к [3]. Такие области называют фазами, из которых практический интерес представляют Кулоновская фаза и фаза конфайнмента.

Для правильного описания свойств частиц параметры р и к следует выбирать вблизи линий раздела фаз [1, 7]. На этих линиях решеточные корреляционные функции имеют определенные сингулярные свойства.

110

ЛЕСНОЙ ВЕСТНИК 2/2007

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

2. Метод гибридного Монте-Карло 2.1. Формула для корреляционных функций

Метод гибридного Монте-Карло предназначен для приближенных, но относительно быстрых вычислений корреляционных функций моделей на решетке путем статистического выбора конечного числа конфигураций калибровочного поля U и вспомогательных полей с весом в подходе динамических фермионов. При этом число поколений фермионов Nf должно быть четным. В данном методе вводят вспомогательные поля P и х таким образом, что формула (4) с учетом свойства (3) для фермионной матрицы M[U] переходит в следующее соотношение [4, 5]

O = Z j O[U] exp(-H[U, P, х]) [dU dPdх], (8)

где H[U,P,x] - гамильтониан, определенный формулой i

H [U, P, х] = Sq [U ] + - X P1 +

^ x,^.

Nf /2

+ Zхt (Mt [U]M [U ])х f . (9)

f=i

где Px - вещественные компоненты вспомогательного поля P, называемого сопряженным импульсом;

Xf - вспомогательные поля с 4V комплексными составляющими, называемые псевдофермионами.

Поля U, P и х необходимо генерировать случайным образом с весом exp^HUP^]). Тогда формула (8) переходит в следующее соотношение, по которому выполняют численные расчеты

1 N

O = lim-£ O[U(а)], (10)

' ' N N ,

1 ’ а=1

где O[U-a)] - значение функции O[U] при конфигурации U = £/а), состоящей из 4 V комплексных составляющих поля Ц/а) ;

а - номер конфигурации поля: а = =1,2,..., N; N - число генерируемых конфигураций поля

up

При расчетах число N является конечным. С увеличением N точность вычислений значений (O) возрастает.

2.2. Генерация полей

Каждую конфигурацию с номером а полей U, P и х с весом распределения, указанным к формуле (10), получают по следующей схеме [4, 5].

Сначала выбирают конфигурацию а сопряженного импульса P, генерируя случайным образом компоненты Px^ с весом exp(-P2x^ /2),

а также выбирают конфигурацию а псевдофермионов х по формуле Xf = М[Ц]ц, где - вектор с 4V комплексными составляющими, генерируемый случайным образом с весом exp(-nt/Hf) (f = 1, • ••, Nf /2). При этом конфигурация а поля U известна.

Такую генерацию Pxц и хf осуществляют методом стохастического вектора. В данном методе генерируют случайные комплексные числа z с весом exp(-C|z|2) по формуле

z = tJ C-1 ln E-1 exp(^),

где E и ф - вещественные случайные числа, равномерно распределенные в интервалах [0,1] и (-п, п].

Каждую компоненту Pxц сопряженного импульса P полагают равной вещественной части Re z при C = 1/2, а каждая из 4V комплексных составляющих вектора nf равна z при C = 1.

Затем выполняют процедуру молекулярной динамики. Эта процедура состоит в получении промежуточной U и вспомогательной P' конфигураций при известной конфигурации а полей U, P и а по следующим формулам

P(1/2) = P + — F [U];

x,n x,n 2 x,|H J’

Ax

U (k) = U (k-1)

x, Ц x, Ц

I P(k+1/2) = P(k-1/2)

(N,-1)

U' = U4"'

x,n x,|+

exp (/AtP*-1/2)),

At Fx^U (k) ];

exp (/AxPx(,N'-1/2)),

P' = P '

x, x,

(N,-1/2)

+

At

+TFx,, [U•].

Здесь At - шаг «времени», N, - число шагов «времени»; k = 1,... Р-1;

t/0) = U ; F [U]

x, ц x;|o/ x,|a l j

dH

dA

x,

Рассматриваемая схема заканчивается шагом Метрополиса для гамильтониана H. Этот шаг заключается в выборе поля U для конфигурации а+1 из поля U, найденного в процедуре молекулярной динамики, и поля U конфигурации а. Поле U принимают за новое (а+1) с вероятностью, равной

wJU ' ,U] = mln{1,exp(H[U,P,х]-H[U 'P ',х])>.

Для этого генерируют случайное число z, равномерно распределенное в интервале [0,1]. Если z < w [U,U], то полученное поле U принимают за U в конфигурации а +1. В противном случае за U в конфигурации а +1 принимают U прежней (а) конфигурации.

ЛЕСНОЙ ВЕСТНИК 2/2007

111

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

В результате выполнения этой схемы для всех N конфигураций U, P и х эти поля оказываются распределенными с весом exp(-#[Up,x]), при котором становится справедливым соотношение (10).

2.3. «Шахматное» разложение решетки

В процедуре молекулярной динамики и в шаге Метрополиса для гамильтониана необходимо умножать псевдофермионы xf на обратные матрицы М-1 и [МД-1, где с учетом свойства (3) матрицу [МД-1 выражают через Mr1. Это умножение выполняют путем решения систем уравнений вида ME, = п с заданным п и неизвестным Е векторами, каждый из которых имеет 4V комплексные составляющие. Для ускорения решения таких систем используют, согласно Т. де Гранду [8], «шахматное» разложение конечной решетки на четные и нечетные узлы с соответствующим разложением величин М, Е и п.

Для этого разложения на узлах решетки х вводят функцию с(х), принимающую значения +1 или -1 и удовлетворяющую соотношению с(х± Ц ) = —с(х) для всех узлов х и всех направлений ц. Такая функция с точностью до знака равна

Е Xv

c(х) = (—1)v .

Узел решетки х называют четным, если

сумма

Е х,

является четным числом и с(х) = 1, и нечетным, если сумма

Е х

V

нечетна и с(х) = — 1.

При четном числе узлов решетки N хотя бы для одного направления ц общее количество четных узлов равно количеству нечетных узлов. Это позволяет разложить фермионную матрицу М по (2) размера 4 V х 4 V на подматрицы одинакового размера 2 V х 2 V в подпространствах с четными и нечетными узлами следующим образом [8]

М =

1

e

М

oe

(11)

Здесь 1е - единичная матрица в подпространствах с четными или нечетными узлами; Мо и Мое - квадратные матрицы, первая из которых отображает подпространство с нечетными узлами решетки, обозначаемое индексом о, в подпространство с четными узлами, обозначаемое как e, а вторая отображает подпространство e в подпро-

странство о. Далее, разлагают векторы п и Е с 4V компонентами на векторы с 2V компонентами в подпространствах с четными и нечетными узлами

|Ч ^

п =

ч ^ )

л =

л

О )

Используя разложение этих векторов, а также разложение (11) фермионной матрицы М, решают следующую систему из 4V уравнений для неизвестных векторов E и Е , эквивалентную МЕ, = п

QE = п — М п ,

o o oe e

Е =п — М Е .

e e e

(12)

Здесь матрица Q размера 2Vx2V, действующая в подпространстве с нечетными узлами, определена по формуле

Q = 1 — М М . (13)

•s-' e oe eo v '

В силу «шахматного» разложения (11)

матрицы М справедливо тождество

detQ = detM (14)

Матрица Q обладает также свойством у5-

эрмитовости, которое следует из (3)

Q} = Y5QY5. (15)

Из системы (12) с учетом (13) видно, что для нахождения вектора Е достаточно решить линейную систему из 2V уравнений, верхнюю в (12) с неизвестным вектором Ео, в подпространстве с нечетными узлами. Вектор Ee элементарно выражается через Ео по системе из 2 V уравнений, нижней в (12).

2.4. Метод сопряженных градиентов

Нахождение вектора Ео в верхней системе уравнений (12) осуществляют методом сопряженных градиентов [7, 9]. Для использования этого метода верхнюю часть системы (12) путем ее умножения на эрмитово сопряженную матрицу Q^ приводят к следующей эквивалентной системе

QtQ^ = Q4 — МЧ».

Данная система имеет структуру в виде

Ах = у,

где А - эрмитова положительно-определенная матрица (в данном случае размера 2 Vx2V); у - заданный вектор (с 2V компонентами); х - неизвестный вектор (с 2V компонентами). Метод сопряженных градиентов представляет собой следующую итеративную схему построения векторов х , сходящихся к х

= х + a g , r = у - Ах , а =

П П~>пР П * пР п

g\ ASn

= r - а Ag , В =

n n <->n7 > n

^, gn+1 = rn+1 + Pngn. (16)

2

r

n

х

2

r

r

n

112

ЛЕСНОЙ ВЕСТНИК 2/2007

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

Здесь n = 0,1,...; ||a|| = у/а]а - норма вектора a; rn - векторы, называемые остатками; gn - векторы, называемые градиентами. Начальные векторы принимают равными g0 = r0 = y - Ax где вектор x0 генерируют методом стохастического вектора с весом exp(-||x0||2). Векторы-остатки rn являются попарно ортогональными.

Выполнение схемы (16) прекращают при достижении условия

||rjl < 5 (17)

где 5 - заданная малая величина: 5 << 1.

Эта величина для остановки метода сопряженных градиентов в процедуре молекулярной динамики равна 5 = 5md, а для остановки этого метода в шаге Метрополиса для гамильтониана 5 = 5acc. Согласно Р. Гупте и др. [10], не обходкмо выполнение условия 5 << 5 . а

Для векторов-остатков rn справедливо следующее неравенство

2pn ILII _ 4Z-1 - X,

r < -

r

2 n II 0|

P =

z=-

1 + p л/с + 1 X min

где X и X . - максимальное и минимальное собс-

' 1 msiY min

твенные значения положительно-определенной матрицы A (в данном случае Q]Q);

Z - число обусловленности матрицы A.

Из указанного неравенства для ||rn|| в случае Z >> 1 получим следующую оценку для числа итераций NCG) схемы метода сопряженных градиентов (16) до достижения условия остановки (17)

N(cg)-^Z In2-. (18)

- 2 5

3. Двухшаговый мультибозонный алгоритм 3.1. Формула для корреляционных функций

Двухшаговый мультибозонный алгоритм является статистическим методом приближенных и относительно быстрых вычислений корреляционных функций моделей на решетке в подходе динамических фермионов при произвольном значении числа поколений фермионов N В данный алгоритм введены вспомогательные поля - мультибозоны, мультибозонное действие и аппроксимирующие многочлены. При этом формула (4) преобразуется с учетом свойств (14) и (15) к следующему виду [6]

= i0[U] det-1P4 (Q\U]Q[U])signdetNf Q[U])ц

" (det-1P4 (Q\U]Q[U])signdetNf Q[U^ .

Здесь обозначено

(A) 12 = Z j A[U] exP( Tg [U] - (20)

-Та[Ф,U])det-1 P2 Q[U]Q[U])[dUdФ]; где SG[U] - действие калибровочного поля U по (1);

Q - матрица размером 2Vx2V, определенная согласно (13);

Та[Ф,K] - мультибозонное действие

n л л

Та [Ф, U] = — Ф] M] [U]Mj [U]Ф .;

j=1

Ф - вспомогательные поля с 4Vкомплексны-

j

ми составляющими, называемые мультибозонами;

Mj [U] - комплексные матрицы размером 4Vx4V, определяемые по формуле л Г 1 M

M j = 8 80

_Y5M« (Y5 -Р j )1e

матрицы 1e, Meo иMoe определены в (11); рj - комплексные корни вещественного многочлена P1(x2): P1(pj2) = 0, удовлетворяющие при j = 1,.,n1/2 условиям в виде pj+n/2 = -р* и Imp. > 0; P1(x), P2(x) и P4(x) - вещественные положительные многочлены степеней n1, n2 и n4, которые используют для грубой, промежуточной и точной аппроксимаций функции x Nf 2 на интервале xe[s, X]

P(x) * x-Nf/2, p(x)P2(x)-x~Nf/2,

P( x)P2( x)PA( x) = x-Nf /2 (21)

Интервал [s, X] содержит минимальное (X .) и максимальное (X ) средние собственные значения матрицы Q]Q при их усреднении по всем конфигурациям калибровочного поля U с весом exp(-TG[ U]) detNf M[ U],

Поля U и Ф необходимо генерировать случайным образом с весом

expbTJUbTa^UDdet-^mUMU]).

Тогда формула (20) переходит в следующее соотношение, по которому выполняются численные расчеты

1N

(A)12 = lim - — A[U(a)], (22)

' ' 12 N^« N 1

1 ’ a=1

где U(a), a и N - те же, что и в формуле (10).

3.2. Генерация полей

Каждую конфигурацию с номером a полей U и Ф с весом распределения, указанным к формуле (22), получают путем выполнения двух шагов.

Для этого вводят матрицы П+ и П_ - проекторы на четное и нечетное подпространства узлов решетки при ее «шахматном» разложении

ЛЕСНОЙ ВЕСТНИК 2/2007

113

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

Г1 0 1 Г0 0e1

е е П = е

А 0е _ ? _ А 1 _

п+ =

где 0е - нулевая матрица в подпространствах с четными или нечетными узлами.

Эти проекторы имеют следующие элементы в пространстве узлов решетки

Е ^

1 ± (-1)

(п±! =■

-5 .

'xy 2 С помощью данных проекторов матрицы

Мj принимают вид

М j = П+М + П_ч,М -р/П_ .

Затем действие УДФ,^] представляют в следующем виде:

SB [Ф, U ] = Я.

V.

Ф +jr-

+ S„

где S - слагаемое, не зависящее от величины

Ф

j,x,r’

Vjxr =kE(U* - ([Рj,xY5(1 + Y-)_А,Ф - +

x _Ц,Д j,x_-,s

+[Y5(1 + Y- )]rs X - ) + Ux,- ([Р j,xY5(1 _

j, x _-,-, s

_Y- ) _ 2LФ - + [Y 5(1 _Y- )L X - A

j,x+-,s j ,x+-,-,s

A = 1 + 4 x 4k2 +

j, x,r

Е x-

+

X.

1 _ (_1) 2

{ Рj I2 _ 2[Y5]rr Re Рj }ф

kE{[Y5(1 _Yv)]„Ux, ф

+

v ,s

+[Y5(1 + Yv )]„U* - Ф - },

x_v, v J ,x_v,s

Е+ Е+

1 + (_1)Д + 1 _ (_1)Д .

p =-------------p +------------p .

Fj,x 2 j 2 j

Здесь j = 1,...,^; x - узел решетки; r,s - номера строк и столбцов матриц Дирака: r,s = 1,2,3,4; матрица y5 диагональна по r и s, и все ее ненулевые элементы [Y5]rr = ±1.

Далее сумму действий SG[U] + +SB^,U] представляют в виде

SG[U] + ^[ОД = _Re(Ux/x,-) + SGB0 ,

где SGB0 - слагаемое, не зависящее от Ux д;

f =ру (U - U* - U* + U\ „ U\ U ) _

x,Д v x,v - - - - /

x+-,v x+v^ x_v+u,v x_v^ x_v,v

+

j, x+ ,s Д

_2kE {[Y5(1 _ Y-)]„ x (pj,xфj,x,rф

j ,r ,s

+ Ф,x,rX . + - +X',x,Д,ГФ . + - ) _ 25rsФ*,x,rФ . + - }.

J’ ’ j,x+ -,s j,x+ ,s J’ ’ j,x+ ,s

Д Д Д

2

j, x,r

j, x+,,s

v

В первом шаге каждую новую (а+1) конфигурацию поля Ф и промежуточную конфигурацию U получают из предыдущих (а) конфигураций U и Ф путем выполнения следующих четырех последовательных процедур методами «тепловой бани» (HB) и верхней релаксации (OR):

1. Метод HB для мультибозонов с получением 4V«1 составляющих поля Ф' из Ф для каждой совокупности индексов (j,x,r) по формулам

Ф' = —

V.

- + Ё

j,x,r ’

j ,x,r

где V^ = Р[Ф,Ц];

Ё. - комплексное случайное число, распре-

*7,x,r

деленное с весом exp(_A. |ё . |2) и гене-

j,x,r j,x,r

рируемое методом стохастического вектора (п. 2.2).

Найденное поле Ф' принимают за Ф.

2. Метод OR для мультибозонов с получением 4V«1 составляющих поля Ф' из найденного в процедуре 1 поля Ф для каждой совокупности (j,x,r) по формулам

V.

Ф' =_Ф _ 2-^

j ,x ,r !,x,r .

A

j, x,r

Найденное поле Ф' принимают за поле Ф в а+1 конфигурации.

3. Метод HB для калибровочного поля с получением 4V составляющих поля U' из U для каждой совокупности узлов и направлений решетки (x,p) при поле Ф в а+1 конфигурации по формулам

U = (F / |F |)_1exp(in ),

x,- v x,- 1 x,-17 1 v lx,-7’

где Fx - = F[Ф,U];

nx- - вещественное случайное число в интервале (_п,п], распределенное с весом

exP(|FJco?

Генерацию таких чисел осуществляют с использованием вспомогательного веса и метода Метрополиса [7]. Найденное поле U принимают за U.

4. Метод OR для калибровочного поля с получением 4V составляющих поля U' из поля U, найденного в процедуре 3, для каждой совокупности (x,p) при поле Ф в а+1 конфигурации по формулам

U = U* (F / IF |)_2.

Далее выполняют второй шаг алгоритма - метод Метрополиса для детерминантов [6]. Этот шаг заключается в выборе поля U для конфигурации а+1 из поля U, найденного в процедуре 4 первого шага, и поля Uконфигурации а. Поле U принимают за новое (а+1) с вероятностью, равной

114

ЛЕСНОЙ ВЕСТНИК 2/2007

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

^ [U ',U ]

(

= min 1,

V

det-1р (&[U']Q[U'])" det-1P2 (Q][U]Q[U])y

, (23)

где, согласно И. Монтвею [6], отношение детерминантов определяют по формуле

det-1P2 (Q\U']Q[U']) =

det-1P2 (Q][U]Q[U]) =

= (exp(-^fP((Qf[U' ]Q[U' ])% + Vn))n. (24)

Здесь обозначено

iB) n = Z J B[n] exp(-^n) [d n]; (25)

где % - вектор с 2 V комплексными составляющими;

% = P,mUQ[U])n;

P3(x) - вещественный положительный многочлен степени n3, аппроксимирующий с высокой точностью функцию P2-1/2(x) на интервале xe [в,Я]

P з(х) = P2-1/2(x). (26)

При расчетах вектор n генерируют случайным образом с весом exp(-^n) методом стохастического вектора (см. п. 3.2), затем вычисляют вектор % и отношение детерминантов по (24). Далее генерируют случайное число z, равномерно распределенное в интервале [0,1]. Полученное поле U принимают за U в конфигурации а, если z < wacc[U' ,U]. В противном случае за U в конфигурации а+1 принимают U прежней (а) конфигурации.

В результате выполнения первого и второго шагов алгоритма для всех N конфигураций U и Ф эти поля оказываются распределенными с требуемым весом

exp(-^[Uj-^,UDdet-1P((Qt[U]Q[U]).

При вычислениях величин (O) по (19) и (22) значения det-1P4(QQ) находят, согласно И. Монтвею [6], по формуле

det-1P 4QQ) = {exp(-n]PJ,Q]Q)n + Пп))„, (27)

где усреднение выполнено по (25).

3.3. Аппроксимирующие многочлены

Апппроксимирующие многочлены P1(x) с заданными степенями n где k = =1,2,3,4, вычисляют, согласно И. Монтвею [11], по следующей схеме.

Выбирают интервал аппроксимаций [в, Я] так, чтобы он содержал средние минимальное (Я . ) и максимальное (Я ) собственные значения матрицы QQ: [(Ятт),(Ят£к)] с [в,Я]. Эти значения (Ятт) и (Я ) вычисляют по формулам (19) и (22),

где в (19) величина O заменяется на Я или Я .

Далее, интегральным методом наименьших квадратов находят коэффициенты аппроксимирующих многочленов P1(x), P2(x) и P4(x) при заданных степенях n n2 и n4

N,/2

x ’

inm, in,

xN/2 P(x)P2(x1-1

-—min,

О фиксирован,

xNf,2P ( p)P*2 (p)P4 (p) -1 -

>min,

Р1, Р2 фиксированы, (28)

где |[/(x)|| - интегральная норма функции /x):

Я V/(

/ (x) =

— J dx | /(x) |(

Я - В ;

s у

(29)

Минимизацию каждого интеграла (29) в (28) выполняют по коэффициентам многочлена P1(x), P2(x) и P4(x) соответственно.

Нахождение коэффициентов многочлена P3(x) заданной степени n3 осуществляют итерационным методом, аналогичным методу Ньютона касательных

P3(j+1) (x) = 1 (p3( 2) (x) + P3( 2) (x)), j = 0,1,....

Здесь Pjx) и P3(j) (x) - многочлены степени n3 для всех j. Многочлен P3(j) (x) аппроксимирует функцию (P2(x) Pjx))-1. Его коэффициенты находят интегральным методом наименьших квадратов

||p, (x)P3( j> (x)P3( j> (x) -1 -- - min,

Р2, P3j фиксированы, (30)

Минимизацию интеграла (29) для (30) выполняют по коэффициентам многочлена P3(j) (x). В пределе j — да многочлен P3j)(x) — P3(x).

После вычисления коэффициентов многочленов находят методом Лежандра корни р. многочлена P1(x2), необходимые для первого шага данного алгоритма.

Вычисление матриц Pk(QQ), где k = 2,3,4, осуществляют, согласно И. Монтвею [11], по следующим формулам, в которых x заменяется на QQ

nk

Рк (x) = £ cj.k Of>(x), k = 2,3,4;

j=0

O+1 (x) = (x + p jk> )O(k} (x) + Y ™O™ (x).

Здесь O(k)(x) - многочлены степени j, удовлетворяющие условию ортогональности

Я

J dx p(k> (x) O\k> (x) Ojk) (x) = 0, l Ф j,

s

ЛЕСНОЙ ВЕСТНИК 2/2007

115

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

где

O * Ч X) = 1, Y -1 = 0,

р(2)( X) = X^Pfi X),

Р(3)( X) = P,2( x)P,2( X),

p(4)( X) = XNfPl1( X) P2( X).

3.4. Явно повторяемый метод Ланцоша

Для выбора границ интервала [в, Я] необходимо вычислять для каждой конфигурации поля U минимальное X и максимальное X

min max

собственные значения матрицы QQ. Эти собственные значения находят явно повторяемым методом Ланцоша [12].

Этот метод в общем случае эрмитовой матрицы A представляет собой итерации по следующей схеме:

1. При первой итерации n = 1 принимают единичный вектор e0(n) произвольного направления. Задают целое число M> 2.

2. Для данного номера итераций n вычисляют M единичных векторов ejn) (j = 0,1,.. .M"— 1) по формуле

Ae(n)

J

(n) (n) , (n) (n)

~ : + a e

j+1 j j

P j e

+p

( n ) e( n ) j—1 J—1 ’

где

a

(n) =

(e (n))tAe

(n),

II, e—1(n) = 0.

j j j

e (n) = ||Ae (")—a we n_В (”)e n I

j j j j j—1 j—1

3. Строят эрмитову трехдиагональную матрицу T(n) размером M х M с элементами a(n) (j = 0,1,...M—1) на главной диагонали и Р,(и) (j = 0,.. .M—2) на соседних диагоналях.

4. Вычисляют методом Ритца либо минимальное, либо максимальное собственное значение X(n) матрицы T(n) и соответствующий единичный собственный вектор-столбец s(n): T(n)s(n) = X(n)s(n).

5. Проверяют выполнение следующего условия для номера итераций n>1:

|X(n) —X(n—1)| < 5, (31)

где 5 - заданное малое число.

Если это условие выполнено, то вычисления останавливают и принимают полученное X(n) либо за X , либо за X собственное значение

min max

матрицы A. В противном случае переходят к следующему шагу.

6. Строят новый исходный вектор e0(n+1) для следующей (n+1) итерации:

( n +1)

e

0

's(я )e(n >,

I j j ’

M —1

j=0

где s(n) - компоненты вектора-столбца s(n).

7. Возвращаются к шагу 2 при n ^ n+1 и

e (n+1) e (n)

e0 e0 .

В случае данного алгоритма матрица A = QQ

3.5. Вероятность принятия новых калибровочных полей

Средняя вероятность (wacc) принятия новых калибровочных полей U в двухшаговом мультибозонном алгоритме равна отношению числа полей U, принятых за новые, к полному числу генерируемых полей. Применяя методы А. Борелли и М. Кройца, с учетом выражений (23) и (24), получим при ((AE)2) < 1 следующее соотношение

<»•„>=1ЧtS (ae )!) • (32)

где

<(AE)2) - V^2(x) —1II2, (33)

где |[/(x)|| - интегральная норма по (29);

AE = ifP2 Q [U' ]Q[U']}^ — п2П;

S = P3 (Q2[U]Q[U])q;

C =1 j C[U' , Ф ' ,U, Ф, q] X

х exp (— Se [U ] — SB [Ф, U ] — n4n)X

X p12[U' , Ф ', U, Ф] [dU d Ф dU ' d Ф ' dq]; p12[U' ,Ф' ,U,0] - плотность вероятности перехода от полей U и Ф конфигурации а к полям U и Ф ' конфигурации а+1.

4. Исходные параметры и производительность

двух алгоритмов для U(1) модели

Для исследования U(1) модели каждым из двух описанных алгоритмов и их сравнения необходимо обоснованно выбрать значения исходных параметров этих алгоритмов. В противном случае время вычислений существенно увеличится и результаты расчетов могут быть неправильными. Ниже получим соотношения для значений исходных параметров и для производительности рассматриваемых алгоритмов при одинаковом заданном числе поколений фермионов N/ = 2.

4.1. Параметры метода гибридного Монте-Карло

Для расчетов методом гибридного Монте-Карло должны быть заданы следующие его исходные параметры: шаг «времени» Ат, число шагов «времени» N малые параметры в условии (17) остановки метода сопряженных градиентов 5md в молекулярной динамике и 5acc в шаге Метрополиса для гамильтониана.

Для получения необходимых соотношений рассмотрим положительно-определенную

116

ЛЕСНОЙ ВЕСТНИК 2/2007

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

матрицу Q размером 4Vx4V, задаваемую формулой

Q =

( д2H Л' ^ д2 и Л

д2 H

1дА J

1дА J

-11/2

где H - гамильтониан по (9);

А - векторный потенциал поля U; усреднение выполнено по полям U, Р и х с весом

exp(-H[UP,x]).

Оценим норму этой матрицы ||Q|| при Nf = 2. Для этого используем свойства симметрии U(1) модели относительно дискретных сдвигов и поворотов, формулу для средних значений псевдофермионов в виде (х'х) = MM и формулу для собственных значений X. матрицы MM в виде X. (MM) = ХД(ММ)1/2]. Кроме того, сумму по собственным значениям X. заменим интегрированием, а также учтем, что для матриц MM и Q}Q справедливы, согласно де Гранду [14], оценки в виде

(х_(м’м)) ~ z (x,je'Q))

(xmn(M'M)}~ Z (x„„(e'e^'

X _(M 'м )~1. (34)

Здесь Z - число обусловленности матрицы QQ. В результате получим следующую оценку:

I 1 \ (Xmax) d x

И~(т>-(м,м)-|>~(1Д~ J -Хт~£, (35)

' 2 j ' (Xmin)

где X , X и X - собственные значения матрицы

у min max г

(MM)1/2.

При исследовании модели свободного бозонного поля на решетке методом гибридного Монте-Карло А. Кеннеди и Б. Пендлетон [13] получили следующие соотношения, обеспечивающие высокие точность и производительность вычислений при средней вероятности принятия новых полей (w ) ~ 1

V||Q||2 (Ах)4 < 1^j^|NxAx ~1, (36)

где

Q =

д2 H дА2

- матрица, не зависящая от свободного бозонного поля A;

H - гамильтониан бозонного поля;

V - объем решетки.

Мы применим к U(1) модели соотношения (36), в которых величину ||Q|| возьмем по (35). Тогда получим [14]

Ат<-

(V Z)1

, N > V1

(37)

1

Для значений параметров 5md и 5acc принимаем, согласно Р. Гупте и др. [10], следующие оценки

5, ~—,5 ~^, (38)

md v ~ acc v 2 ' 4 '

Итак, получены соотношения (37) и (38) для исходных параметров метода гибридного Монте-Карло применительно к U(1) модели при N = 2.

Из соотношений (18) и (38) имеем следующую оценку для значения среднего числа итераций (N(CG)) в методе сопряженных градиентов

(N (cg ])~^Z ln V. (39)

4.2. Параметры двухшагового мультибозонного алгоритма

Для расчетов двухшаговым мультибозонным алгоритмом должны быть заданы следующие его исходные параметры: границы в и X интервала аппроксимаций функции x Nf 2 по (21) и степени и n2, n3 и n4 аппроксимирующих многочленов.

Границы в и X интервала аппроксимаций выберем следующим образом [14]

В = WmmX Х = (12 - 1.4)(Xmax), (40)

где (Xmin) и (Xmax) - средние минимальное и максимальное собственные значения матрицы QQ, усредненные, как в (19).

Такой выбор учитывает необходимость нахождения крайних значений (X ) - сХ и

(Xmax) + СХ внутри интервала [в, X], где сХ и

Gx - среднеквадратичные отклонения.

max

Согласно И. Монтвею [6], высокие точность и производительность данного алгоритма обеспечиваются при условиях

(wacc) - 1, (det-^ee)) - 1,

где усреднение вероятности wacc выполнено так же, как в (32), а усреднение величины det-1P4(Q'Q) - как в (19). При этих условиях из (27), (32) и (33) получим следующие оценки для многочленов P2(x) и P4(x)

||P2(х) - 1||~-р||Рл(X) -1||~ -1, (41)

где норма |fx)|| определена по (29).

Точные аппроксимации (21) и (26) обеспечиваются, согласно И. Монтвею [6], при условиях в виде

Р2(х)Р2(х) -1

V

2

1

xNf/2 Р( х) Р2( X) Р4 (х) -1

1

V ■

(42)

ЛЕСНОЙ ВЕСТНИК 2/2007

117

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

Используя (41) и (42), получим следующие оценки

1

y[V ’

xNf /2р(x) -1

xNf/2P1(x)P2(x) -1

1

V'

(43)

Для оценки значений степеней nk многочленов Pk(x) введем, удобные для анализа многочлены Pk (x) степеней Пк, составленные из многочленов Чебышева Tn +1 (x), где к = 1,2,3,4. Многочлены Pk (x) при к = 1,2,4 аппроксимируют функцию x f на отрезке [s',Я'] с одинаковым небольшим отклонением от 1 функций xNf 2Pk (x), а многочлен P3(x) аппроксимирует функцию P2(x)P32(x)

P(x) * P(x) * x‘Nf/2,

P2(x) * p(x)P, (x)-x-Nf/2,

P53(x) * P> (x)P32(x) = 1,

P,(x) * P (x)P2 (x)P4 (x) = x~

Здесь

(44)

P (x) = 1

x

f

1-k

T.

T /2 *-( X+s')V

+4 V-s' /

(45)

/ V+s'\

, ‘nk+1VV-s7/

где TA +j (x)=cos[(Пк +1)arccos x].

Согласно нашим расчетам при Nf = 2, многочлены p (x) и P1(x) практически одинаково аппроксимируют функцию x Nf /2 при выполнении следующих условий

П = n1, s' -yfsV, V ' = V. (46)

Для многочленов Чебышева справедлива оценка в виде

' V + s''

2( Пk +1)V s' / V

T„ +11 -V' ' I ~ . (47)

V -s _

Из определения числа обусловленности Z по (34) с учетом (40) следует оценка в виде

V ~ z .

s

(48)

Из (43), первого соотношения (44) и (45) - (48) получим

X . .4 -1

0, |П,*-1/4 1

-2(n1+1)Z 1

[ V' + s')

xNf р(x) -1 ~ T 1 I

nk У V-s')

W’

откуда имеем оценку степени n1 [14]

n1 - Z1/4lnV. (49)

Расчеты при Nf = 2 показали, что многочлены P2 (x) и P1(x)P2(x) практически одинаково аппроксимируют функцию x Nf /2 при условиях

n2 = n1 + n2, s '~ s, V ' = V. (50)

Из (43), второго соотношения (44), а также (45), (47), (48) и (50) получим

N./2

[ V + s'I

xNf P2( x) -1 ~ T n +1

1 V -s )

'е

-2(«1 +n2 +1)Z

1/2

1

V ’

откуда имеем оценку для степени n2 [14]

n2 - Z1/2lnV. (51)

Результаты аналогичных расчетов значений многочленов P3 (x), P4 (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. В расчетах корреляционных функций (О) используются только аппроксимирующие многочлены Pk(x), степени которых nk выбирают по (49, 51, 52).

4.3. Производительность алгоритмов

Производительностью алгоритма вычислений конкретной корреляционной функции (O) модели на решетке назовем величину P [14], равную

P = (^рХ^ (53)

где (N ) - среднее число компьютерных операций, необходимых для перехода от одной конфигурации полей к последующей; (iint) - интегральное время автокорреляций, равное среднему расстоянию между ближайшими статистически независимыми выборочными значениями 0[Ц-а)] этой корреляционной функции [7].

По способу суммирования [9] это время вычисляют по формуле

= 1+1 у

2 W t? Го (0)

(54)

где W - некоторое обрезающее число вычисле-

I

int

ний: 1 << W << N ;

Г0(а) - автокорреляционная функция

1 N -а / 1 N-а \( i N

Г0 (а) = O )[0+а- — У о

N-а j.=! ( N-а k=1 А N-ак=а+1

Здесь Оа = 0[Ц^а)]; остальные обозначения те же, что и в формуле (10).

Время вычисления величины (О) определяется формулой

t = t (N ) N = ,

c \ °р/ p ’

где N. - число статистически независимых значений О : N. = N/ (т. .);

tc - время выполнения одной компьютерной операции.

118

ЛЕСНОЙ ВЕСТНИК 2/2007

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

Лучшим алгоритмом является такой, у которого большая производительность Р при одинаковых параметрах модели N V, в, к и числе N генерируемых конфигураций полей.

В случае метода гибридного Монте-Карло (HMC) на модели свободного бозонного поля на решетке А. Кеннеди и Б. Пендлетон [13] получили следующую оценку времени автокорреляций: <TintHMC) ~ (N Ах)-2. (55)

Мы так же, как и в п. 4.1, применим это соотношение к U(1) модели, в которое подставим наши оценки (37).

В случае двухшагового мультибозонного алгоритма (TSMB), согласно М. Пиардону, К. Янсену и А. Бориси, для времени автокорреляций справедливо соотношение

(Tmt)~ ^mln |Im р)\ j , (56)

где р2 - j-й корень аппроксимирующего многочлена РДх): РДр.2) = 0.

Для рассматриваемой U(1) модели при N = 2 заменим, согласно (44), многочлен РДх) на Р (х) и с учетом (45) найдем корни р. многочлена Р (х), при которых р (pj2) =0

2 /т ' . ■ 2 nj ■ ГГГ~1 • 2nj

р . = (к +8 ) sin---+ 1у1 к 8 sin-,

П + 1 П + 1

j = 1,2,..., n,.

Подставляя это выражение в (56) и учитывая (40, 46, 48 и 34), получим следующую оценку ( ------------------ 2п ^ n (к'

к 8 sin -

n +1J

к

VkT

V 8 J

л f f

к V 8

y\.

nZ

(57)

Для каждого алгоритма вычислений в U(1) модели получим следующие оценки значений среднего числа операций (NopHMC) и (NopTSMB), необходимых для перехода от одной конфигурации полей к последующей

(NHPMC) = 4 VN (1 + 2(NCG))+

+8V(NCo) ~ V{Nco)Nt , (58)

N™B) = 2V(128n, + n2 + n3)~ Vn2. (59)

т

int

1/2

n

Из (53) с учетом (37), (39), (55) и (58) для метода HMC, и с учетом (49), (51), (57) и (59) для алгоритма TSMB получим следующие оценки производительностей этих алгоритмов применительно к U(1) модели при N, = 2

Р

HMC ZV5/4 In V

Р

’ TSMB

1

ZVIn2 V

(60)

1

Отсюда для отношения производительностей этих алгоритмов вычислений в U(1) модели имеем оценку [14]

Р /Р

TSMB / HI

V1

In V

(61)

Видно, что при увеличении объема решетки V U(1) модели двухшаговый мультибозонный алгоритм становится лучшим по производительности, чем метод гибридного Монте-Карло.

5. Результаты вычислений и сравнение двух алгоритмов

5.1. Параметры модели и параметры алгоритмов

Нами были разработаны, отлажены и применены для расчетов программы вычислений средних значений корреляционных функций рассматриваемой U(1) модели по Вильсону на четырехмерной решетке в подходе динамических фермионов как методом гибридного Монте-Карло (HMC), так и двухшаговым мультибозонным алгоритмом (TSMB). По этим программам были вычислены указанные в разделе 2 величины (EG), (f f) и (П) данной U(1) модели в областях Кулоновской фазы и фазы конфайнмента. Вычисления проведены на суперкомпьютере VPP-300.

Были заданы следующие параметры U(1) модели, одинаковые при вычислениях каждым алгоритмом [14]: число поколений фермионов N, = 2; объем решетки V = 63х12; параметры в = 2 и к = 0.130 для Кулоновской фазы, в = 0 и к = 0.238 для фазы конфайнмента; число генерируемых полей U-a) N = 10000; периодические граничные условия по пространственным направлениям v = 1,2,3 и антипериодические по направлению времени v = 4 для фермионных полей.

Выбранные значения в и к находятся вблизи линии раздела фаз, что необходимо для правильного описания свойств частиц. Кроме того, интересно узнать характеристики алгоритмов при «трудном» для них значении в = 0.

Для расчета значений исходных параметров каждого алгоритма в каждой фазе были вычислены средние собственные значения (ктт) и (ктах) матрицы Q}Q в приближении статических фермионов путем усреднения значений к и к по калибровочным полям U с весом exp(-£0[и|), как в формуле (22). Для каждой конфигурации поля U значения к и к рассчитаны явно повторяемым методом Ланцоша (п. 3.4) при значениях M = 10 и в формуле (31) 5 = 10-6. Генерация

ЛЕСНОЙ ВЕСТНИК 2/2007

119

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

конфигураций полей U с указанным весом выполнена методом «тепловой бани» при Ф = 0 (п. 3.2).

По найденным (к . ) и (X ) в приближении статических фермионов были вычислены для каждой фазы модели по (34), (37) и (38) предварительные значения числа обусловленности Z и следующих исходных параметров метода HMC: N, Ах, 5 и 5 .

Затем методом HMC были сгенерированы конфигурации калибровочных полей U, распределенные с весом exp(-SG[U])det fM[U]. Для каждой из таких конфигураций были снова найдены явно повторяемым методом Ланцоша значения к и к . Эти значения были усреднены

min max

по Uс новым весом exp(-SG[U]) det f M[U], как в формуле (10), и в результате для них были получены следующие значения, одинаковые в каждом алгоритме [14]: (Xmin) = 0.13(1) и (Xmax) = 1.63(1) в Кулоновской фазе, и 0,0005(1) и 6,59(1) в фазе конфайнмента.

Используя эти окончательные значения (Xmin) и (Xmax) в подходе динамических фермионов, по соотношениям (34, 37, 38) для метода HMC, и по (40, 49, 51, 52) для алгоритма TSMB были рассчитаны значения исходных параметров каждого алгоритма. Эти значения приведены в табл.

1 [14].

В случае алгоритма TSMB интегральным методом наименьших квадратов (п. 3.3) найдены коэффициенты аппроксимирующих многочленов Pk(x), к = 1,2,3,4, а также найдены корни P1(x). Расчетные значения данных многочленов удовлетворяют оценкам (41-43).

5.2. Значения корреляционных функций и характеристики алгоритмов

Каждым из расматриваемых алгоритмов с использованием формул (5-7, 10, 19, 20, 22) были рассчитаны корреляционные функции (EG), ( у у) и (П) данной U(1) модели в двух ее фазах. Кроме того, были рассчитаны значения величин (w ), (т ) способом суммирования по (54), PTSMB / PHMC и значения среднего компьютерного времени (tCPU) перехода от одной конфигурации полей к последующей, равного

(tCPU) = tc(Np). (62)

Значения отношения Ptsmb / Phmc были рассчитаны по формуле, следующей из (53) и (62)

P /Р =

* TSMB I * HMC

HMC

HMC

t

т

CPU

int

TSMB

TSMB

t

т

CPU

int

Т а б л и ц а 1

Значения исходных параметров двух алгоритмов для двух фаз

Фаза HMC

Ат N 5 d 5

Кулоновская конфайнмента 0,025 40 10-3 10-7

0,01 10 10-3 10-7

TSMB

8 к n П3 П4

Кулоновская конфайнмента 0,025 2,5 6 36 48 200

0,000225 9 50 360 450 500

Т а б л и ц а 2

Значения корреляционных функций (O) в двух фазах U(1) модели, вычисленные двумя

алгоритмами, и характеристики алгоритмов

(Eg) (у у) (П)

Кулоновская фаза

(OHMC) 0,1332(1) 0,9381(1) 1,378(1)

(QTSMB) 0,1331(1) 0,9379(1) 1,376(1)

(w HMC) 0,94(1)

(w TSMB) 0,48(1)

(TinHMC) 3,2(3) 2,0(2) 25(4)

(y™) 3,0(3) 2,8(2) 50(8)

(tcPuHMC), с 15,1(2)

(tCPrTSMB), с 8,96(2)

P / P 1,8(2) 1,2(2) 0,8(2)

фаза конфайнмента

(OHMC) 0,939(1) 0,95(1) 13,9(2)

(QTSMB') 0,938(1) 0,96(1) 13,7(2)

(w HMC) 4 acc ' 0,72(1)

(w TSMB) acc 0,68(1)

(т inlHMC) 65(7) 60(7) 35(5)

(т nTSMB) 120(20) 125(15) 45(5)

(tcPUHMC), с 76(1)

tf TSMB\ n \1cpu /, с 69(1)

P / P TSMB ± HMC 0,6(1) 0,5(1) 0,9(1)

Результаты выполненных расчетов приведены в табл. 2 [14]. Значения каждой из трех корреляционных функций в каждой из двух фаз модели, полученные с применением двух разных алгоритмов, оказались одинаковыми. Это свидетельствует о пригодности каждого из рассматриваемых алгоритмов, отличающихся многими сложными процедурами и формулами, для вычислений с существенным сокращением времени в U(1) модели по Вильсону на четырехмерной решетке. Кроме того, это подтверждает правильность установленных соотношений для исходных параметров каждого алгоритма. Особо отметим совпадение результатов вычислений каждым алгоритмом и пригодность этих алгоритмов для «трудного» случая в = 0 в фазе конфайнмента.

120

ЛЕСНОЙ ВЕСТНИК 2/2007

ФУНДАМЕНТАЛЬНЫЕ И ПРИКЛАДНЫЕ ИССЛЕДОВАНИЯ

Полученные средние значения (wacc) для каждого алгоритма находятся в пределах 0,48 - 0.94, что соответствует принятой оценке (wacc) ~ 1. Значения (xmt) для каждой корреляционной функции оказались одного порядка величины для разных алгоритмов, причем для величин (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) модели элементарных частиц на четырехмерной решетке пространства-времени: метод гибридного Монте-Карло и двухшаговый мультибозонный алгоритм. Каждый из этих алгоритмов является статистическим методом получения калибровочных полей с весом, включающим фермионный детерминант. Такие поля необходимы для вычисления корреляционных функций модели. В каждый алгоритм входит метод стохастического вектора, и используется «шахматное» разложение узлов решетки на четные и нечетные узлы.

Вместе с тем, данные алгоритмы отличаются многими входящими в их состав методами, процедурами и подходами. В метод гибридного Монте-Карло входят использование вспомогательных полей - псевдофермионов и сопряженного импульса, процедура молекулярной динамики,

метод сопряженных градиентов и шаг Метрополиса для гамильтониана. Двухшаговый мультибозонный алгоритм включает использование вспомогательных полей - мультибозонов, методы «тепловой бани» и верхней релаксации, применение аппроксимирующих многочленов, явно повторяемый метод Ланцоша и шаг Метрополиса для детерминанта.

Установлены и применены в расчетах соотношения для параметров и характеристик этих алгоритмов применительно к рассматриваемой U(1) модели. Это зависимости для шага «времени», числа шагов «времени» и параметров остановки метода сопряженых градиентов в методе гибридного Монте-Карло, соотношения для границ интервала аппроксимации и для степеней четырех аппроксимирующих многочленов в случае двухшагового мультибозонного алгоритма, а также зависимости для производительности каждого алгоритма. Данные соотношения необходимы для обеспечения достаточно точных и быстрых вычислений.

Разработаны, отлажены и реализованы программы вычислений корреляционных функций U(1) модели каждым из указаных алгоритмов.

С помощью каждого алгоритма в областях Кулоновской фазы при р = 2 и фазы конфайнмента при «трудном» значении р = 0 вычислены значения следующих корреляционных функций: средней «энергии» калибровочного поля, скалярного конденсата и «пионной» нормы. Значения корреляционных функций в каждой фазе, рассчитанные разными алгоритмами, оказались одинаковыми. Это свидетельствует о пригодности каждого из рассмотренных альтернативных алгоритмов для исследований данной U(1) модели при существенном сокращении времени вычислений, а также подтверждает правильность полученных соотношений для параметров алгоритмов.

Установлено аналитически и подтверждено численными расчетами, что для выбранной решетки 63х12 время вычислений корреляционных функций с помощью каждого алгоритма является примерно одинаковым. При увеличении объема решетки время вычислений каждым алгоритмом возрастает, но алгоритмом с меньшим временем вычислений становится двухшаговый мультибозонный алгоритм.

Указанные алгоритмы с параметрами по установленным зависимостям целесообразно ис-

ЛЕСНОЙ ВЕСТНИК 2/2007

121