Problem Summary

Problem 825 asks for the partial sum

$$T(N)=\sum_{n=2}^{N} s_n,\qquad N=10^{14},$$

where the rational sequence \(s_n\) comes from the chasing game and is generated by two linear recurrences. A direct summation over \(10^{14}-1\) terms is impossible, so the solution extracts the exact dominant asymptotic term, sums that term in closed form, and keeps only a short rapidly convergent correction.

Mathematical Approach

Write

$$s_n=\frac{a_n}{b_n},$$

with initial values

$$a_2=7,\quad a_3=35,\quad a_4=121,$$

$$b_2=11,\quad b_3=73,\quad b_4=395,\quad b_5=1933,$$

and recurrences

$$a_n=3a_{n-1}+3a_{n-2}-a_{n-3}\qquad (n\ge 5),$$

$$b_n=8b_{n-1}-18b_{n-2}+8b_{n-3}-b_{n-4}\qquad (n\ge 6).$$

The implementations are built around the structure hidden in these recurrences.

Step 1: Factor the characteristic polynomials

The numerator recurrence has characteristic polynomial

$$r^3-3r^2-3r+1=(r+1)(r^2-4r+1).$$

The denominator recurrence has characteristic polynomial

$$r^4-8r^3+18r^2-8r+1=(r^2-4r+1)^2.$$

So the key roots are

$$\lambda=2+\sqrt3,\qquad \mu=2-\sqrt3=\lambda^{-1}.$$

The repeated factor in the denominator is the reason an extra linear factor in \(n\) appears there, and that is exactly what turns the ratio \(a_n/b_n\) into a shifted harmonic term.

Step 2: Solve the recurrences explicitly

Using the initial values, the closed forms become

$$A=\frac{3-\sqrt3}{2},\qquad B=\frac{3+\sqrt3}{2},$$

$$a_n=A\lambda^n+B\mu^n-2(-1)^n,$$

$$b_n=\left(An-\frac12\right)\lambda^n+\left(Bn-\frac12\right)\mu^n.$$

These formulas satisfy both the recurrences and the listed starting values. They also make the asymptotic behavior transparent: \(\lambda>1\) dominates, while \(\mu<1\) contributes only exponentially tiny corrections.

Step 3: Extract the shifted harmonic main term

Divide numerator and denominator by \(\lambda^n\):

$$s_n=\frac{A+B\mu^{2n}-2(-1)^n\lambda^{-n}}{An-\frac12+\left(Bn-\frac12\right)\mu^{2n}}.$$

Since \(\mu^{2n}=\lambda^{-2n}\), all non-dominant pieces decay exponentially. Therefore

$$s_n=\frac{A+O(\lambda^{-n})}{An-\frac12+O(n\lambda^{-2n})}.$$

Because

$$\frac{1}{2A}=\frac{3+\sqrt3}{6},$$

it is natural to define

$$c=\frac{3+\sqrt3}{6}.$$

Then

$$An-\frac12=A(n-c),$$

and hence

$$s_n=\frac{1}{n-c}+O(\lambda^{-n}).$$

The crucial point is that the correction is exponentially small, not merely polynomially small.

Step 4: Split the huge sum into a closed form and a short correction

Now decompose

$$T(N)=\sum_{n=2}^{N}\frac{1}{n-c}+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right).$$

The second sum converges very quickly because the summand is exponentially small. If we truncate it at a modest cutoff \(K\), we get

$$E_K=\sum_{n=2}^{K}\left(s_n-\frac{1}{n-c}\right),$$

and the omitted tail is negligible for the requested \(8\)-decimal output. The implementations choose

$$K=60.$$

Step 5: Sum the main term with the digamma identity

The digamma function satisfies

$$\psi(x+1)-\psi(x)=\frac1x.$$

Applying this telescoping identity to the shifted harmonic part gives

$$\sum_{n=2}^{N}\frac{1}{n-c}=\psi(N+1-c)-\psi(2-c).$$

So the desired value is

$$T(N)=\psi(N+1-c)-\psi(2-c)+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right),$$

and numerically it is enough to use

$$T(N)\approx \psi(N+1-c)-\psi(2-c)+E_{60}.$$

Worked Example: the checkpoint \(T(10)\)

The first few terms are

$$s_2=\frac{7}{11}\approx 0.63636364,\qquad s_3=\frac{35}{73}\approx 0.47945205,\qquad s_4=\frac{121}{395}\approx 0.30632911.$$

With

$$c=\frac{3+\sqrt3}{6}\approx 0.78867513,$$

the shifted harmonic part up to \(10\) is

$$\sum_{n=2}^{10}\frac{1}{n-c}\approx 2.54851473.$$

The correction over the same range is

$$\sum_{n=2}^{10}\left(s_n-\frac{1}{n-c}\right)\approx -0.16616191.$$

Therefore

$$T(10)\approx 2.54851473-0.16616191=2.38235282,$$

which matches the numerical checkpoint used by the implementations.

How the Code Works

The C++, Python, and Java implementations build the rational terms only up to \(n=60\). They generate the numerator and denominator sequences from the recurrences using exact integer arithmetic first, and only convert to floating point when forming each ratio \(s_n\) inside the correction sum. This keeps the short precomputation stable and accurate.

After the correction \(E_{60}\) is accumulated, the remaining task is to evaluate two digamma values: one at \(N+1-c\) and one at \(2-c\). The implementation does this numerically in two stages. First, it repeatedly uses the recurrence \(\psi(x)=\psi(x+1)-1/x\) until the argument is safely large. Then it applies the asymptotic expansion

$$\psi(x)=\log x-\frac{1}{2x}-\frac{1}{12x^2}+\frac{1}{120x^4}-\frac{1}{252x^6}+\frac{1}{240x^8}-\frac{1}{132x^{10}}+\frac{691}{32760x^{12}}+\cdots.$$

Finally the program combines

$$\psi(N+1-c)-\psi(2-c)+E_{60}$$

and prints the result to \(8\) decimal places.

Complexity Analysis

If the cutoff is denoted by \(K\), generating all needed recurrence values and summing the correction costs \(O(K)\) time and \(O(K)\) memory. The two digamma evaluations use only constant-time arithmetic and a fixed number of asymptotic-series terms. Since the implementation fixes \(K=60\), the whole computation is effectively constant-time even though \(N=10^{14}\) is enormous.

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=825
  2. Digamma function: Wikipedia — Digamma function
  3. Linear recurrence with constant coefficients: Wikipedia — Linear recurrence with constant coefficients
  4. Asymptotic expansion: Wikipedia — Asymptotic expansion
  5. Bernoulli numbers in special-function expansions: Wikipedia — Bernoulli number

