Формула Харди–Рамануджана

Формула Харди–Рамануджана — это асимптотическое, а в уточнённой форме точное выражение для вычисления функции разбиений \( p(n) \), разработанное в 1918 году английскими математиками Годфри Харолдом Харди и Сринивасой Рамануджаном. Данная формула стала первым значительным достижением аналитической теории чисел в исследовании аддитивных функций и положила начало методу кругов (методу Харди–Литлвуда).

Определение

Формула Харди–Рамануджана даёт асимптотическое представление для числа неупорядоченных разбиений натурального числа \( n \), обозначаемого \( p(n) \). В первоначальном виде она записывается как:

\[ p(n) \sim \frac{1}{4\pi\sqrt{3}} \exp\left( \pi \sqrt{\frac{2n}{3}} \right) \quad \text{при} \quad n \to \infty \]

В дальнейшем Харди и Рамануджан вывели уточнённую версию в форме бесконечного ряда:

\[ p(n) = \frac{1}{\pi \sqrt{2}} \sum_{k=1}^\infty A_k(n) \sqrt{k} \cdot \frac{d}{dn} \left( \frac{ \sinh\left( \frac{\pi}{k} \sqrt{ \frac{2}{3} \left( n - \frac{1}{24} \right) } \right) }{ \sqrt{ n - \frac{1}{24} } } \right) \]

которая является точной для всех целых \( n \geq 1 \).

Исторический контекст

Исследования разбиений чисел берут начало в работах Готфрида Вильгельма Лейбница и Леонарда Эйлера. Последний в 1748 году определил производящую функцию:

\[ \sum_{n=0}^\infty p(n) x^n = \prod_{k=1}^\infty \frac{1}{1 - x^k} \]

Длительное время эффективные способы вычисления \( p(n) \) для больших значений \( n \) отсутствовали. В 1914 году Рамануджан приехал в Кембридж и начал совместную работу с Харди. В период с 1916 по 1917 год они создали метод кругов, который позволяет анализировать производящие функции на границе круга сходимости посредством локального анализа вблизи корней из единицы.

В 1918 году учёные опубликовали труд «Asymptotic formulae in combinatory analysis», где впервые была представлена асимптотическая формула для \( p(n) \), а затем её усовершенствованный вариант.

Компоненты формулы

Для любого целого \( n \geq 1 \) точная формула имеет вид:

\[ p(n) = \frac{1}{\pi \sqrt{2}} \sum_{k=1}^\infty A_k(n) \sqrt{k} \cdot \frac{d}{dn} \left( \frac{ \sinh\left( \frac{\pi}{k} \sqrt{ \frac{2}{3} \left( n - \frac{1}{24} \right) } \right) }{ \sqrt{ n - \frac{1}{24} } } \right) \]

где используются следующие компоненты.

Функция разбиений \( p(n) \)

Функция \( p(n) \) определяет количество способов представления числа \( n \) в виде суммы натуральных чисел без учёта порядка слагаемых. По определению \( p(0) = 1 \).

Примеры значений:

• \( p(1) = 1 \)
• \( p(2) = 2 \)
• \( p(3) = 3 \)
• \( p(4) = 5 \)
• \( p(5) = 7 \)

Арифметическая сумма \( A_k(n) \)

\[ A_k(n) = \sum_{\substack{0 \leqslant m < k \\ \gcd(m,k)=1}} \exp\left( \pi i s(m,k) - 2\pi i \frac{nm}{k} \right) \]

• Условие \( \gcd(m,k)=1 \) означает, что суммирование проводится по всем вычетам \( m \) по модулю \( k \), взаимно простым с \( k \).
• Количество слагаемых равно значению функции Эйлера \( \varphi(k) \).
• \( A_k(n) \) является алгебраическим числом, причём \( |A_k(n)| \leq \varphi(k) \).

Сумма Дедекинда \( s(m,k) \)

