Problem Summary

For this problem, the drone-delivery expectation reduces to evaluating

$$E(n)=n\sum_{1\le i\le j\le n}\frac{1}{ij}=\frac n2\left(H_n^2+H_n^{(2)}\right),$$

at the target value \(n=10^8\), where \(H_n=\sum_{k=1}^n \frac1k\) and \(H_n^{(2)}=\sum_{k=1}^n \frac1{k^2}\). The implementations exploit this exact identity and then switch to Euler-Maclaurin expansions, because summing \(10^8\) terms directly is unnecessary.

Mathematical Approach

Once the underlying probability model is simplified, everything depends only on harmonic sums. That replaces a large discrete process with a short analytic formula.

Step 1: Write the exact closed form

Introduce the ordinary harmonic number

$$H_n=\sum_{k=1}^n\frac1k$$

and the second-order harmonic number

$$H_n^{(2)}=\sum_{k=1}^n\frac1{k^2}.$$

The reduced expectation is then

$$E(n)=\frac n2\left(H_n^2+H_n^{(2)}\right).$$

This identity is exact, not asymptotic. For moderate \(n\), one may therefore compute \(E(n)\) directly from the defining sums.

Step 2: See why the double sum becomes harmonic numbers

Expanding the square gives

$$H_n^2=\left(\sum_{i=1}^n\frac1i\right)\left(\sum_{j=1}^n\frac1j\right)=\sum_{i=1}^n\frac1{i^2}+2\sum_{1\le i<j\le n}\frac1{ij}.$$

Adding \(H_n^{(2)}\) duplicates the diagonal terms, so

$$H_n^2+H_n^{(2)}=2\sum_{1\le i\le j\le n}\frac1{ij}.$$

Hence

$$E(n)=n\sum_{1\le i\le j\le n}\frac1{ij}.$$

This triangular-sum viewpoint explains why the closed form contains both \(H_n^2\) and \(H_n^{(2)}\): the off-diagonal pairs come from the square, while the diagonal \(i=j\) contributes the second-order harmonic sum.

Step 3: Use Euler-Maclaurin for large \(n\)

For the target input, direct summation is wasteful. The implementations use the standard asymptotic expansions

$$H_n=\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}+O\left(n^{-10}\right),$$

$$H_n^{(2)}=\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}+O\left(n^{-11}\right).$$

Here \(\gamma\) is the Euler-Mascheroni constant, and \(\pi^2/6\) is the limit of the second-order harmonic sum.

Step 4: Substitute the expansions into the exact formula

Plugging those series into the identity for \(E(n)\) yields an approximation that is effectively exact for \(n=10^8\):

$$E(n)\approx \frac n2\left[\left(\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}\right)^2+\left(\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}\right)\right].$$

The dominant growth is

$$E(n)=\frac n2\left((\log n+\gamma)^2+\frac{\pi^2}{6}\right)+O(\log n).$$

So the quantity grows on the order of \(n\log^2 n\), and the omitted tail of the asymptotic series is far too small to affect the nearest integer when \(n=10^8\).

Step 5: Why rounding to the nearest integer is safe

After multiplication by \(n/2\), the first neglected terms are still tiny: the harmonic-number remainder contributes on the order of \(n^{-9}\), and the second-order remainder contributes on the order of \(n^{-10}\). At \(n=10^8\), those analytic errors are astronomically below \(1/2\).

That is why the implementations can evaluate the truncated series with high precision and then round the result to the nearest integer.

Worked Example: \(n=5\)

For \(n=5\), direct summation is easy:

$$H_5=1+\frac12+\frac13+\frac14+\frac15=\frac{137}{60},$$

$$H_5^{(2)}=1+\frac14+\frac19+\frac1{16}+\frac1{25}=\frac{5269}{3600}.$$

Therefore

$$E(5)=\frac52\left[\left(\frac{137}{60}\right)^2+\frac{5269}{3600}\right]=\frac{12019}{720}\approx 16.6930555556.$$

This matches the exact checkpoint used by the implementations and confirms the closed form before any asymptotic approximation is applied.

How the Code Works

The C++, Python, and Java implementations all evaluate the same formula. They first prepare the constants \(\gamma\), \(\pi\), and the target input \(n=10^8\). Then they build the inverse powers \(n^{-1},n^{-2},\dots,n^{-9}\), because the Euler-Maclaurin expansions reuse them many times.

The C++ implementation also checks a few small values by direct summation, which verifies the exact harmonic identity numerically before the large-\(n\) path is used. For the target computation, the implementations form the truncated series for \(H_n\) and \(H_n^{(2)}\), substitute them into \(\frac n2(H_n^2+H_n^{(2)})\), and round to the nearest integer.

The Python and Java implementations use 50-digit decimal arithmetic. The Java version stores the natural logarithm of the target input as a high-precision constant, while the C++ version computes the logarithm directly in extended floating-point arithmetic.

Complexity Analysis

If one summed the defining series directly, the cost would be \(O(n)\) time and \(O(1)\) extra memory. The asymptotic method used for the real input evaluates only a fixed number of arithmetic operations on high-precision numbers, so with fixed precision it is \(O(1)\) time and \(O(1)\) memory. That is the point of replacing the huge sums with Euler-Maclaurin expansions.

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=724
  2. Harmonic number: Wikipedia - Harmonic number
  3. Generalized harmonic number: Wikipedia - Generalized harmonic number
  4. Euler-Maclaurin formula: Wikipedia - Euler-Maclaurin formula
  5. Euler's constant: Wikipedia - Euler's constant
  6. Basel problem: Wikipedia - Basel problem

Problemzusammenfassung

Für dieses Problem reduziert sich der Erwartungswert des Zustellmodells auf die Auswertung von

$$E(n)=n\sum_{1\le i\le j\le n}\frac{1}{ij}=\frac n2\left(H_n^2+H_n^{(2)}\right),$$