Problemzusammenfassung

Problem 825 verlangt die Partialsumme

$$T(N)=\sum_{n=2}^{N} s_n,\qquad N=10^{14},$$

wobei die rationale Folge \(s_n\) aus dem Chasing Game stammt und durch zwei lineare Rekurrenzen erzeugt wird. Eine direkte Summation über \(10^{14}-1\) Glieder ist ausgeschlossen, daher isoliert die Lösung den exakten dominanten Asymptotikterm, summiert diesen in geschlossener Form und ergänzt nur noch eine kurze, schnell konvergente Korrektur.

Mathematischer Ansatz

Wir schreiben

$$s_n=\frac{a_n}{b_n},$$

mit Anfangswerten

$$a_2=7,\quad a_3=35,\quad a_4=121,$$

$$b_2=11,\quad b_3=73,\quad b_4=395,\quad b_5=1933,$$

und Rekurrenzen

$$a_n=3a_{n-1}+3a_{n-2}-a_{n-3}\qquad (n\ge 5),$$

$$b_n=8b_{n-1}-18b_{n-2}+8b_{n-3}-b_{n-4}\qquad (n\ge 6).$$

Der gesamte Rechenweg beruht auf der algebraischen Struktur dieser beiden Folgen.

Schritt 1: Charakteristische Polynome zerlegen

Für den Zähler lautet das charakteristische Polynom

$$r^3-3r^2-3r+1=(r+1)(r^2-4r+1).$$

Für den Nenner gilt

$$r^4-8r^3+18r^2-8r+1=(r^2-4r+1)^2.$$

Die entscheidenden Wurzeln sind also

$$\lambda=2+\sqrt3,\qquad \mu=2-\sqrt3=\lambda^{-1}.$$

Dass der Nenner einen doppelten dominanten Faktor besitzt, erklärt den zusätzlichen linearen Faktor \(n\), der später zur verschobenen harmonischen Form führt.

Schritt 2: Rekurrenzen explizit lösen

Aus den Startwerten ergeben sich die geschlossenen Formen

$$A=\frac{3-\sqrt3}{2},\qquad B=\frac{3+\sqrt3}{2},$$

$$a_n=A\lambda^n+B\mu^n-2(-1)^n,$$

$$b_n=\left(An-\frac12\right)\lambda^n+\left(Bn-\frac12\right)\mu^n.$$

Diese Darstellungen erfüllen sowohl die Rekurrenzen als auch die Anfangswerte. Gleichzeitig zeigen sie direkt das asymptotische Verhalten: \(\lambda>1\) dominiert, während \(\mu<1\) nur exponentiell kleine Beiträge liefert.

Schritt 3: Den verschobenen harmonischen Hauptterm isolieren

Teilt man Zähler und Nenner durch \(\lambda^n\), erhält man

$$s_n=\frac{A+B\mu^{2n}-2(-1)^n\lambda^{-n}}{An-\frac12+\left(Bn-\frac12\right)\mu^{2n}}.$$

Da \(\mu^{2n}=\lambda^{-2n}\), verschwinden alle nichtdominanten Anteile exponentiell. Also gilt

$$s_n=\frac{A+O(\lambda^{-n})}{An-\frac12+O(n\lambda^{-2n})}.$$

Weil

$$\frac{1}{2A}=\frac{3+\sqrt3}{6},$$

setzen wir

$$c=\frac{3+\sqrt3}{6}.$$

Dann ist

$$An-\frac12=A(n-c),$$

und daraus folgt

$$s_n=\frac{1}{n-c}+O(\lambda^{-n}).$$

Die Korrektur ist also exponentiell klein und nicht nur von einer Potenzordnung in \(n\).

Schritt 4: Die riesige Summe zerlegen

Nun schreiben wir

$$T(N)=\sum_{n=2}^{N}\frac{1}{n-c}+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right).$$

Die zweite Summe konvergiert sehr schnell. Trunkieren wir sie bei einer kleinen Schranke \(K\), dann ist

$$E_K=\sum_{n=2}^{K}\left(s_n-\frac{1}{n-c}\right)$$

bereits eine ausgezeichnete Näherung, und der Rest ist für \(8\) Nachkommastellen vernachlässigbar. Die Implementierungen wählen

$$K=60.$$

Schritt 5: Den Hauptterm mit der Digamma-Funktion summieren

Für die Digamma-Funktion gilt die Identität

$$\psi(x+1)-\psi(x)=\frac1x.$$

Damit ergibt sich für den verschobenen harmonischen Anteil

$$\sum_{n=2}^{N}\frac{1}{n-c}=\psi(N+1-c)-\psi(2-c).$$

Somit ist

$$T(N)=\psi(N+1-c)-\psi(2-c)+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right),$$

und numerisch reicht

$$T(N)\approx \psi(N+1-c)-\psi(2-c)+E_{60}.$$

Durchgerechnetes Beispiel: der Kontrollwert \(T(10)\)

Die ersten Glieder sind

$$s_2=\frac{7}{11}\approx 0.63636364,\qquad s_3=\frac{35}{73}\approx 0.47945205,\qquad s_4=\frac{121}{395}\approx 0.30632911.$$

Mit

$$c=\frac{3+\sqrt3}{6}\approx 0.78867513$$

ist der verschobene harmonische Anteil bis \(10\)

$$\sum_{n=2}^{10}\frac{1}{n-c}\approx 2.54851473.$$

Die zugehörige Korrektur lautet

$$\sum_{n=2}^{10}\left(s_n-\frac{1}{n-c}\right)\approx -0.16616191.$$

Daher

$$T(10)\approx 2.54851473-0.16616191=2.38235282,$$

genau der Kontrollwert, den die Implementierungen verwenden.

Wie der Code arbeitet

Die C++-, Python- und Java-Implementierungen erzeugen die rationalen Glieder nur bis \(n=60\). Zunächst werden Zähler- und Nennerfolgen mit exakter Ganzzahlarithmetik aus den Rekurrenzen aufgebaut; erst beim Bilden jedes Quotienten \(s_n\) für die Korrektursumme wird in Gleitkommaarithmetik gewechselt. So bleibt die kurze Vorberechnung numerisch stabil.

Anschließend werden zwei Digamma-Werte berechnet: bei \(N+1-c\) und bei \(2-c\). Dazu nutzt die Implementierung zuerst wiederholt die Rekursion \(\psi(x)=\psi(x+1)-1/x\), bis das Argument groß genug ist. Danach wird die asymptotische Entwicklung

$$\psi(x)=\log x-\frac{1}{2x}-\frac{1}{12x^2}+\frac{1}{120x^4}-\frac{1}{252x^6}+\frac{1}{240x^8}-\frac{1}{132x^{10}}+\frac{691}{32760x^{12}}+\cdots$$