\[ s(m,k) = \sum_{r=1}^{k-1} \left( \left( \frac{r}{k} \right) \right) \left( \left( \frac{mr}{k} \right) \right), \quad ((x)) = \begin{cases} x - \lfloor x \rfloor - \frac{1}{2}, & x \notin \mathbb{Z} \\ 0, & x \in \mathbb{Z} \end{cases} \]

Свойства:

• \( s(m,k) \in \mathbb{Q} \)
• \( s(-m,k) = -s(m,k) \)
• \( s(m+k,k) = s(m,k) \)
• При \( \gcd(m,k) = 1 \): \( s(m,k) = s(m^{-1}, k) \), где \( m^{-1} \) — обратный элемент по модулю \( k \).

Аналитическая часть

• Величина \( u(n) = n - \frac{1}{24} \) представляет собой модулярную поправку, связанную с функцией эта Дедекинда:

 \[
 \eta(\tau) = q^{1/24} \prod_{m=1}^{\infty} (1 - q^m), \quad q = e^{2\pi i \tau}
 \]

• Аргумент гиперболического синуса: \( B_k(n) = \frac{\pi}{k} \sqrt{ \frac{2}{3} u(n) } \)
• Производная по \( n \) возникает при применении формулы Коши и дифференцировании параметрического интеграла.

Вывод формулы

1. **Производящая функция**:

  \[
  P(q) = \sum_{n=0}^\infty p(n) q^n = \frac{1}{(q; q)_\infty} = \frac{q^{1/24}}{\eta(\tau)}, \quad q = e^{2\pi i \tau}
  \]

2. **Интеграл Коши**:

  \[
  p(n) = \frac{1}{2\pi i} \oint_{|q|=r} \frac{P(q)}{q^{n+1}} \, dq
  \]

3. **Метод кругов**:

  * Окружность разбивается на мажорные дуги вблизи точек \( q = e^{2\pi i m/k} \) и минорные дуги между ними.
  * Вклад минорных дуг пренебрежимо мал.
  * Вклад мажорной дуги около \( m/k \) выражается через модулярное преобразование \( \eta(-1/\tau) = \sqrt{-i\tau} \, \eta(\tau) \).

4. **Асимптотика интеграла** даёт слагаемое, зависящее от \( k \), \( m \) и \( n \). 5. **Суммирование** по всем \( m \), взаимно простым с \( k \), приводит к \( A_k(n) \). 6. **Суммирование** по всем \( k \) даёт итоговую формулу.

Формула Радемахера

В 1937 году Ханс Радемахер упростил формулу, получив абсолютно сходящийся ряд без производной:

\[ p(n) = \frac{1}{\pi \sqrt{2}} \sum_{k=1}^{\infty} A_k(n) \sqrt{k} \cdot \frac{d}{dn} \left( \frac{ I_{3/2}\left( \dfrac{\pi}{k} \sqrt{ \dfrac{2}{3} \left( n - \dfrac{1}{24} \right) } \right) }{ \sqrt{ n - \dfrac{1}{24} } } \right) \]

Используя соотношение \( I_{3/2}(z) = \sqrt{\dfrac{2}{\pi z}} \sinh(z) \), можно показать эквивалентность с формулой Харди–Рамануджана.

Свойства и сходимость

• Ряд сходится абсолютно для всех \( n \geq 1 \).
• Для вычисления \( p(n) \) с погрешностью менее 0.5 достаточно взять \( k \leq C \sqrt{n} \), где \( C \approx 1 \).
• При \( k=1 \): \( A_1(n) = 1 \), и соответствующее слагаемое даёт асимптотическую формулу.
• Формула является точной: при достаточном количестве слагаемых результат представляет собой целое число.

Уточнение приближения с ростом числа слагаемых

Точность формулы Харди–Рамануджана (или эквивалентной формулы Радемахера) значительно повышается с увеличением количества слагаемых \( K \) в сумме. Это иллюстрирует следующая таблица для трёх значений \( n \).