für den Zielwert \(n=10^8\). Dabei ist \(H_n=\sum_{k=1}^n \frac1k\) die harmonische Zahl und \(H_n^{(2)}=\sum_{k=1}^n \frac1{k^2}\) die harmonische Zahl zweiter Ordnung. Die Implementierungen nutzen diese exakte Identität und wechseln dann zu Euler-Maclaurin-Entwicklungen, weil eine direkte Summe über \(10^8\) Terme unnötig wäre.

Mathematischer Ansatz

Nach der Vereinfachung des zugrunde liegenden Wahrscheinlichkeitsmodells hängt alles nur noch von harmonischen Summen ab. Damit wird ein großes diskretes Problem auf eine kurze analytische Formel zurückgeführt.

Schritt 1: Die exakte geschlossene Form aufschreiben

Wir definieren die gewöhnliche harmonische Zahl

$$H_n=\sum_{k=1}^n\frac1k$$

und die harmonische Zahl zweiter Ordnung

$$H_n^{(2)}=\sum_{k=1}^n\frac1{k^2}.$$

Dann lautet der reduzierte Erwartungswert

$$E(n)=\frac n2\left(H_n^2+H_n^{(2)}\right).$$

Diese Identität ist exakt und nicht nur asymptotisch. Für moderate \(n\) kann man \(E(n)\) daher direkt aus den definierenden Summen berechnen.

Schritt 2: Warum die Doppelsumme zu harmonischen Zahlen wird

Das Ausmultiplizieren des Quadrats ergibt

$$H_n^2=\left(\sum_{i=1}^n\frac1i\right)\left(\sum_{j=1}^n\frac1j\right)=\sum_{i=1}^n\frac1{i^2}+2\sum_{1\le i<j\le n}\frac1{ij}.$$

Durch Addition von \(H_n^{(2)}\) werden die Diagonalterme verdoppelt, also

$$H_n^2+H_n^{(2)}=2\sum_{1\le i\le j\le n}\frac1{ij}.$$

Folglich gilt

$$E(n)=n\sum_{1\le i\le j\le n}\frac1{ij}.$$

Diese Dreieckssumme erklärt, warum in der geschlossenen Form sowohl \(H_n^2\) als auch \(H_n^{(2)}\) auftauchen: Der nichtdiagonale Teil kommt aus dem Quadrat, die Diagonale \(i=j\) liefert die Summe zweiter Ordnung.

Schritt 3: Euler-Maclaurin für große \(n\)

Für den Zielwert ist direktes Summieren unnötig. Stattdessen verwenden die Implementierungen die Standardentwicklungen

$$H_n=\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}+O\left(n^{-10}\right),$$

$$H_n^{(2)}=\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}+O\left(n^{-11}\right).$$

Dabei ist \(\gamma\) die Euler-Mascheroni-Konstante, und \(\pi^2/6\) ist der Grenzwert der harmonischen Summe zweiter Ordnung.

Schritt 4: Die Entwicklungen in die exakte Formel einsetzen

Setzt man diese Reihen in die Identität für \(E(n)\) ein, erhält man für \(n=10^8\) eine praktisch exakte Näherung:

$$E(n)\approx \frac n2\left[\left(\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}\right)^2+\left(\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}\right)\right].$$

Das dominierende Wachstum ist

$$E(n)=\frac n2\left((\log n+\gamma)^2+\frac{\pi^2}{6}\right)+O(\log n).$$

Die Größe wächst also in der Ordnung \(n\log^2 n\), und der abgeschnittene Rest der asymptotischen Reihen ist für \(n=10^8\) viel zu klein, um die nächste ganze Zahl zu verändern.

Schritt 5: Warum das Runden sicher ist

Auch nach der Multiplikation mit \(n/2\) bleiben die ersten ausgelassenen Terme winzig: Der Rest der harmonischen Zahl trägt von der Größenordnung \(n^{-9}\) bei, der Rest von \(H_n^{(2)}\) von der Größenordnung \(n^{-10}\). Für \(n=10^8\) liegen diese analytischen Fehler astronomisch weit unter \(1/2\).

Deshalb können die Implementierungen die abgeschnittenen Reihen mit hoher Präzision auswerten und anschließend gefahrlos auf die nächste ganze Zahl runden.

Durchgerechnetes Beispiel: \(n=5\)

Für \(n=5\) ist direkte Summation leicht:

$$H_5=1+\frac12+\frac13+\frac14+\frac15=\frac{137}{60},$$

$$H_5^{(2)}=1+\frac14+\frac19+\frac1{16}+\frac1{25}=\frac{5269}{3600}.$$

Damit folgt

$$E(5)=\frac52\left[\left(\frac{137}{60}\right)^2+\frac{5269}{3600}\right]=\frac{12019}{720}\approx 16.6930555556.$$

Das stimmt mit dem exakten Kontrollwert der Implementierungen überein und bestätigt die geschlossene Form noch vor jeder asymptotischen Näherung.

Wie der Code funktioniert

Die Implementierungen in C++, Python und Java werten alle dieselbe Formel aus. Zuerst werden die Konstanten \(\gamma\), \(\pi\) und der Zielwert \(n=10^8\) vorbereitet. Danach werden die inversen Potenzen \(n^{-1},n^{-2},\dots,n^{-9}\) gebildet, weil sie in den Euler-Maclaurin-Reihen mehrfach gebraucht werden.

Die C++-Implementierung prüft zusätzlich einige kleine Werte durch direkte Summation und verifiziert damit die exakte harmonische Identität numerisch, bevor der Groß-\(n\)-Pfad verwendet wird. Für die eigentliche Zielberechnung werden die abgeschnittenen Reihen für \(H_n\) und \(H_n^{(2)}\) aufgebaut, in \(\frac n2(H_n^2+H_n^{(2)})\) eingesetzt und anschließend auf die nächste ganze Zahl gerundet.

Die Python- und Java-Implementierungen arbeiten mit 50 Dezimalstellen Präzision. In Java wird der natürliche Logarithmus des Zieleingangs als hochpräzise Konstante hinterlegt, während C++ den Logarithmus direkt in erweiterter Gleitkommaarithmetik berechnet.