ausgewertet. Abschließend wird

$$\psi(N+1-c)-\psi(2-c)+E_{60}$$

kombiniert und mit \(8\) Nachkommastellen ausgegeben.

Komplexitätsanalyse

Bezeichnet \(K\) die Abschneidegrenze, dann kosten das Erzeugen der Rekurrenzwerte und das Summieren der Korrektur \(O(K)\) Zeit und \(O(K)\) Speicher. Die beiden Digamma-Auswertungen benötigen nur konstante Arbeit und eine feste Anzahl von Termen der asymptotischen Reihe. Da hier \(K=60\) fest gewählt wird, ist die gesamte Berechnung praktisch konstant, obwohl \(N=10^{14}\) riesig ist.

Fußnoten und Quellen

  1. Problemseite: https://projecteuler.net/problem=825
  2. Digamma-Funktion: Wikipedia — Digamma function
  3. Lineare Rekurrenzen mit konstanten Koeffizienten: Wikipedia — Linear recurrence with constant coefficients
  4. Asymptotische Entwicklung: Wikipedia — Asymptotic expansion
  5. Bernoulli-Zahlen in Spezialfunktionsreihen: Wikipedia — Bernoulli number

Problem Özeti

Problem 825,

$$T(N)=\sum_{n=2}^{N} s_n,\qquad N=10^{14}$$

toplamını istemektedir. Chasing Game’den gelen \(s_n\) rasyonel dizisi iki lineer bağıntıyla üretilir. \(10^{14}-1\) terimi tek tek toplamak mümkün olmadığı için çözüm, dizinin baskın asimptotik terimini tam olarak ayırır, bu ana terimi kapalı formda toplar ve geriye yalnızca kısa, çok hızlı yakınsayan bir düzeltme bırakır.

Matematiksel Yaklaşım

Şöyle yazalım:

$$s_n=\frac{a_n}{b_n},$$

başlangıç değerleri

$$a_2=7,\quad a_3=35,\quad a_4=121,$$

$$b_2=11,\quad b_3=73,\quad b_4=395,\quad b_5=1933,$$

ve bağıntılar

$$a_n=3a_{n-1}+3a_{n-2}-a_{n-3}\qquad (n\ge 5),$$

$$b_n=8b_{n-1}-18b_{n-2}+8b_{n-3}-b_{n-4}\qquad (n\ge 6).$$

Uygulamanın bütün fikri bu iki reküransın cebirsel yapısından gelir.

Adım 1: Karakteristik polinomları çarpanlara ayır

Pay için karakteristik polinom

$$r^3-3r^2-3r+1=(r+1)(r^2-4r+1)$$

olur. Payda için ise

$$r^4-8r^3+18r^2-8r+1=(r^2-4r+1)^2$$

elde edilir. Temel kökler bu yüzden

$$\lambda=2+\sqrt3,\qquad \mu=2-\sqrt3=\lambda^{-1}$$

şeklindedir. Özellikle paydada çift katlı baskın çarpan bulunması, daha sonra \(n\) ile çarpılan bir ana terim oluşturur; oranı da kaydırılmış harmonik biçime dönüştüren tam olarak budur.

Adım 2: Reküransları kapalı biçimde çöz

Başlangıç değerleri yerleştirildiğinde

$$A=\frac{3-\sqrt3}{2},\qquad B=\frac{3+\sqrt3}{2}$$

ve

$$a_n=A\lambda^n+B\mu^n-2(-1)^n,$$

$$b_n=\left(An-\frac12\right)\lambda^n+\left(Bn-\frac12\right)\mu^n$$

elde edilir. Bu ifadeler hem verilen başlangıç değerlerini hem de reküransları sağlar. Ayrıca asimptotik davranışı açıkça gösterir: \(\lambda>1\) baskındır, \(\mu<1\) içeren terimler ise üstel olarak küçülür.

Adım 3: Kaydırılmış harmonik ana terimi ayıkla

Payı ve paydayı \(\lambda^n\)’ye bölersek

$$s_n=\frac{A+B\mu^{2n}-2(-1)^n\lambda^{-n}}{An-\frac12+\left(Bn-\frac12\right)\mu^{2n}}$$

olur. Çünkü \(\mu^{2n}=\lambda^{-2n}\), baskın olmayan bütün parçalar üstel hızla kaybolur. Dolayısıyla

$$s_n=\frac{A+O(\lambda^{-n})}{An-\frac12+O(n\lambda^{-2n})}.$$

Şimdi

$$\frac{1}{2A}=\frac{3+\sqrt3}{6}$$

olduğundan

$$c=\frac{3+\sqrt3}{6}$$

tanımını yaparız. Böylece

$$An-\frac12=A(n-c)$$

ve sonuç olarak

$$s_n=\frac{1}{n-c}+O(\lambda^{-n})$$

elde edilir. Yani düzeltme yalnızca polinomik değil, üstel derecede küçüktür.

Adım 4: Dev toplamı kapalı form ve kısa düzeltme olarak ayır

Şimdi

$$T(N)=\sum_{n=2}^{N}\frac{1}{n-c}+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right)$$

yazarız. İkinci toplam, terimleri üstel küçük olduğu için çok hızlı yakınsar. Küçük bir \(K\) sınırında kesersek

$$E_K=\sum_{n=2}^{K}\left(s_n-\frac{1}{n-c}\right)$$

çok iyi bir yaklaşım verir; kalan kuyruk \(8\) ondalık basamak için ihmal edilebilir. Uygulamalar

$$K=60$$

değerini seçer.

Adım 5: Ana terimi digamma özdeşliğiyle topla

Digamma fonksiyonu

$$\psi(x+1)-\psi(x)=\frac1x$$

özelliğini sağlar. Bu yüzden kaydırılmış harmonik kısım için

$$\sum_{n=2}^{N}\frac{1}{n-c}=\psi(N+1-c)-\psi(2-c)$$

olur. Dolayısıyla

$$T(N)=\psi(N+1-c)-\psi(2-c)+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right),$$

ve sayısal olarak

$$T(N)\approx \psi(N+1-c)-\psi(2-c)+E_{60}$$

yeterlidir.

Çözümlü Örnek: kontrol değeri \(T(10)\)

İlk birkaç terim

$$s_2=\frac{7}{11}\approx 0.63636364,\qquad s_3=\frac{35}{73}\approx 0.47945205,\qquad s_4=\frac{121}{395}\approx 0.30632911$$

şeklindedir. Ayrıca

$$c=\frac{3+\sqrt3}{6}\approx 0.78867513$$

