Zeta function

Сегодняшняя заметка посвящена Роминой задаче. Решить ее с наскоку не получилось (хотя я уже близок и идеологически все хорошо понимаю, и вообще ускорение сходимости, например, с помощью преобразования Бореля, для физиков это everyday life, то есть быт. Спойлеров в комментариях буду строго карать.) Зато с наскоку получилось вывести формулу, которая хоть и выглядит монструозно, но поставленную задачу решает.

UPD. Исправил опечатки, добавил выражение через числа Бернулли и сравнение с формулой Эйлера-Маклорена.

А задача такова:

Просуммировать N членов некого ряда и получить численное значение, совпадающее с дзета функцией в N десятичных знаках.

То есть наш ряд должен давать экспоненциальную точность. В силу врожденной тупости и недостатка времени для серьезного вдумывания, я решил задачу просто в лоб, т.е. пересуммировал все степенные поправки для функции ζ(α) − H_α(N), где H_α(N) это просто гармонические числа:

Hα(N) =  N ∑ n = 1   1 nα

(1)

Давайте начнем с выражения для остатка в интегральном виде:

ζ(α) − Hα(N) =  1     ∞ ∫ 0   dx  xα-1 e−N x  . Γ(α) ex − 1

(2)

Для вычисления интеграла (2) давайте воспользуемся следующей формулой:

1  =  +∞ ∑ n = −∞   (−1)n  =  1  +  +∞ ∑ n = 0   2y (−1)n  . sinh y y − iπn y y2 + (πn)2

(3)

Собственно, эвристический вывод этой формулы уже написан  —  мы просто взяли и просуммировали по всем полюсам с учетом вычетов. Такие формулы я научился выводить, когда преподавал “Квантовую физику” в НГУ. Например, в теории сверхпроводимости БКШ для определения параметра щели Δ(T) в пределе когда температура стремиться к критической T → T_c, необходимо разрешить следующее выражение в пределе Δ → 0:

ln  Tc  =    ∞ ∫ 0   dξ  tanh (ξ/2T)  −  tanh [(ξ2 + Δ2)1/2/2T]  . T ξ (ξ2 + Δ2)1/2

(4)

Для этого проще всего воспользоваться выражением:

tanh (πx)  =  2x +∞ ∑ n = 0   1  . π x2 + (n + 1/2)2

(5)

После этого выражение (4) легко раскладывается в ряд по Δ, интегрируется и решается.

На самом деле легко дать эвристическому рассмотрению формул (3) и (5) строгое математическое доказательство, воспользовавшись цепочкой рассуждений называемой приемом Герглотца. Формулу (5) удобнее доказывать для мнимых аргументов, т.е. сделать замену x → ix. После этого давайте левую часть равенства (5) обозначим за f(x) а правую за g(x). Прием Герглотца состоит в доказательстве цепочки утверждений:

1. f(x) и g(x) определены и непрерывны во всех точках за исключением полуцелых точек.
2. f(x) и g(x) периодические с периодом единица.
3. f(x) и g(x) удовлетворяют одному и тому же функциональному уравнению f(x/2)+f(x/2+1/2) = 2f(x+1/2).
4. h(x) = f(x)−g(x) дополненная условием h(n+1/2) = 0 непрерывна везде и обладает всеми вышеперечисленными свойствами. Полюса мы все уже вычли, а ноль легко доказать прямо из выражения (5) и предела (для которого надо дважды применить правило Лопиталя)
   lim_x→1/2[π tan (πx) − (x − 1/2)⁻¹] = 0.
5. Пусть h(x₀+1/2) = m максимум функции h на интервале [1/2, 3/2], тогда согласно принципу Дирихле h(x₀/2) = h(x₀/2+1/2) = m. Продолжая процесс итеративно, получим h(x₀/2ⁿ+1/2) = m или h(1/2) = m стало быть m = 0.

Вернемся к нашим баранам. Используя формулу (3) получаем следующее выражение:

∞ ∫ 0   dx  xα-1 e−N x  =  Γ(α − 1)  + +∞ ∑ n = 0 2(−1)n     ∞ ∫ 0   dx  xα-1 e− (N+1/2) x  . ex − 1 (N + 1/2)α − 1 x2 + (2πn)2

(6)

Поскольку в интеграле существенна только область x ∼ 1/N разлагаем в ряд по x

1  =  1 +∞ ∑ k = 0 (−x2)k  . x2 + (2πn)2 (2πn)2 (2πn)2k

(7)

Интегрируя по x и суммируя по n получаем:

ζ(α) = Hα(N) +  1  −  2 +∞ ∑ k = 0   (−1)k Γ(α + 2k + 1)(1 − 2−2k−1)  ζ(2k + 2)  . (α − 1)(N + 1/2)α − 1 Γ(α) (2π)2k+2 (N + 1/2)α + 2k + 1

(8)