Komplexitätsanalyse

Würde man die definierenden Reihen direkt summieren, lägen die Kosten bei \(O(n)\) Zeit und \(O(1)\) zusätzlichem Speicher. Die für den echten Eingabewert verwendete asymptotische Methode braucht dagegen nur eine feste Zahl hochpräziser Rechenoperationen. Bei fester Präzision ist sie also \(O(1)\) in Zeit und \(O(1)\) im Speicher. Genau darum werden die riesigen Summen durch Euler-Maclaurin-Entwicklungen ersetzt.

Fußnoten und Quellen

  1. Problemseite: https://projecteuler.net/problem=724
  2. Harmonische Zahl: Wikipedia - Harmonic number
  3. Verallgemeinerte harmonische Zahl: Wikipedia - Generalized harmonic number
  4. Euler-Maclaurin-Formel: Wikipedia - Euler-Maclaurin formula
  5. Eulersche Konstante: Wikipedia - Euler's constant
  6. Basel-Problem: Wikipedia - Basel problem

Problem Özeti

Bu problemde drone teslimat modelinden çıkan beklenen değer

$$E(n)=n\sum_{1\le i\le j\le n}\frac{1}{ij}=\frac n2\left(H_n^2+H_n^{(2)}\right)$$

ifadesini \(n=10^8\) için hesaplamaya indirgenir. Burada \(H_n=\sum_{k=1}^n \frac1k\) harmonik sayı, \(H_n^{(2)}=\sum_{k=1}^n \frac1{k^2}\) ise ikinci mertebeden harmonik toplamdır. Uygulamalar bu tam özdeşliği kullanır ve ardından Euler-Maclaurin açılımlarına geçer; çünkü \(10^8\) terimi doğrudan toplamak gereksizdir.

Matematiksel Yaklaşım

Alttaki olasılık modeli sadeleştirildiğinde geriye yalnızca harmonik toplamlar kalır. Böylece büyük bir ayrık süreç, kısa bir analitik formüle dönüşür.

Adım 1: Tam kapalı formu yaz

Önce sıradan harmonik sayıyı

$$H_n=\sum_{k=1}^n\frac1k$$

ve ikinci mertebeden harmonik sayıyı

$$H_n^{(2)}=\sum_{k=1}^n\frac1{k^2}$$

tanımlayalım. İndirgenmiş beklenen değer

$$E(n)=\frac n2\left(H_n^2+H_n^{(2)}\right)$$

olur. Bu formül yaklaşık değil, tamdır. Dolayısıyla \(n\) çok büyük değilse \(E(n)\) doğrudan bu iki toplamdan hesaplanabilir.

Adım 2: Çift toplamın neden harmonik sayılara indiğini gör

Karenin açılımı

$$H_n^2=\left(\sum_{i=1}^n\frac1i\right)\left(\sum_{j=1}^n\frac1j\right)=\sum_{i=1}^n\frac1{i^2}+2\sum_{1\le i<j\le n}\frac1{ij}$$

şeklindedir. Buna \(H_n^{(2)}\) eklenince köşegen terimler ikinci kez gelir ve

$$H_n^2+H_n^{(2)}=2\sum_{1\le i\le j\le n}\frac1{ij}$$

elde edilir. Dolayısıyla

$$E(n)=n\sum_{1\le i\le j\le n}\frac1{ij}.$$

Bu üçgensel toplam yorumu, kapalı formda neden hem \(H_n^2\) hem de \(H_n^{(2)}\) bulunduğunu açıklar: köşegen dışı çiftler kareden gelir, \(i=j\) olan köşegen ise ikinci mertebeden harmonik toplamı üretir.

Adım 3: Büyük \(n\) için Euler-Maclaurin kullan

Hedef girdi için doğrudan toplama yapmak gereksizdir. Bunun yerine uygulamalar standart asimptotik açılımları kullanır:

$$H_n=\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}+O\left(n^{-10}\right),$$

$$H_n^{(2)}=\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}+O\left(n^{-11}\right).$$

Burada \(\gamma\) Euler-Mascheroni sabitidir; \(\pi^2/6\) ise ikinci mertebeden harmonik toplamın limitidir.

Adım 4: Açılımları tam formüle yerleştir

Bu serileri \(E(n)\) özdeşliğine koyunca, \(n=10^8\) için fiilen tam olan bir yaklaşım elde edilir:

$$E(n)\approx \frac n2\left[\left(\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}\right)^2+\left(\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}\right)\right].$$

Baskın büyüme terimi

$$E(n)=\frac n2\left((\log n+\gamma)^2+\frac{\pi^2}{6}\right)+O(\log n)$$

şeklindedir. Yani büyüklük mertebesi \(n\log^2 n\) olur ve seri kuyruğunda atılan terimler \(n=10^8\) için en yakın tam sayıyı değiştiremeyecek kadar küçüktür.

Adım 5: En yakın tam sayıya yuvarlama neden güvenli

\(n/2\) ile çarpıldıktan sonra bile ilk atlanan terimler çok küçüktür: harmonik sayı kalanının katkısı \(n^{-9}\) mertebesindedir, ikinci mertebe kalan ise \(n^{-10}\) mertebesindedir. \(n=10^8\) için bu analitik hata \(1/2\)'den astronomik ölçüde küçüktür.

Bu yüzden uygulamalar kesilmiş serileri yüksek hassasiyetle hesaplayıp sonucu güvenle en yakın tam sayıya yuvarlayabilir.

Çözümlü Örnek: \(n=5\)

\(n=5\) için doğrudan toplamak kolaydır:

$$H_5=1+\frac12+\frac13+\frac14+\frac15=\frac{137}{60},$$

$$H_5^{(2)}=1+\frac14+\frac19+\frac1{16}+\frac1{25}=\frac{5269}{3600}.$$

Böylece

$$E(5)=\frac52\left[\left(\frac{137}{60}\right)^2+\frac{5269}{3600}\right]=\frac{12019}{720}\approx 16.6930555556.$$