olduğundan, \(10\)’a kadar kaydırılmış harmonik toplam

$$\sum_{n=2}^{10}\frac{1}{n-c}\approx 2.54851473$$

ve aynı aralıktaki düzeltme

$$\sum_{n=2}^{10}\left(s_n-\frac{1}{n-c}\right)\approx -0.16616191$$

olur. Böylece

$$T(10)\approx 2.54851473-0.16616191=2.38235282,$$

ki bu da uygulamalardaki sayısal kontrol değeriyle tam uyumludur.

Kod Nasıl Çalışır

C++, Python ve Java implementasyonları rasyonel terimleri yalnızca \(n=60\)’a kadar üretir. Önce pay ve payda dizileri reküranslardan tam sayı aritmetiğiyle oluşturulur; düzeltme toplamındaki her \(s_n\) oranı alınırken kayan nokta aritmetiğine geçilir. Böylece kısa ön hesaplama hem kararlı hem de hassas kalır.

Daha sonra iki digamma değeri hesaplanır: biri \(N+1-c\) noktasında, diğeri \(2-c\) noktasında. Bunun için uygulama önce \(\psi(x)=\psi(x+1)-1/x\) bağıntısını tekrar tekrar kullanarak girdiyi yeterince büyük yapar. Ardından

$$\psi(x)=\log x-\frac{1}{2x}-\frac{1}{12x^2}+\frac{1}{120x^4}-\frac{1}{252x^6}+\frac{1}{240x^8}-\frac{1}{132x^{10}}+\frac{691}{32760x^{12}}+\cdots$$

asimptotik serisini uygular. Son adımda

$$\psi(N+1-c)-\psi(2-c)+E_{60}$$

birleştirilir ve sonuç \(8\) ondalık basamakla yazdırılır.

Karmaşıklık Analizi

Kesme noktası \(K\) ile gösterilirse, gerekli rekürans değerlerini üretmek ve düzeltmeyi toplamak \(O(K)\) zaman ve \(O(K)\) bellek gerektirir. İki digamma hesabı yalnızca sabit sayıda aritmetik işlem ve sabit uzunlukta bir asimptotik seri ister. Uygulamada \(K=60\) sabit seçildiği için, \(N=10^{14}\) gibi dev bir değer olsa da toplam maliyet pratikte sabit zamandır.

Dipnotlar ve Kaynakça

  1. Problem sayfası: https://projecteuler.net/problem=825
  2. Digamma fonksiyonu: Wikipedia — Digamma function
  3. Sabit katsayılı lineer reküranslar: Wikipedia — Linear recurrence with constant coefficients
  4. Asimptotik açılım: Wikipedia — Asymptotic expansion
  5. Özel fonksiyon açılımlarında Bernoulli sayıları: Wikipedia — Bernoulli number

Resumen del Problema

El problema 825 pide la suma parcial

$$T(N)=\sum_{n=2}^{N} s_n,\qquad N=10^{14},$$

donde la sucesión racional \(s_n\) procede del Chasing Game y se genera mediante dos recurrencias lineales. Sumar directamente \(10^{14}-1\) términos es imposible, así que la solución separa el término asintótico dominante exacto, suma ese término en forma cerrada y deja solo una corrección corta que converge muy deprisa.

Enfoque Matemático

Escribimos

$$s_n=\frac{a_n}{b_n},$$

con valores iniciales

$$a_2=7,\quad a_3=35,\quad a_4=121,$$

$$b_2=11,\quad b_3=73,\quad b_4=395,\quad b_5=1933,$$

y recurrencias

$$a_n=3a_{n-1}+3a_{n-2}-a_{n-3}\qquad (n\ge 5),$$

$$b_n=8b_{n-1}-18b_{n-2}+8b_{n-3}-b_{n-4}\qquad (n\ge 6).$$

La solución aprovecha la estructura exacta escondida en estas dos relaciones.

Paso 1: Factorizar los polinomios característicos

La recurrencia del numerador tiene polinomio característico

$$r^3-3r^2-3r+1=(r+1)(r^2-4r+1).$$

La del denominador tiene

$$r^4-8r^3+18r^2-8r+1=(r^2-4r+1)^2.$$

Por tanto, las raíces clave son

$$\lambda=2+\sqrt3,\qquad \mu=2-\sqrt3=\lambda^{-1}.$$

El hecho de que el denominador tenga un factor dominante repetido explica por qué aparece un factor lineal en \(n\), y eso es precisamente lo que transforma el cociente en un término armónico desplazado.

Paso 2: Resolver las recurrencias en forma cerrada

A partir de los valores iniciales se obtiene

$$A=\frac{3-\sqrt3}{2},\qquad B=\frac{3+\sqrt3}{2},$$

$$a_n=A\lambda^n+B\mu^n-2(-1)^n,$$

$$b_n=\left(An-\frac12\right)\lambda^n+\left(Bn-\frac12\right)\mu^n.$$

Estas fórmulas satisfacen tanto las recurrencias como los datos iniciales. Además dejan claro el comportamiento asintótico: \(\lambda>1\) domina, mientras que los términos con \(\mu<1\) son exponencialmente pequeños.

Paso 3: Extraer el término armónico desplazado

Si dividimos numerador y denominador por \(\lambda^n\), obtenemos

$$s_n=\frac{A+B\mu^{2n}-2(-1)^n\lambda^{-n}}{An-\frac12+\left(Bn-\frac12\right)\mu^{2n}}.$$

Como \(\mu^{2n}=\lambda^{-2n}\), todas las partes no dominantes decaen exponencialmente. Entonces

$$s_n=\frac{A+O(\lambda^{-n})}{An-\frac12+O(n\lambda^{-2n})}.$$

Además,

$$\frac{1}{2A}=\frac{3+\sqrt3}{6},$$

de modo que definimos

$$c=\frac{3+\sqrt3}{6}.$$

Con eso,

$$An-\frac12=A(n-c),$$

y por tanto

$$s_n=\frac{1}{n-c}+O(\lambda^{-n}).$$

La corrección es exponencialmente pequeña, no solo de orden polinómico.

Paso 4: Separar la suma enorme en una parte cerrada y una corrección corta

Descomponemos

$$T(N)=\sum_{n=2}^{N}\frac{1}{n-c}+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right).$$

La segunda suma converge muy rápido. Si la truncamos en un valor modesto \(K\), definimos

$$E_K=\sum_{n=2}^{K}\left(s_n-\frac{1}{n-c}\right),$$

y la cola omitida ya es despreciable para una salida de \(8\) decimales. Las implementaciones toman

$$K=60.$$

Paso 5: Sumar la parte principal con la identidad de la digamma