“Какого черта?!”— спросите вы меня, —“Начинали с дзета-функции и под знаком суммы опять получили дзета-функции.” Ну что ж, как говориться, “такова селя ви”. Единственное, что я могу придумать в свое оправдание, так это то, что слева стоит дзета-функция произвольного нецелого аргумента, а справа стоят дзета функции в четных точках, которые есть ничто иное как четные степени π деленные на целые числа. Из выражения (8) (по‐моему) видно, что утверждение Ромы, о том, что остаток раскладывается в ряд по 1/N неверно. В силу того, что Γ-функция растет очень быстро, ряд в выражении (6) расходиться, то есть является асимптотическим. В общем, это было понятно с самого начала — Рома утверждал, что имеются экспоненциальные поправки, а функция exp(−N) в ряд по 1/N не раскладывается, поскольку убывает, быстрее чем любая степень.

Асимптотические ряды  —  это то с чем физики живут в обнимку всю свою жизнь. Практически все ряды теории возмущений являются асимптотическими. (Далее автор ударяется в мутные пояснения) Рассмотрим глубокую локализованную яму и основное состояние волновой функции в ней. Давайте добавим еще одну такую яму далеко за радиусом локализации первой. Поскольку волновая функция основного состояния торчит из первой ямы совсем немного, то вторую яму можно попытаться рассмотреть как возмущение и вычислять поправки к уровню энергии основного состояния просто по теории возмущений. Однако продолжать эти вычисления до бесконечности у вас не получиться  —  ряд теории возмущений начнет расходиться. Физическая причина этого ясна: дело в том, что в присутствии новой ямы волновая функция основного состояния совсем не похожа на прежнюю. За счет туннелирования волновая будет поровну разделена по этим двум ямам. Вместо одного уровня, у нас появится два  —  и все что мы считали по теории возмущений это систематический сдвиг этих обоих уровней. Когда точность теории возмущений начнет сравниваться с величиной расщепления  —  ряд теории возмущений начнет расходиться, поскольку величина расщепления экспоненциальна мала (понятно по какому параметру) и в ряд разложена быть не может. Это типичное проявление непертурбативных эффектов (процесс туннелирования между различными вакуумами КХД называется инстантоном).

Исходя из этой идеологии, становиться понятно, что мы должны оборвать ряд (8) на том члене, когда поправки становятся экспоненциально малы. Пососав палец и вспомнив формулу Стирлинга, оборвем ряд на k = N−1. В итоге получим:

ζ(α, N) =  1  + N ∑ n = 1   1  −  (−1)n   2Γ(α + 2n − 1)(1 − 21−2n)  ζ(2n)  . (α − 1)(N + 1/2)α − 1 nα (2π)2n Γ(α)(N + 1/2)α+2n−1

(9)

Формула (9) решает сформулированную выше задачу, как говорят в этих ваших интернетах, чуть более чем полностью. При суммировании N + 1 члена ряда мы получим 2N + 2 десятичных цифр. Формулу (9) можно еще чуть чуть упростить, выразив дзета-функции четного аргумента через числа Бернулли.

ζ(α) =  (− 1)n − 1(2π)2n   B2n . 2Γ(2n + 1)

(10)

В этом случае формула (9) будет выглядеть так:

ζ(α, N) =  1  + N ∑ n = 1   1  −  Γ(α + 2n − 1)(1 − 21−2n)   B2n  . (α − 1)(N + 1/2)α − 1 nα Γ(2n + 1)Γ(α) (N + 1/2)α+2n−1

(11)

Естественно, мы могли бы сразу воспользоваться определением чисел Бернулли:

x  =  N ∑ n = 1   xn  Bn , ex − 1 n!

(12)

и тогда сразу бы получили выражение

ζ(α, N) =  −  1  + N ∑ n = 1   1  −  Γ(α + 2n − 1)   B2n  , 2Nα − 1 nα Γ(2n + 1)Γ(α) Nα+2n−1

(13)

но тогда у меня не было бы повода рассказать про прием Герглотца, а между тем это и есть самый простой способ получить связь (10) между числами Бернулли и дзета-функцией четного аргумента.

В комментариях Рома резонно заметил, что формула (9) до боли напоминает формулу Эйлера-Маклорена, которая, будучи применена к остатку

∞ ∑ n = 0   1  , (N + 1 + n)α

(14)

дает следующее приближение:

ζEM(α, N) = Hα(N) +  1  +  1  + N − 1 ∑ n = 1   Γ(α + 2n − 1)    B(2n)  . (α − 1)(N + 1)α − 1 2α(N + 1)α Γ(2k+1)Γ(α) (N + 1)α+2n−1

(15)

Я специально оборвал ряд на N − 1 члене, чтобы общее число слагаемых в этих двух формулах совпадало. Заметим, что аналога члена (N + 1)^−α в формуле (9) нет.

С точки зрения экспоненциальной точности формулы (9), (13) и (15) эквивалентны. Давайте определим следующую функцию относительной ошибки:

Err(α, N) = 2    ζ(α) − ζ(α, N)   ζ(α) + ζ(α, N)

(16)

Ниже на рисунке она показана в логарифмическом масштабе для α = 1.1: красным формула (9), зеленым (13), а синим Эйлер с Маклореном (15)