Bu değer uygulamalardaki tam kontrol sonucuyla aynıdır ve asimptotik yaklaşıma geçmeden önce kapalı formülü doğrular.

Kod Nasıl Çalışır

C++, Python ve Java uygulamalarının hepsi aynı formülü değerlendirir. Önce \(\gamma\), \(\pi\) ve hedef girdi \(n=10^8\) hazırlanır. Sonra \(n^{-1},n^{-2},\dots,n^{-9}\) ters kuvvetleri hesaplanır; çünkü Euler-Maclaurin serilerinde bu terimler tekrar tekrar kullanılır.

C++ uygulaması ayrıca birkaç küçük değeri doğrudan toplayarak tam harmonik özdeşliği sayısal olarak sınar ve ancak ondan sonra büyük-\(n\) yoluna geçer. Hedef hesapta ise uygulamalar \(H_n\) ve \(H_n^{(2)}\) için kesilmiş serileri kurar, bunları \(\frac n2(H_n^2+H_n^{(2)})\) içine yerleştirir ve sonucu en yakın tam sayıya yuvarlar.

Python ve Java uygulamaları 50 basamaklı ondalık aritmetik kullanır. Java sürümü, hedef girdinin doğal logaritmasını yüksek hassasiyetli sabit olarak saklar; C++ sürümü ise logaritmayı genişletilmiş kayan nokta aritmetiğiyle doğrudan hesaplar.

Karmaşıklık Analizi

Tanımdaki seriler doğrudan toplansaydı süre maliyeti \(O(n)\), ek bellek maliyeti \(O(1)\) olurdu. Gerçek girdi için kullanılan asimptotik yöntem ise yalnızca sabit sayıda yüksek hassasiyetli aritmetik işlem yapar. Hassasiyet sabit kabul edilirse zaman karmaşıklığı \(O(1)\), bellek karmaşıklığı \(O(1)\) olur. Büyük toplamları Euler-Maclaurin açılımıyla değiştirme fikrinin amacı tam olarak budur.

Dipnotlar ve Kaynakça

  1. Problem sayfası: https://projecteuler.net/problem=724
  2. Harmonik sayı: Wikipedia - Harmonic number
  3. Genelleştirilmiş harmonik sayı: Wikipedia - Generalized harmonic number
  4. Euler-Maclaurin formülü: Wikipedia - Euler-Maclaurin formula
  5. Euler sabiti: Wikipedia - Euler's constant
  6. Basel problemi: Wikipedia - Basel problem

Resumen del Problema

En este problema, la esperanza asociada al modelo de entrega con drones se reduce a evaluar

$$E(n)=n\sum_{1\le i\le j\le n}\frac{1}{ij}=\frac n2\left(H_n^2+H_n^{(2)}\right),$$

para el valor objetivo \(n=10^8\), donde \(H_n=\sum_{k=1}^n \frac1k\) y \(H_n^{(2)}=\sum_{k=1}^n \frac1{k^2}\). Las implementaciones aprovechan esta identidad exacta y luego pasan a expansiones de Euler-Maclaurin, porque sumar \(10^8\) términos de forma directa no tiene sentido.

Enfoque Matemático

Una vez simplificada la parte probabilística del problema, todo queda expresado en términos de sumas armónicas. Así, un proceso discreto grande se convierte en una fórmula analítica corta.

Paso 1: Escribir la forma cerrada exacta

Definimos el número armónico ordinario

$$H_n=\sum_{k=1}^n\frac1k$$

y el número armónico de segundo orden

$$H_n^{(2)}=\sum_{k=1}^n\frac1{k^2}.$$

Entonces la esperanza reducida es

$$E(n)=\frac n2\left(H_n^2+H_n^{(2)}\right).$$

Esta identidad es exacta, no asintótica. Por eso, para valores moderados de \(n\), \(E(n)\) puede calcularse directamente a partir de las sumas definitorias.

Paso 2: Ver por qué la suma doble se convierte en números armónicos

Al expandir el cuadrado se obtiene

$$H_n^2=\left(\sum_{i=1}^n\frac1i\right)\left(\sum_{j=1}^n\frac1j\right)=\sum_{i=1}^n\frac1{i^2}+2\sum_{1\le i<j\le n}\frac1{ij}.$$

Al sumar \(H_n^{(2)}\), los términos diagonales quedan duplicados, de modo que

$$H_n^2+H_n^{(2)}=2\sum_{1\le i\le j\le n}\frac1{ij}.$$

Por tanto

$$E(n)=n\sum_{1\le i\le j\le n}\frac1{ij}.$$

Esta forma triangular explica por qué aparecen simultáneamente \(H_n^2\) y \(H_n^{(2)}\): la parte fuera de la diagonal procede del cuadrado, mientras que la diagonal \(i=j\) aporta la suma armónica de segundo orden.

Paso 3: Usar Euler-Maclaurin para \(n\) grande

Para la entrada objetivo, la suma directa es innecesaria. Las implementaciones usan las expansiones estándar

$$H_n=\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}+O\left(n^{-10}\right),$$

$$H_n^{(2)}=\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}+O\left(n^{-11}\right).$$

Aquí \(\gamma\) es la constante de Euler-Mascheroni, y \(\pi^2/6\) es el límite de la suma armónica de segundo orden.

Paso 4: Sustituir las expansiones en la fórmula exacta

Al sustituir estas series en la identidad de \(E(n)\), se obtiene una aproximación que para \(n=10^8\) es, en la práctica, exacta:

$$E(n)\approx \frac n2\left[\left(\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}\right)^2+\left(\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}\right)\right].$$

El crecimiento dominante es

$$E(n)=\frac n2\left((\log n+\gamma)^2+\frac{\pi^2}{6}\right)+O(\log n).$$

Así, la magnitud es del orden de \(n\log^2 n\), y la cola omitida de la serie asintótica es demasiado pequeña para modificar el entero más cercano cuando \(n=10^8\).