La función digamma cumple

$$\psi(x+1)-\psi(x)=\frac1x.$$

Aplicando esta identidad al término armónico desplazado se obtiene

$$\sum_{n=2}^{N}\frac{1}{n-c}=\psi(N+1-c)-\psi(2-c).$$

Así,

$$T(N)=\psi(N+1-c)-\psi(2-c)+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right),$$

y numéricamente basta usar

$$T(N)\approx \psi(N+1-c)-\psi(2-c)+E_{60}.$$

Ejemplo trabajado: el punto de control \(T(10)\)

Los primeros términos son

$$s_2=\frac{7}{11}\approx 0.63636364,\qquad s_3=\frac{35}{73}\approx 0.47945205,\qquad s_4=\frac{121}{395}\approx 0.30632911.$$

Con

$$c=\frac{3+\sqrt3}{6}\approx 0.78867513,$$

la parte armónica desplazada hasta \(10\) vale

$$\sum_{n=2}^{10}\frac{1}{n-c}\approx 2.54851473.$$

La corrección sobre el mismo rango es

$$\sum_{n=2}^{10}\left(s_n-\frac{1}{n-c}\right)\approx -0.16616191.$$

Por tanto,

$$T(10)\approx 2.54851473-0.16616191=2.38235282,$$

exactamente el valor de comprobación numérica usado por las implementaciones.

Cómo Funciona el Código

Las implementaciones en C++, Python y Java generan los términos racionales solo hasta \(n=60\). Primero construyen las sucesiones del numerador y del denominador con aritmética entera exacta; solo convierten a coma flotante al formar cada cociente \(s_n\) dentro de la suma de corrección. Así la precomputación corta permanece estable y precisa.

Después se evalúan dos valores de la función digamma: uno en \(N+1-c\) y otro en \(2-c\). Para ello, la implementación usa primero repetidamente la relación \(\psi(x)=\psi(x+1)-1/x\) hasta que el argumento es suficientemente grande. Luego aplica la expansión asintótica

$$\psi(x)=\log x-\frac{1}{2x}-\frac{1}{12x^2}+\frac{1}{120x^4}-\frac{1}{252x^6}+\frac{1}{240x^8}-\frac{1}{132x^{10}}+\frac{691}{32760x^{12}}+\cdots.$$

Al final combina

$$\psi(N+1-c)-\psi(2-c)+E_{60}$$

y muestra el resultado con \(8\) decimales.

Análisis de Complejidad

Si denotamos por \(K\) el punto de corte, generar todos los valores de las recurrencias y sumar la corrección cuesta \(O(K)\) tiempo y \(O(K)\) memoria. Las dos evaluaciones de digamma usan solo trabajo constante y un número fijo de términos de la serie asintótica. Como la implementación fija \(K=60\), todo el procedimiento es efectivamente de tiempo constante aunque \(N=10^{14}\) sea gigantesco.

Notas y Referencias

  1. Página del problema: https://projecteuler.net/problem=825
  2. Función digamma: Wikipedia — Digamma function
  3. Recurrencias lineales con coeficientes constantes: Wikipedia — Linear recurrence with constant coefficients
  4. Expansión asintótica: Wikipedia — Asymptotic expansion
  5. Números de Bernoulli en expansiones de funciones especiales: Wikipedia — Bernoulli number

问题概述

第 825 题要求计算部分和

$$T(N)=\sum_{n=2}^{N} s_n,\qquad N=10^{14},$$

其中有理数列 \(s_n\) 来自 Chasing Game,并由两个线性递推定义。对 \(10^{14}-1\) 个项直接求和显然不可行,所以正确做法是先找出这个数列的精确主导渐近项,把主项用闭式求和,再只保留一段很短但收敛极快的修正和。

数学方法

$$s_n=\frac{a_n}{b_n},$$

其初值为

$$a_2=7,\quad a_3=35,\quad a_4=121,$$

$$b_2=11,\quad b_3=73,\quad b_4=395,\quad b_5=1933,$$

递推关系为

$$a_n=3a_{n-1}+3a_{n-2}-a_{n-3}\qquad (n\ge 5),$$

$$b_n=8b_{n-1}-18b_{n-2}+8b_{n-3}-b_{n-4}\qquad (n\ge 6).$$

整个解法的关键,就是把这两个递推背后的代数结构完全展开。

步骤 1:分解特征多项式

分子递推的特征多项式是

$$r^3-3r^2-3r+1=(r+1)(r^2-4r+1).$$

分母递推的特征多项式是

$$r^4-8r^3+18r^2-8r+1=(r^2-4r+1)^2.$$

因此最重要的两个根为

$$\lambda=2+\sqrt3,\qquad \mu=2-\sqrt3=\lambda^{-1}.$$

这里最值得注意的是:分母的主导因子是重根,所以分母的主项不是简单的常数乘 \(\lambda^n\),而是会多出一个线性因子 \(n\)。正是这个结构,使得 \(a_n/b_n\) 的主项自然变成“带位移的调和项”。

步骤 2:写出递推的闭式表达

把初值代入求解,可以得到

$$A=\frac{3-\sqrt3}{2},\qquad B=\frac{3+\sqrt3}{2},$$

$$a_n=A\lambda^n+B\mu^n-2(-1)^n,$$

$$b_n=\left(An-\frac12\right)\lambda^n+\left(Bn-\frac12\right)\mu^n.$$

这两条公式同时满足给定的初值和递推。更重要的是,它们直接暴露了渐近行为:因为 \(\lambda>1\),\(\lambda^n\) 对应的部分会主导增长;因为 \(\mu<1\),涉及 \(\mu^n\) 的项会指数级衰减;额外的 \((-1)^n\) 项在除以 \(\lambda^n\) 后也会迅速消失。

步骤 3:提取主项 \(\frac{1}{n-c}\)

把分子和分母同时除以 \(\lambda^n\),得到

$$s_n=\frac{A+B\mu^{2n}-2(-1)^n\lambda^{-n}}{An-\frac12+\left(Bn-\frac12\right)\mu^{2n}}.$$

由于 \(\mu=\lambda^{-1}\),所以 \(\mu^{2n}=\lambda^{-2n}\)。这说明除主项以外的所有部分都在指数级缩小,因此

$$s_n=\frac{A+O(\lambda^{-n})}{An-\frac12+O(n\lambda^{-2n})}.$$

再注意到

$$\frac{1}{2A}=\frac{3+\sqrt3}{6},$$

于是定义

$$c=\frac{3+\sqrt3}{6}.$$

这样就有

$$An-\frac12=A(n-c),$$

从而得到

$$s_n=\frac{1}{n-c}+O(\lambda^{-n}).$$