\( n \)	\( K \) (число слагаемых)	Приближение \( p(n) \)	Фактическое \( p(n) \)	Погрешность (%)
10	1	48	42	14,533
	2	43		2,519
	3	42		0,050
50	1	217 590	204 226	6,544
	2	205 239		0,496
	3	204 257		0,015
100	1	199 280 893	190 569 292	4,571
	2	191 032 426		0,243
	3	190 619 341		0,026

Как видно из таблицы, даже добавление второго и третьего слагаемых (\( K=2, 3 \)) существенно снижает погрешность, что подтверждает: чем больше \( K \), тем точнее приближение.

Применения

• **Комбинаторика**: точный подсчёт разбиений чисел.
• **Теория чисел**: пример применения метода кругов.
• **Физика**:

 * В **теории струн** \( p(n) \) описывает спектр состояний.
 * В **термодинамике чёрных дыр** (энтропия Бекенштейна–Хокинга) энтропия пропорциональна \( \log p(n) \).

Реализация формулы на языке Python

Ниже приведён код на Python для вычисления \( p(n) \) с использованием формулы Харди–Рамануджана. Параметр `n` задаёт число для разбиения, а `t` — количество знаков после запятой. Например, при \( n = 1\,000\,000 \) и \( t = 10 \) вывод будет «1.471684986e+1107», а при \( t = 10\,000 \) результат будет точным без использования степеней.

```python import mpmath import math

mpmath.mp.dps = 80

j = mpmath.j pi = mpmath.pi

def u(n):

   return n - mpmath.mpf(1)/24

def B(k, n):

   return (pi / k) * mpmath.sqrt(mpmath.mpf(2)/3 * u(n))

def ProizB(k, n):

   return (pi / (2 * k)) * mpmath.sqrt(2 / (3 * u(n)))

def Proizddn(k, n):

   Bval = B(k, n)
   term1 = mpmath.cosh(Bval) * ProizB(k, n) / mpmath.sqrt(u(n))
   term2 = mpmath.sinh(Bval) / (2 * (u(n) ** (3/2)))
   return term1 - term2

def s(m, k):

   if k == 1:
       return mpmath.mpf('0')
   total = mpmath.mpf('0')
   for r in range(1, k):
       term1 = mpmath.mpf(r)/k - mpmath.mpf('0.5')
       mr_mod_k = (m * r) % k
       term2 = mpmath.mpf(mr_mod_k)/k - mpmath.mpf('0.5')
       total += term1 * term2
   return total

def A(k, n):

   total = mpmath.mpc('0')
   for m in range(0, k):
       if math.gcd(m, k) == 1:  
           exponent = pi * j * s(m, k) - 2 * pi * j * (mpmath.mpf(n) * m / k)
           total += mpmath.exp(exponent)
   return total

def p(n):

   n = mpmath.mpf(n)
   itog = mpmath.mpc('0')
   for k in range(1, 5):
       akn = A(k, n)
       sqrt_k = mpmath.sqrt(k)
       deriv = Proizddn(k, n)
       itog += akn * sqrt_k * deriv

       if abs(akn * sqrt_k * deriv) < abs(itog) * 1e-20 and k > 1:
           break
   result = itog / (pi * mpmath.sqrt(2))
   if abs(result.imag) < abs(result.real) * 1e-10:
       return result.real
   return result

result = p(n) print(mpmath.nstr(result, t)) ```

Литература

1. Hardy, G. H.; Ramanujan, S. (1918). «Asymptotic formulae in combinatory analysis». Proceedings of the London Mathematical Society.
2. Rademacher, H. (1937). «A convergent series for the partition function p(n)». Proceedings of the National Academy of Sciences.
3. Apostol, T. M. (1990). Modular Functions and Dirichlet Series in Number Theory. Springer.
4. Andrews, G. E. (1998). The Theory of Partitions. Cambridge University Press.