Paso 5: Por qué el redondeo es seguro

Incluso después de multiplicar por \(n/2\), los primeros términos omitidos siguen siendo minúsculos: el resto de \(H_n\) aporta del orden de \(n^{-9}\), y el resto de \(H_n^{(2)}\) aporta del orden de \(n^{-10}\). Para \(n=10^8\), esos errores analíticos quedan astronómicamente por debajo de \(1/2\).

Por eso las implementaciones pueden evaluar la serie truncada con alta precisión y redondear al entero más cercano sin riesgo.

Ejemplo trabajado: \(n=5\)

Para \(n=5\), la suma directa es sencilla:

$$H_5=1+\frac12+\frac13+\frac14+\frac15=\frac{137}{60},$$

$$H_5^{(2)}=1+\frac14+\frac19+\frac1{16}+\frac1{25}=\frac{5269}{3600}.$$

Entonces

$$E(5)=\frac52\left[\left(\frac{137}{60}\right)^2+\frac{5269}{3600}\right]=\frac{12019}{720}\approx 16.6930555556.$$

Esto coincide con el punto de comprobación exacto usado por las implementaciones y valida la forma cerrada antes de aplicar cualquier aproximación asintótica.

Cómo Funciona el Código

Las implementaciones en C++, Python y Java evalúan exactamente la misma fórmula. Primero preparan las constantes \(\gamma\), \(\pi\) y la entrada objetivo \(n=10^8\). Después construyen las potencias inversas \(n^{-1},n^{-2},\dots,n^{-9}\), porque las expansiones de Euler-Maclaurin reutilizan esos términos muchas veces.

La implementación en C++ además verifica algunos valores pequeños mediante suma directa, confirmando numéricamente la identidad armónica exacta antes de usar la ruta para \(n\) grande. Para el cálculo objetivo, las implementaciones forman las series truncadas de \(H_n\) y \(H_n^{(2)}\), las sustituyen en \(\frac n2(H_n^2+H_n^{(2)})\) y redondean al entero más cercano.

Las implementaciones en Python y Java usan aritmética decimal de 50 cifras. La versión en Java almacena el logaritmo natural de la entrada objetivo como una constante de alta precisión, mientras que la versión en C++ calcula el logaritmo directamente con aritmética de punto flotante extendida.

Análisis de Complejidad

Si se sumaran directamente las series definitorias, el coste sería \(O(n)\) en tiempo y \(O(1)\) en memoria extra. El método asintótico usado para la entrada real solo requiere un número fijo de operaciones aritméticas de alta precisión, así que con precisión fija su coste es \(O(1)\) en tiempo y \(O(1)\) en memoria. Ese es precisamente el beneficio de reemplazar sumas enormes por expansiones de Euler-Maclaurin.

Notas y Referencias

  1. Página del problema: https://projecteuler.net/problem=724
  2. Número armónico: Wikipedia - Harmonic number
  3. Número armónico generalizado: Wikipedia - Generalized harmonic number
  4. Fórmula de Euler-Maclaurin: Wikipedia - Euler-Maclaurin formula
  5. Constante de Euler: Wikipedia - Euler's constant
  6. Problema de Basilea: Wikipedia - Basel problem

问题概述

在这道题里,和无人机配送过程相关的期望值最终可化为

$$E(n)=n\sum_{1\le i\le j\le n}\frac{1}{ij}=\frac n2\left(H_n^2+H_n^{(2)}\right),$$

并要求在目标输入 \(n=10^8\) 时进行计算。这里 \(H_n=\sum_{k=1}^n \frac1k\) 是调和数,\(H_n^{(2)}=\sum_{k=1}^n \frac1{k^2}\) 是二阶调和和。实现采用这条精确恒等式,再配合 Euler-Maclaurin 展开,因为直接累加 \(10^8\) 项没有必要。

数学方法

一旦把题目的概率部分化简完,剩下的内容就只和调和求和有关。这样一个规模很大的离散过程,就被压缩成了少量可直接分析的公式。

步骤 1:写出精确闭式

先定义普通调和数

$$H_n=\sum_{k=1}^n\frac1k$$

以及二阶调和数

$$H_n^{(2)}=\sum_{k=1}^n\frac1{k^2}.$$

化简后的期望值就是

$$E(n)=\frac n2\left(H_n^2+H_n^{(2)}\right).$$

这不是近似式,而是精确公式。因此当 \(n\) 不大时,完全可以直接按定义求出两个和,再代入得到 \(E(n)\)。

步骤 2:说明双重求和为何化为调和数

把平方展开:

$$H_n^2=\left(\sum_{i=1}^n\frac1i\right)\left(\sum_{j=1}^n\frac1j\right)=\sum_{i=1}^n\frac1{i^2}+2\sum_{1\le i<j\le n}\frac1{ij}.$$

再加上 \(H_n^{(2)}\) 后,对角线项会再出现一次,因此

$$H_n^2+H_n^{(2)}=2\sum_{1\le i\le j\le n}\frac1{ij}.$$

于是

$$E(n)=n\sum_{1\le i\le j\le n}\frac1{ij}.$$

这个三角区域求和的形式解释了为什么闭式里同时出现 \(H_n^2\) 和 \(H_n^{(2)}\):非对角项来自平方展开,而对角线 \(i=j\) 恰好对应二阶调和和。

步骤 3:对大 \(n\) 使用 Euler-Maclaurin 展开

对于目标输入,直接求和不划算。实现使用的是标准渐近展开:

$$H_n=\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}+O\left(n^{-10}\right),$$

$$H_n^{(2)}=\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}+O\left(n^{-11}\right).$$

其中 \(\gamma\) 是 Euler-Mascheroni 常数,而 \(\pi^2/6\) 则是二阶调和和的极限值。

步骤 4:把展开式代回精确公式

将上面的级数代入 \(E(n)\) 的恒等式,可得对 \(n=10^8\) 来说几乎就是精确值的表达式:

$$E(n)\approx \frac n2\left[\left(\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}\right)^2+\left(\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}\right)\right].$$

它的主导增长项是

$$E(n)=\frac n2\left((\log n+\gamma)^2+\frac{\pi^2}{6}\right)+O(\log n).$$

也就是说,这个量的规模大约是 \(n\log^2 n\)。而且在 \(n=10^8\) 时,渐近展开尾项被截掉之后的误差小到不可能影响最近整数。

步骤 5:为什么可以放心取最近整数

即使乘上 \(n/2\) 之后,被省略的首项仍然非常小:调和数展开的余项量级是 \(n^{-9}\),二阶调和和展开的余项量级是 \(n^{-10}\)。当 \(n=10^8\) 时,这些解析误差与 \(1/2\) 相比小得不可同日而语。

因此,只要用足够高的数值精度计算截断后的公式,再四舍五入到最近整数即可。

示例:\(n=5\)

当 \(n=5\) 时,直接求和很容易:

$$H_5=1+\frac12+\frac13+\frac14+\frac15=\frac{137}{60},$$

$$H_5^{(2)}=1+\frac14+\frac19+\frac1{16}+\frac1{25}=\frac{5269}{3600}.$$

因此

$$E(5)=\frac52\left[\left(\frac{137}{60}\right)^2+\frac{5269}{3600}\right]=\frac{12019}{720}\approx 16.6930555556.$$

这个结果与实现中使用的精确校验点一致,说明在引入渐近近似之前,闭式本身就是正确的。

代码如何工作

C++、Python 和 Java 三个实现都在计算同一个公式。它们先准备 \(\gamma\)、\(\pi\) 以及目标输入 \(n=10^8\),然后一次性构造 \(n^{-1},n^{-2},\dots,n^{-9}\) 这些倒数幂,因为 Euler-Maclaurin 展开会重复使用它们。

C++ 实现还会用直接求和的方式检查几个较小的数值,从而在进入大 \(n\) 路径之前先对精确调和恒等式做数值验证。对于真正的目标计算,三个实现都会先构造 \(H_n\) 与 \(H_n^{(2)}\) 的截断展开,再代入 \(\frac n2(H_n^2+H_n^{(2)})\),最后取最近整数。

Python 与 Java 使用 50 位十进制高精度运算。Java 版本由于标准高精度十进制算术没有现成的自然对数函数,所以把目标输入的 \(\log(10^8)\) 直接作为高精度常数写入;C++ 版本则用扩展浮点运算直接求对数。

复杂度分析

如果直接按定义累加,时间复杂度是 \(O(n)\),额外空间复杂度是 \(O(1)\)。而实际采用的渐近方法只需要固定次数的高精度加减乘除以及一个对数运算,因此在精度固定的前提下,时间复杂度为 \(O(1)\),空间复杂度也为 \(O(1)\)。这正是用 Euler-Maclaurin 展开替代巨型求和的意义所在。

注释与参考资料

  1. 题目页面:https://projecteuler.net/problem=724
  2. 调和数:Wikipedia - Harmonic number
  3. 广义调和数:Wikipedia - Generalized harmonic number
  4. Euler-Maclaurin 公式:Wikipedia - Euler-Maclaurin formula
  5. Euler 常数:Wikipedia - Euler's constant
  6. 巴塞尔问题:Wikipedia - Basel problem

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

В этой задаче ожидание, возникающее в модели доставки дронами, сводится к вычислению

$$E(n)=n\sum_{1\le i\le j\le n}\frac{1}{ij}=\frac n2\left(H_n^2+H_n^{(2)}\right),$$

для целевого значения \(n=10^8\). Здесь \(H_n=\sum_{k=1}^n \frac1k\) — гармоническое число, а \(H_n^{(2)}=\sum_{k=1}^n \frac1{k^2}\) — гармоническая сумма второго порядка. Реализации используют эту точную формулу и затем переходят к разложениям Эйлера-Маклорена, потому что напрямую суммировать \(10^8\) членов нет смысла.

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

После упрощения вероятностной части задачи все сводится к гармоническим суммам. Тем самым большой дискретный процесс заменяется короткой аналитической формулой.

Шаг 1: Записать точную замкнутую формулу

Введем обычное гармоническое число

$$H_n=\sum_{k=1}^n\frac1k$$

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

$$H_n^{(2)}=\sum_{k=1}^n\frac1{k^2}.$$

Тогда искомое выражение равно

$$E(n)=\frac n2\left(H_n^2+H_n^{(2)}\right).$$

Это точное тождество, а не асимптотика. Поэтому при умеренных \(n\) значение \(E(n)\) можно получать непосредственно из определяющих сумм.

Шаг 2: Почему двойная сумма сворачивается в гармонические числа

Раскрывая квадрат, получаем

$$H_n^2=\left(\sum_{i=1}^n\frac1i\right)\left(\sum_{j=1}^n\frac1j\right)=\sum_{i=1}^n\frac1{i^2}+2\sum_{1\le i<j\le n}\frac1{ij}.$$

Если затем прибавить \(H_n^{(2)}\), диагональные слагаемые удваиваются, и потому

$$H_n^2+H_n^{(2)}=2\sum_{1\le i\le j\le n}\frac1{ij}.$$

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

$$E(n)=n\sum_{1\le i\le j\le n}\frac1{ij}.$$

Эта треугольная сумма объясняет появление сразу двух объектов \(H_n^2\) и \(H_n^{(2)}\): внедиагональная часть приходит из квадрата, а диагональ \(i=j\) дает сумму второго порядка.

Шаг 3: Формула Эйлера-Маклорена для больших \(n\)

Для целевого входа прямое суммирование избыточно. Поэтому реализации используют стандартные асимптотические разложения

$$H_n=\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}+O\left(n^{-10}\right),$$

$$H_n^{(2)}=\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}+O\left(n^{-11}\right).$$

Здесь \(\gamma\) — константа Эйлера-Маскерони, а \(\pi^2/6\) — предел гармонической суммы второго порядка.