这里最核心的结论是:减掉 \(\frac{1}{n-c}\) 以后,剩下的误差不是慢慢按幂次衰减,而是按指数速度衰减,所以只需要很短的一段修正和就足够了。

步骤 4:把超长求和拆成闭式部分和快速收敛修正

因此可以写成

$$T(N)=\sum_{n=2}^{N}\frac{1}{n-c}+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right).$$

第二部分的求和项指数级变小,所以收敛非常快。把它截断到一个不大的 \(K\),定义

$$E_K=\sum_{n=2}^{K}\left(s_n-\frac{1}{n-c}\right),$$

则后面的尾项对最终 \(8\) 位小数输出已经几乎没有影响。实际实现中取的是

$$K=60.$$

步骤 5:用 digamma 恒等式求和主项

digamma 函数满足

$$\psi(x+1)-\psi(x)=\frac1x.$$

把这条恒等式用于移位调和和,就得到

$$\sum_{n=2}^{N}\frac{1}{n-c}=\psi(N+1-c)-\psi(2-c).$$

所以目标值可以写成

$$T(N)=\psi(N+1-c)-\psi(2-c)+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right).$$

在数值计算时,只需使用

$$T(N)\approx \psi(N+1-c)-\psi(2-c)+E_{60}$$

即可达到所需精度。

算例:检查点 \(T(10)\)

前几项为

$$s_2=\frac{7}{11}\approx 0.63636364,\qquad s_3=\frac{35}{73}\approx 0.47945205,\qquad s_4=\frac{121}{395}\approx 0.30632911.$$

同时

$$c=\frac{3+\sqrt3}{6}\approx 0.78867513.$$

那么到 \(10\) 为止的主项和为

$$\sum_{n=2}^{10}\frac{1}{n-c}\approx 2.54851473,$$

同一区间上的修正和为

$$\sum_{n=2}^{10}\left(s_n-\frac{1}{n-c}\right)\approx -0.16616191.$$

两者相加得到

$$T(10)\approx 2.54851473-0.16616191=2.38235282,$$

这正好与实现中使用的数值检查点一致。

代码如何实现

C++、Python 和 Java 实现都只把有理项显式生成到 \(n=60\)。它们先用精确整数算术根据递推式构造分子和分母,再在累加修正和时把每一项 \(s_n\) 转成浮点数。这种做法既保留了递推阶段的准确性,也避免了无意义的大规模逐项求和。

接下来只需要计算两个 digamma 值:\(\psi(N+1-c)\) 和 \(\psi(2-c)\)。实现采用两步法。第一步,反复使用递推关系 \(\psi(x)=\psi(x+1)-1/x\),把自变量移动到足够大的区域。第二步,在大参数区间使用渐近展开

$$\psi(x)=\log x-\frac{1}{2x}-\frac{1}{12x^2}+\frac{1}{120x^4}-\frac{1}{252x^6}+\frac{1}{240x^8}-\frac{1}{132x^{10}}+\frac{691}{32760x^{12}}+\cdots.$$

最后组合

$$\psi(N+1-c)-\psi(2-c)+E_{60}$$

并输出 \(8\) 位小数结果。

复杂度分析

如果把截断点记为 \(K\),那么生成递推值并累加修正和的复杂度是 \(O(K)\) 时间、\(O(K)\) 空间。两次 digamma 计算只需要常数次算术操作以及固定项数的渐近级数。由于实现中 \(K=60\) 是固定常数,所以即使 \(N=10^{14}\) 非常大,整体计算量仍然近似常数。

参考资料

  1. 题目页面:https://projecteuler.net/problem=825
  2. Digamma 函数:Wikipedia — Digamma function
  3. 常系数线性递推:Wikipedia — Linear recurrence with constant coefficients
  4. 渐近展开:Wikipedia — Asymptotic expansion
  5. 特殊函数展开中的 Bernoulli 数:Wikipedia — Bernoulli number

Краткое описание задачи

В задаче 825 требуется вычислить частичную сумму

$$T(N)=\sum_{n=2}^{N} s_n,\qquad N=10^{14},$$

где рациональная последовательность \(s_n\) возникает в Chasing Game и задается двумя линейными рекуррентными соотношениями. Напрямую складывать \(10^{14}-1\) слагаемых невозможно, поэтому решение выделяет точный главный асимптотический член, суммирует его в закрытом виде и оставляет лишь короткую быстро сходящуюся поправку.

Математический подход

Запишем

$$s_n=\frac{a_n}{b_n},$$

где начальные значения равны

$$a_2=7,\quad a_3=35,\quad a_4=121,$$

$$b_2=11,\quad b_3=73,\quad b_4=395,\quad b_5=1933,$$

а рекуррентные формулы имеют вид

$$a_n=3a_{n-1}+3a_{n-2}-a_{n-3}\qquad (n\ge 5),$$

$$b_n=8b_{n-1}-18b_{n-2}+8b_{n-3}-b_{n-4}\qquad (n\ge 6).$$

Вся математика решения строится на точном анализе этих двух рекуррентных последовательностей.

Шаг 1: Разложить характеристические многочлены

Для числителя характеристический многочлен равен

$$r^3-3r^2-3r+1=(r+1)(r^2-4r+1).$$

Для знаменателя получаем

$$r^4-8r^3+18r^2-8r+1=(r^2-4r+1)^2.$$

Следовательно, ключевые корни таковы:

$$\lambda=2+\sqrt3,\qquad \mu=2-\sqrt3=\lambda^{-1}.$$

У знаменателя доминирующий множитель кратный, поэтому в ведущем члене появляется дополнительный множитель \(n\). Именно это и превращает отношение \(a_n/b_n\) в сдвинутый гармонический член.

Шаг 2: Выписать замкнутые формулы

Из начальных условий получаются константы

$$A=\frac{3-\sqrt3}{2},\qquad B=\frac{3+\sqrt3}{2},$$

и выражения

$$a_n=A\lambda^n+B\mu^n-2(-1)^n,$$

$$b_n=\left(An-\frac12\right)\lambda^n+\left(Bn-\frac12\right)\mu^n.$$

Эти формулы удовлетворяют и рекуррентным соотношениям, и заданным стартовым значениям. Кроме того, из них сразу видна асимптотика: часть с \(\lambda>1\) доминирует, а слагаемые с \(\mu<1\) убывают экспоненциально.

Шаг 3: Выделить главный член \(\frac{1}{n-c}\)

Разделим числитель и знаменатель на \(\lambda^n\):

$$s_n=\frac{A+B\mu^{2n}-2(-1)^n\lambda^{-n}}{An-\frac12+\left(Bn-\frac12\right)\mu^{2n}}.$$

Так как \(\mu^{2n}=\lambda^{-2n}\), все недоминирующие части исчезают экспоненциально быстро. Поэтому

$$s_n=\frac{A+O(\lambda^{-n})}{An-\frac12+O(n\lambda^{-2n})}.$$

Поскольку

$$\frac{1}{2A}=\frac{3+\sqrt3}{6},$$

удобно ввести

$$c=\frac{3+\sqrt3}{6}.$$

Тогда

$$An-\frac12=A(n-c),$$

и отсюда следует

$$s_n=\frac{1}{n-c}+O(\lambda^{-n}).$$

То есть после вычитания главного члена остается не медленный степенной хвост, а экспоненциально малая поправка.

Шаг 4: Разбить огромную сумму на закрытую часть и короткую поправку

Теперь можно записать

$$T(N)=\sum_{n=2}^{N}\frac{1}{n-c}+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right).$$

Вторая сумма сходится очень быстро. Если обрезать ее на умеренном \(K\), получим

$$E_K=\sum_{n=2}^{K}\left(s_n-\frac{1}{n-c}\right),$$

а оставшийся хвост уже не влияет на ответ с точностью до \(8\) знаков после запятой. В реализации берется

$$K=60.$$

Шаг 5: Просуммировать главный член через дигамма-функцию

Для дигамма-функции верно тождество

$$\psi(x+1)-\psi(x)=\frac1x.$$

Применяя его к сдвинутой гармонической сумме, получаем

$$\sum_{n=2}^{N}\frac{1}{n-c}=\psi(N+1-c)-\psi(2-c).$$

Следовательно,

$$T(N)=\psi(N+1-c)-\psi(2-c)+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right),$$

и численно достаточно использовать

$$T(N)\approx \psi(N+1-c)-\psi(2-c)+E_{60}.$$

Разобранный пример: контрольное значение \(T(10)\)

Первые члены таковы:

$$s_2=\frac{7}{11}\approx 0.63636364,\qquad s_3=\frac{35}{73}\approx 0.47945205,\qquad s_4=\frac{121}{395}\approx 0.30632911.$$

При

$$c=\frac{3+\sqrt3}{6}\approx 0.78867513$$

сдвинутая гармоническая сумма до \(10\) равна

$$\sum_{n=2}^{10}\frac{1}{n-c}\approx 2.54851473.$$

Поправка на том же диапазоне равна

$$\sum_{n=2}^{10}\left(s_n-\frac{1}{n-c}\right)\approx -0.16616191.$$

Отсюда

$$T(10)\approx 2.54851473-0.16616191=2.38235282,$$

что совпадает с численной контрольной точкой в реализациях.

Как Работает Код

Реализации на C++, Python и Java строят рациональные члены только до \(n=60\). Сначала последовательности числителей и знаменателей вычисляются точно в целых числах, а в вещественную арифметику переход осуществляется лишь при формировании каждого отношения \(s_n\) в сумме поправки. Это делает короткую предвычислительную часть устойчивой и точной.

После этого нужно вычислить два значения дигамма-функции: в точках \(N+1-c\) и \(2-c\). Для этого реализация сначала многократно применяет рекуррентную формулу \(\psi(x)=\psi(x+1)-1/x\), пока аргумент не станет достаточно большим. Затем используется асимптотическое разложение

$$\psi(x)=\log x-\frac{1}{2x}-\frac{1}{12x^2}+\frac{1}{120x^4}-\frac{1}{252x^6}+\frac{1}{240x^8}-\frac{1}{132x^{10}}+\frac{691}{32760x^{12}}+\cdots.$$

В конце программа объединяет

$$\psi(N+1-c)-\psi(2-c)+E_{60}$$

и печатает ответ с \(8\) знаками после запятой.

Анализ Сложности

Если обозначить точку обрезания через \(K\), то вычисление значений рекуррентных последовательностей и суммирование поправки занимают \(O(K)\) времени и \(O(K)\) памяти. Две оценки дигамма-функции требуют только константного числа арифметических операций и фиксированного количества членов асимптотического ряда. Поскольку в реализации зафиксировано \(K=60\), вся процедура фактически работает за постоянное время, хотя \(N=10^{14}\) чрезвычайно велико.

Ссылки и Литература

  1. Страница задачи: https://projecteuler.net/problem=825
  2. Дигамма-функция: Wikipedia — Digamma function
  3. Линейные рекуррентные соотношения с постоянными коэффициентами: Wikipedia — Linear recurrence with constant coefficients
  4. Асимптотическое разложение: Wikipedia — Asymptotic expansion
  5. Числа Бернулли в разложениях специальных функций: Wikipedia — Bernoulli number

ملخص المسألة

تطلب المسألة 825 حساب المجموع الجزئي

$$T(N)=\sum_{n=2}^{N} s_n,\qquad N=10^{14},$$

حيث إن المتتالية النسبية \(s_n\) ناشئة من Chasing Game وتُعطى بواسطة علاقتين خطيتين تراجعيّتين. لا يمكن جمع \(10^{14}-1\) حدًا مباشرة، لذلك تعتمد الفكرة الصحيحة على استخراج الحد التقاربي الرئيسي بدقة، وجمعه بصيغة مغلقة، ثم الإبقاء فقط على تصحيح قصير سريع التقارب.

المنهج الرياضي

نكتب

$$s_n=\frac{a_n}{b_n},$$

مع القيم الابتدائية

$$a_2=7,\quad a_3=35,\quad a_4=121,$$

$$b_2=11,\quad b_3=73,\quad b_4=395,\quad b_5=1933,$$

وعلاقتي التراجع

$$a_n=3a_{n-1}+3a_{n-2}-a_{n-3}\qquad (n\ge 5),$$

$$b_n=8b_{n-1}-18b_{n-2}+8b_{n-3}-b_{n-4}\qquad (n\ge 6).$$

جوهر الحل هو تحليل البنية الجبرية الدقيقة الكامنة خلف هاتين العلاقتين.

الخطوة 1: تحليل كثيرات الحدود المميِّزة

لعلاقة البسط يكون كثير الحدود المميز

$$r^3-3r^2-3r+1=(r+1)(r^2-4r+1).$$

أما علاقة المقام فتعطي

$$r^4-8r^3+18r^2-8r+1=(r^2-4r+1)^2.$$

ومن ثم فالجذران الحاسمان هما

$$\lambda=2+\sqrt3,\qquad \mu=2-\sqrt3=\lambda^{-1}.$$

وجود العامل المسيطر مكررًا في المقام هو السبب في ظهور عامل خطي في \(n\)، وهذا بالضبط ما يحوّل النسبة \(a_n/b_n\) إلى حد توافقي مزاح.