Шаг 4: Подстановка разложений в точную формулу

Подставляя эти ряды в тождество для \(E(n)\), получаем приближение, которое для \(n=10^8\) практически совпадает с точным значением:

$$E(n)\approx \frac n2\left[\left(\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}\right)^2+\left(\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}\right)\right].$$

Главный рост задается формулой

$$E(n)=\frac n2\left((\log n+\gamma)^2+\frac{\pi^2}{6}\right)+O(\log n).$$

То есть величина растет как \(n\log^2 n\), а отброшенный хвост асимптотических рядов при \(n=10^8\) слишком мал, чтобы повлиять на ближайшее целое число.

Шаг 5: Почему округление безопасно

Даже после умножения на \(n/2\) первые отброшенные члены остаются ничтожными: остаток для \(H_n\) дает вклад порядка \(n^{-9}\), а остаток для \(H_n^{(2)}\) — порядка \(n^{-10}\). При \(n=10^8\) эти аналитические ошибки астрономически меньше \(1/2\).

Поэтому реализации могут вычислить усеченные ряды с высокой точностью и затем безопасно округлить результат до ближайшего целого.

Разобранный пример: \(n=5\)

Для \(n=5\) прямое суммирование совсем просто:

$$H_5=1+\frac12+\frac13+\frac14+\frac15=\frac{137}{60},$$

$$H_5^{(2)}=1+\frac14+\frac19+\frac1{16}+\frac1{25}=\frac{5269}{3600}.$$

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

$$E(5)=\frac52\left[\left(\frac{137}{60}\right)^2+\frac{5269}{3600}\right]=\frac{12019}{720}\approx 16.6930555556.$$

Это совпадает с точной контрольной точкой, используемой в реализациях, и подтверждает верность замкнутой формулы до перехода к асимптотике.

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

Реализации на C++, Python и Java вычисляют одну и ту же формулу. Сначала подготавливаются константы \(\gamma\), \(\pi\) и целевой вход \(n=10^8\). Затем один раз строятся обратные степени \(n^{-1},n^{-2},\dots,n^{-9}\), потому что они многократно используются в разложениях Эйлера-Маклорена.

Реализация на C++ дополнительно проверяет несколько малых значений прямым суммированием и тем самым численно подтверждает точное гармоническое тождество, прежде чем перейти к ветке для большого \(n\). В целевом вычислении реализации строят усеченные ряды для \(H_n\) и \(H_n^{(2)}\), подставляют их в \(\frac n2(H_n^2+H_n^{(2)})\) и округляют ответ до ближайшего целого.

Реализации на Python и Java используют десятичную арифметику с 50 знаками точности. В Java натуральный логарифм целевого входа хранится как заранее подготовленная высокоточная константа, тогда как версия на C++ вычисляет логарифм непосредственно в расширенной плавающей арифметике.

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

Если бы определяющие ряды суммировались напрямую, это потребовало бы \(O(n)\) времени и \(O(1)\) дополнительной памяти. Асимптотический метод, применяемый для реального входа, использует лишь фиксированное число высокоточных арифметических операций. При фиксированной точности это означает \(O(1)\) по времени и \(O(1)\) по памяти. Именно в этом и состоит выигрыш от замены огромных сумм формулой Эйлера-Маклорена.

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

  1. Страница задачи: https://projecteuler.net/problem=724
  2. Гармоническое число: Wikipedia - Harmonic number
  3. Обобщенное гармоническое число: Wikipedia - Generalized harmonic number
  4. Формула Эйлера-Маклорена: Wikipedia - Euler-Maclaurin formula
  5. Постоянная Эйлера: Wikipedia - Euler's constant
  6. Базельская задача: Wikipedia - Basel problem

ملخص المسألة

في هذه المسألة يختزل مقدار التوقع المرتبط بنموذج التوصيل بالطائرات المسيّرة إلى حساب

$$E(n)=n\sum_{1\le i\le j\le n}\frac{1}{ij}=\frac n2\left(H_n^2+H_n^{(2)}\right),$$

عند القيمة الهدف \(n=10^8\). هنا \(H_n=\sum_{k=1}^n \frac1k\) هو العدد التوافقي، و\(H_n^{(2)}=\sum_{k=1}^n \frac1{k^2}\) هو المجموع التوافقي من الرتبة الثانية. تعتمد التطبيقات على هذه الهوية الدقيقة ثم تنتقل إلى متسلسلات Euler-Maclaurin، لأن جمع \(10^8\) حدًا مباشرةً غير ضروري.

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

بعد تبسيط الجزء الاحتمالي من المسألة، يصبح كل شيء معتمدًا فقط على المجاميع التوافقية. وهكذا يتحول نموذج متقطع كبير إلى صيغة تحليلية قصيرة.

الخطوة 1: كتابة الصيغة المغلقة الدقيقة

نعرّف أولًا العدد التوافقي العادي

$$H_n=\sum_{k=1}^n\frac1k$$

والعدد التوافقي من الرتبة الثانية

$$H_n^{(2)}=\sum_{k=1}^n\frac1{k^2}.$$

عندئذ يصبح التوقع المختزل

$$E(n)=\frac n2\left(H_n^2+H_n^{(2)}\right).$$

هذه هوية تامة وليست تقريبًا. لذلك عندما تكون \(n\) متوسطة الحجم يمكن حساب \(E(n)\) مباشرة من مجموعيه التعريفيين.

الخطوة 2: لماذا ينهار المجموع المزدوج إلى أعداد توافقية

بتوسيع المربع نحصل على

$$H_n^2=\left(\sum_{i=1}^n\frac1i\right)\left(\sum_{j=1}^n\frac1j\right)=\sum_{i=1}^n\frac1{i^2}+2\sum_{1\le i<j\le n}\frac1{ij}.$$

وعند إضافة \(H_n^{(2)}\) تتكرر حدود القطر، ولذلك

$$H_n^2+H_n^{(2)}=2\sum_{1\le i\le j\le n}\frac1{ij}.$$