الخطوة 2: كتابة الحلين المغلقين

من الشروط الابتدائية نحصل على

$$A=\frac{3-\sqrt3}{2},\qquad B=\frac{3+\sqrt3}{2},$$

ثم

$$a_n=A\lambda^n+B\mu^n-2(-1)^n,$$

$$b_n=\left(An-\frac12\right)\lambda^n+\left(Bn-\frac12\right)\mu^n.$$

هذه الصيغ تحقق علاقات التراجع والقيم الابتدائية معًا. وهي أيضًا تكشف السلوك التقاربي مباشرة: الجزء الذي يحوي \(\lambda>1\) هو المسيطر، أما الحدود التي تحوي \(\mu<1\) فتتناقص أسيًا.

الخطوة 3: استخراج الحد الرئيسي \(\frac{1}{n-c}\)

عند قسمة البسط والمقام على \(\lambda^n\) نحصل على

$$s_n=\frac{A+B\mu^{2n}-2(-1)^n\lambda^{-n}}{An-\frac12+\left(Bn-\frac12\right)\mu^{2n}}.$$

وبما أن \(\mu^{2n}=\lambda^{-2n}\)، فإن جميع الأجزاء غير المسيطرة تضمحل أسيًا، ومن ثم

$$s_n=\frac{A+O(\lambda^{-n})}{An-\frac12+O(n\lambda^{-2n})}.$$

ولأن

$$\frac{1}{2A}=\frac{3+\sqrt3}{6},$$

نعرّف

$$c=\frac{3+\sqrt3}{6}.$$

وعندئذٍ

$$An-\frac12=A(n-c),$$

فتصبح الصيغة

$$s_n=\frac{1}{n-c}+O(\lambda^{-n}).$$

إذن الفرق بعد طرح الحد الرئيسي ليس مجرد ذيل كثير حدود، بل تصحيح أسي الصغر.

الخطوة 4: تفكيك المجموع الكبير إلى جزء مغلق وتصحيح قصير

نكتب الآن

$$T(N)=\sum_{n=2}^{N}\frac{1}{n-c}+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right).$$

المجموع الثاني يتقارب بسرعة كبيرة لأن حدوده صغيرة أسيًا. وإذا قطعناه عند قيمة معتدلة \(K\)، نعرّف

$$E_K=\sum_{n=2}^{K}\left(s_n-\frac{1}{n-c}\right),$$

فتصبح بقية الذيل مهملة على مستوى الدقة المطلوبة ذي \(8\) منازل عشرية. التنفيذات تستخدم

$$K=60.$$

الخطوة 5: جمع الجزء الرئيسي بهوية دالة digamma

تحقق دالة digamma الهوية

$$\psi(x+1)-\psi(x)=\frac1x.$$

وبتطبيقها على الجزء التوافقي المزاح نحصل على

$$\sum_{n=2}^{N}\frac{1}{n-c}=\psi(N+1-c)-\psi(2-c).$$

إذن القيمة المطلوبة هي

$$T(N)=\psi(N+1-c)-\psi(2-c)+\sum_{n=2}^{N}\left(s_n-\frac{1}{n-c}\right),$$

ومن الناحية العددية يكفي استعمال

$$T(N)\approx \psi(N+1-c)-\psi(2-c)+E_{60}.$$

مثال محلول: نقطة التحقق \(T(10)\)

الحدود الأولى هي

$$s_2=\frac{7}{11}\approx 0.63636364,\qquad s_3=\frac{35}{73}\approx 0.47945205,\qquad s_4=\frac{121}{395}\approx 0.30632911.$$

ومع

$$c=\frac{3+\sqrt3}{6}\approx 0.78867513$$

يكون الجزء التوافقي المزاح حتى \(10\) مساويًا لـ

$$\sum_{n=2}^{10}\frac{1}{n-c}\approx 2.54851473.$$

أما التصحيح على المدى نفسه فهو

$$\sum_{n=2}^{10}\left(s_n-\frac{1}{n-c}\right)\approx -0.16616191.$$

ومن ثم

$$T(10)\approx 2.54851473-0.16616191=2.38235282,$$

وهو نفس المقدار العددي الذي تُثبّت به التنفيذات صحة الحساب.

كيف تعمل الشيفرة

تنفيذات C++ وPython وJava تولد الحدود النسبية فقط حتى \(n=60\). أولًا تُبنى متتاليتا البسط والمقام بدقة صحيحة كاملة، ثم تتحول كل نسبة \(s_n\) إلى عدد عشري فقط عند إضافتها إلى مجموع التصحيح. بهذه الطريقة تبقى مرحلة التمهيد القصيرة دقيقة ومستقرة عدديًا.

بعد ذلك لا يبقى إلا تقييم قيمتين لدالة digamma: عند \(N+1-c\) وعند \(2-c\). ويُنجز ذلك على مرحلتين. في المرحلة الأولى تُستخدم العلاقة \(\psi(x)=\psi(x+1)-1/x\) مرارًا حتى يصبح الوسيط كبيرًا بما يكفي. وفي المرحلة الثانية تُطبَّق المتسلسلة التقاربية

$$\psi(x)=\log x-\frac{1}{2x}-\frac{1}{12x^2}+\frac{1}{120x^4}-\frac{1}{252x^6}+\frac{1}{240x^8}-\frac{1}{132x^{10}}+\frac{691}{32760x^{12}}+\cdots.$$

وأخيرًا تُجمَع الكمية

$$\psi(N+1-c)-\psi(2-c)+E_{60}$$

ويُطبع الناتج حتى \(8\) منازل عشرية.

تحليل التعقيد

إذا رمزنا إلى نقطة القطع بـ \(K\)، فإن توليد قيم علاقات التراجع وجمع التصحيح يحتاجان إلى زمن \(O(K)\) وذاكرة \(O(K)\). أما تقييمتا digamma فتحتاجان إلى عمل ثابت وعدد ثابت من حدود المتسلسلة التقاربية. وبما أن التنفيذ يثبّت \(K=60\)، فإن الكلفة الإجمالية تصبح عمليًا ثابتة حتى عندما يكون \(N=10^{14}\) هائلًا.

مراجع وحواشٍ

  1. صفحة المسألة: https://projecteuler.net/problem=825
  2. دالة digamma: Wikipedia — Digamma function
  3. العلاقات التراجعية الخطية ذات المعاملات الثابتة: Wikipedia — Linear recurrence with constant coefficients
  4. التوسع التقاربي: Wikipedia — Asymptotic expansion
  5. أعداد برنولي في توسعات الدوال الخاصة: Wikipedia — Bernoulli number