ومن ثم

$$E(n)=n\sum_{1\le i\le j\le n}\frac1{ij}.$$

وهذا التفسير على شكل مجموع مثلثي يوضح لماذا تظهر الكميتان \(H_n^2\) و\(H_n^{(2)}\) معًا: الأزواج خارج القطر تأتي من تربيع المجموع، أما القطر \(i=j\) فينتج المجموع التوافقي من الرتبة الثانية.

الخطوة 3: استخدام صيغة Euler-Maclaurin عندما يكون \(n\) كبيرًا

بالنسبة إلى القيمة الهدف، لا حاجة إلى الجمع المباشر. لذلك تستخدم التطبيقات التوسعات القياسية

$$H_n=\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}+O\left(n^{-10}\right),$$

$$H_n^{(2)}=\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}+O\left(n^{-11}\right).$$

هنا \(\gamma\) هو ثابت Euler-Mascheroni، و\(\pi^2/6\) هو النهاية التي يقترب منها المجموع التوافقي من الرتبة الثانية.

الخطوة 4: التعويض في الصيغة الدقيقة

عند تعويض هذه السلاسل في هوية \(E(n)\)، نحصل على تقريب يُعد عمليًا مطابقًا للقيمة الدقيقة عندما \(n=10^8\):

$$E(n)\approx \frac n2\left[\left(\log n+\gamma+\frac{1}{2n}-\frac{1}{12n^2}+\frac{1}{120n^4}-\frac{1}{252n^6}+\frac{1}{240n^8}\right)^2+\left(\frac{\pi^2}{6}-\frac1n+\frac{1}{2n^2}-\frac{1}{6n^3}+\frac{1}{30n^5}-\frac{1}{42n^7}+\frac{1}{30n^9}\right)\right].$$

أما السلوك الغالب فهو

$$E(n)=\frac n2\left((\log n+\gamma)^2+\frac{\pi^2}{6}\right)+O(\log n).$$

أي إن المقدار ينمو من رتبة \(n\log^2 n\)، كما أن الذيل المحذوف من التوسع أصغر بكثير من أن يغيّر أقرب عدد صحيح عند \(n=10^8\).

الخطوة 5: لماذا يكون التقريب إلى أقرب عدد صحيح آمنًا

حتى بعد الضرب في \(n/2\)، تبقى أول الحدود المهملة صغيرة للغاية: فباقي توسع \(H_n\) يساهم بمقدار من رتبة \(n^{-9}\)، وباقي توسع \(H_n^{(2)}\) يساهم بمقدار من رتبة \(n^{-10}\). وعند \(n=10^8\) تكون هذه الأخطاء التحليلية أصغر من \(1/2\) بمسافة هائلة.

ولهذا تستطيع التطبيقات تقييم السلسلة المبتورة بدقة عالية ثم تقريب الناتج إلى أقرب عدد صحيح بثقة.

مثال محلول: \(n=5\)

عندما \(n=5\) يكون الجمع المباشر سهلًا:

$$H_5=1+\frac12+\frac13+\frac14+\frac15=\frac{137}{60},$$

$$H_5^{(2)}=1+\frac14+\frac19+\frac1{16}+\frac1{25}=\frac{5269}{3600}.$$

إذًا

$$E(5)=\frac52\left[\left(\frac{137}{60}\right)^2+\frac{5269}{3600}\right]=\frac{12019}{720}\approx 16.6930555556.$$

وهذا يطابق نقطة التحقق الدقيقة المستخدمة في التطبيقات، ويؤكد صحة الصيغة المغلقة قبل الانتقال إلى التقريب الأسيمبتوطي.

كيف يعمل الكود

تقوم تطبيقات C++ وPython وJava جميعها بحساب الصيغة نفسها. فهي تبدأ بتجهيز الثوابت \(\gamma\) و\(\pi\) وقيمة الإدخال الهدف \(n=10^8\)، ثم تبني القوى العكسية \(n^{-1},n^{-2},\dots,n^{-9}\) مرة واحدة لأن متسلسلات Euler-Maclaurin تعيد استخدامها مرارًا.

وتضيف نسخة C++ أيضًا عدة فحوص لقيم صغيرة باستخدام الجمع المباشر، وبذلك تتحقق عدديًا من الهوية التوافقية الدقيقة قبل استخدام مسار \(n\) الكبير. أما في الحساب النهائي فتكوّن التطبيقات التوسعات المبتورة لكل من \(H_n\) و\(H_n^{(2)}\)، ثم تعوض بها في \(\frac n2(H_n^2+H_n^{(2)})\)، ثم تقرّب الناتج إلى أقرب عدد صحيح.

وتستخدم نسختا Python وJava حسابًا عشريًا عالي الدقة من 50 رقمًا. أما نسخة Java فتخزن اللوغاريتم الطبيعي للإدخال الهدف على هيئة ثابت عالي الدقة، في حين تحسبه نسخة C++ مباشرةً باستخدام حساب الفاصلة العائمة الممتد.

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

لو جُمعت السلاسل التعريفية مباشرة لكان الزمن \(O(n)\) والذاكرة الإضافية \(O(1)\). أما الطريقة الأسيمبتوطية المستعملة في الإدخال الحقيقي فتحتاج فقط إلى عدد ثابت من العمليات الحسابية عالية الدقة، ولذلك تكون عند ثبات الدقة \(O(1)\) في الزمن و\(O(1)\) في الذاكرة. وهذا هو الهدف من استبدال المجاميع الضخمة بتوسعات Euler-Maclaurin.

هوامش ومراجع

  1. صفحة المسألة: https://projecteuler.net/problem=724
  2. العدد التوافقي: Wikipedia - Harmonic number
  3. العدد التوافقي المعمم: Wikipedia - Generalized harmonic number
  4. صيغة Euler-Maclaurin: Wikipedia - Euler-Maclaurin formula
  5. ثابت Euler: Wikipedia - Euler's constant
  6. مسألة بازل: Wikipedia - Basel problem