Problem Summary

The Golomb sequence \(G(n)\) is the unique nondecreasing sequence of positive integers in which the integer \(n\) appears exactly \(G(n)\) times. Its beginning is

$$1,2,2,3,3,4,4,4,5,5,5,6,6,6,6,\dots$$

Project Euler 341 asks for

$$\sum_{1\le n \lt 10^6} G(n^3).$$

The official checkpoints \(G(10^3)=86\), \(G(10^6)=6137\), and \(\sum_{1\le n \lt 10^3}G(n^3)=153506976\) show that the queried positions become enormous very quickly. A naive approach that expands the sequence all the way to the largest queried position \((10^6-1)^3\) would be hopelessly too slow and too memory-hungry.

Mathematical Approach

Self-Describing Structure

By definition, value \(k\) occurs exactly \(G(k)\) times. The implementation generates the sequence with the classical recurrence

$$G(1)=1,\qquad G(n)=1+G\bigl(n-G(G(n-1))\bigr)\quad(n\ge 2).$$

This recurrence is equivalent to the usual Project Euler formulation and allows us to precompute \(G(1),G(2),\dots\) once.

First Prefix Sum: Value Boundaries

Define the first prefix sum

$$S(k)=\sum_{i=1}^{k}G(i).$$

Since the value \(k\) appears \(G(k)\) times, the positions holding \(k\) are exactly

$$S(k-1)+1,\ S(k-1)+2,\ \dots,\ S(k).$$

Therefore

$$\boxed{G(n)=k\iff S(k-1)\lt n\le S(k).}$$

So \(S(k)\) marks where the block of value \(k\) ends.

Second Prefix Sum: Blocking the Indices Themselves

Now look at the indices \(m\) for which \(G(m)=k\). From the previous identity, those indices are

$$m=S(k-1)+1,\dots,S(k),$$

and there are exactly \(G(k)\) of them. Each such index \(m\) contributes the value \(k\), so the total length contributed by the entire \(k\)-block is \(k\cdot G(k)\). This motivates the weighted prefix sum

$$P(k)=\sum_{i=1}^{k} i\,G(i).$$

After finishing all blocks up to \(k\), the occupied sequence positions are precisely \(1,\dots,P(k)\).

Inverting a \(P\)-Block

Suppose \(x\) satisfies

$$P(k-1)\lt x\le P(k).$$

Then \(G(x)\) must lie inside the block made of the values

$$S(k-1)+1,\dots,S(k),$$

where each value is repeated exactly \(k\) times. Let

$$d=x-P(k-1).$$

Every consecutive group of \(k\) positions advances the answer by 1, hence the offset inside the block is

$$\left\lceil\frac{d}{k}\right\rceil.$$

So we get the fast query formula

$$\boxed{G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.}$$

This is the core identity used by the function golomb_at.

Worked Example: \(x=10\)

The first boundaries are

$$S(1)=1,\qquad S(2)=3,\qquad S(3)=5,$$

and the weighted boundaries are

$$P(1)=1,\qquad P(2)=5,\qquad P(3)=11.$$

Since \(P(2)=5\lt 10\le 11=P(3)\), the query lies in the \(k=3\) block. Therefore

$$G(10)=S(2)+\left\lceil\frac{10-P(2)}{3}\right\rceil =3+\left\lceil\frac{5}{3}\right\rceil =5.$$

Indeed, the tenth term of the Golomb sequence is 5.

Why the Monotone Sweep Works

The problem only asks for \(x=n^3\), and these queries arrive in strictly increasing order. Therefore the correct block index \(k\) never moves backward. The code keeps the current \(k\), \(S(k)\), and \(P(k)\) in a reusable state object and only advances them when the next cube crosses a new block boundary. Precomputation stops as soon as \(P(k)\) reaches the largest required cube \((10^6-1)^3\).

How the Code Works

precompute_golomb_until(max_query_position) builds the table golomb[n]=G(n) using the recurrence above. While doing so, it also tracks products, which is exactly the weighted prefix sum \(P(k)\). Thus it stops as soon as the largest needed query position is covered.

QueryState stores the current block index \(k\), the current boundaries \(S(k)\) and \(P(k)\), and the previous boundaries \(S(k-1)\) and \(P(k-1)\). In golomb_at(x, golomb, st), a simple while loop advances to the unique block satisfying \(P(k-1)\lt x\le P(k)\), then evaluates

$$G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.$$

The ceiling is implemented with pure integer arithmetic:

$$\left\lceil\frac{a}{b}\right\rceil=\left\lfloor\frac{a+b-1}{b}\right\rfloor.$$

Finally, sum_golomb_cubes(limit, golomb) loops over \(n=1,2,\dots,\text{limit}-1\), computes \(x=n^3\), queries \(G(x)\), and accumulates the total.

Complexity Analysis

Let \(K\) be the largest block index needed so that \(P(K)\ge (10^6-1)^3\). Precomputing \(G(1),\dots,G(K)\) costs \(O(K)\) time and \(O(K)\) memory. The cube sweep costs \(O(10^6+K)\) amortized, because the pointer over \(k\)-blocks only moves forward. There is no per-query binary search and no attempt to materialize the sequence up to index \(10^{18}\).

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=341
  2. Golomb sequence: Wikipedia - Golomb sequence
  3. Prefix sums: Wikipedia - Prefix sum

Problemzusammenfassung

Die Golomb-Folge \(G(n)\) ist die eindeutige nichtfallende Folge positiver ganzer Zahlen, in der die Zahl \(n\) genau \(G(n)\)-mal vorkommt. Die ersten Werte sind

$$1,2,2,3,3,4,4,4,5,5,5,6,6,6,6,\dots$$

Gesucht ist in Problem 341

$$\sum_{1\le n \lt 10^6} G(n^3).$$

Die Kontrollwerte \(G(10^3)=86\), \(G(10^6)=6137\) und \(\sum_{1\le n \lt 10^3}G(n^3)=153506976\) zeigen bereits, wie schnell die abgefragten Positionen wachsen. Ein direktes Ausrollen der Folge bis zur größten Anfrage \((10^6-1)^3\) wäre praktisch unmöglich.

Mathematischer Ansatz

Selbstbeschreibende Struktur

Per Definition erscheint der Wert \(k\) genau \(G(k)\)-mal. Der Code erzeugt die Folge mit der klassischen Rekursion

$$G(1)=1,\qquad G(n)=1+G\bigl(n-G(G(n-1))\bigr)\quad(n\ge 2).$$

Diese Form ist äquivalent zur üblichen Euler-Definition und erlaubt eine einmalige Vorberechnung von \(G(1),G(2),\dots\).

Erste Präfixsumme: Wertgrenzen

Definiere

$$S(k)=\sum_{i=1}^{k}G(i).$$

Da der Wert \(k\) genau \(G(k)\)-mal vorkommt, liegen seine Positionen genau bei

$$S(k-1)+1,\ S(k-1)+2,\ \dots,\ S(k).$$

Daraus folgt

$$\boxed{G(n)=k\iff S(k-1)\lt n\le S(k).}$$

\(S(k)\) markiert also das Ende des Blocks mit Wert \(k\).

Zweite Präfixsumme: Blockstruktur der Indizes

Betrachte nun die Indizes \(m\) mit \(G(m)=k\). Nach der obigen Charakterisierung sind das

$$m=S(k-1)+1,\dots,S(k),$$

also genau \(G(k)\) viele. Jeder dieser Indizes trägt den Wert \(k\), daher hat der gesamte \(k\)-Block die Länge \(k\cdot G(k)\). Deshalb definieren wir die gewichtete Präfixsumme

$$P(k)=\sum_{i=1}^{k} i\,G(i).$$

Nach allen Blöcken bis \(k\) sind genau die Folgepositionen \(1,\dots,P(k)\) belegt.

Inversion eines \(P\)-Blocks

Sei nun

$$P(k-1)\lt x\le P(k).$$

Dann liegt \(G(x)\) im Block aus den Werten

$$S(k-1)+1,\dots,S(k),$$

wobei jeder dieser Werte genau \(k\)-mal wiederholt wird. Setze

$$d=x-P(k-1).$$

Jede Gruppe von \(k\) aufeinanderfolgenden Positionen erhöht die Antwort um 1, also ist der Blockoffset

$$\left\lceil\frac{d}{k}\right\rceil.$$

Somit erhalten wir die schnelle Formel

$$\boxed{G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.}$$

Genau diese Identität verwendet die Funktion golomb_at.

Beispiel: \(x=10\)

Die ersten Grenzen sind

$$S(1)=1,\qquad S(2)=3,\qquad S(3)=5,$$

und gewichtet

$$P(1)=1,\qquad P(2)=5,\qquad P(3)=11.$$

Da \(P(2)=5\lt 10\le 11=P(3)\), liegt die Anfrage im Block \(k=3\). Also

$$G(10)=S(2)+\left\lceil\frac{10-P(2)}{3}\right\rceil =3+\left\lceil\frac{5}{3}\right\rceil =5.$$

Tatsächlich ist der zehnte Term der Golomb-Folge gleich 5.

Warum der monotone Sweep schnell ist

Gefragt werden nur Werte \(x=n^3\), und diese Anfragen sind streng aufsteigend. Deshalb bewegt sich der richtige Blockindex \(k\) nie rückwärts. Der Code speichert \(k\), \(S(k)\) und \(P(k)\) in einem Zustand und schiebt ihn nur dann weiter, wenn der nächste Würfel eine neue Blockgrenze überschreitet. Die Vorberechnung endet, sobald \(P(k)\) den größten benötigten Würfel \((10^6-1)^3\) erreicht.

Wie der Code funktioniert

precompute_golomb_until(max_query_position) baut die Tabelle golomb[n]=G(n) mit der Rekursion oben auf. Gleichzeitig wird products mitgeführt; das ist genau die gewichtete Präfixsumme \(P(k)\). Dadurch wird nur so weit vorab gerechnet, wie es für die größte Anfrage nötig ist.

QueryState speichert den aktuellen Blockindex \(k\), die aktuellen Grenzen \(S(k)\) und \(P(k)\) sowie die vorherigen Grenzen \(S(k-1)\) und \(P(k-1)\). In golomb_at(x, golomb, st) wird per while-Schleife bis zum eindeutigen Block mit \(P(k-1)\lt x\le P(k)\) vorgerückt; danach wird

$$G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil$$

mit Ganzzahlarithmetik ausgewertet. Die Aufrundung wird als

$$\left\lceil\frac{a}{b}\right\rceil=\left\lfloor\frac{a+b-1}{b}\right\rfloor$$

implementiert. Schließlich summiert sum_golomb_cubes(limit, golomb) alle Werte \(G(n^3)\) für \(1\le n\lt \text{limit}\).

Komplexitätsanalyse

Sei \(K\) der größte benötigte Blockindex mit \(P(K)\ge (10^6-1)^3\). Die Vorberechnung von \(G(1),\dots,G(K)\) kostet \(O(K)\) Zeit und \(O(K)\) Speicher. Die Sweep-Phase über alle Würfel kostet amortisiert \(O(10^6+K)\), weil der Zeiger über die Blöcke nur vorwärts läuft. Weder ist eine binäre Suche pro Anfrage nötig noch eine Materialisierung der Folge bis Index \(10^{18}\).

Fußnoten und Quellen

  1. Problemseite: https://projecteuler.net/problem=341
  2. Golomb-Folge: Wikipedia - Golomb-Folge
  3. Präfixsummen: Wikipedia - Partialsumme

Problem Özeti

Golomb dizisi \(G(n)\), \(n\) sayısının dizide tam olarak \(G(n)\) kez göründüğü tek azalmayan pozitif tamsayı dizisidir. İlk terimleri şöyledir:

$$1,2,2,3,3,4,4,4,5,5,5,6,6,6,6,\dots$$

Project Euler 341'de istenen toplam

$$\sum_{1\le n \lt 10^6} G(n^3)$$

şeklindedir. Resmî kontrol değerleri \(G(10^3)=86\), \(G(10^6)=6137\) ve \(\sum_{1\le n \lt 10^3}G(n^3)=153506976\) olduğundan, sorgu konumlarının çok hızlı büyüdüğü görülür. En büyük sorgu \((10^6-1)^3\) mertebesinde olduğu için diziyi doğrudan o indekse kadar üretmek zaman ve bellek açısından pratik değildir.

Matematiksel Yaklaşım

Kendini Tanımlayan Yapı

Tanım gereği \(k\) değeri dizide tam \(G(k)\) kez görünür. Kod, diziyi şu klasik bağıntı ile üretir:

$$G(1)=1,\qquad G(n)=1+G\bigl(n-G(G(n-1))\bigr)\quad(n\ge 2).$$

Bu bağıntı, Project Euler'daki standart tanımla eşdeğerdir ve \(G(1),G(2),\dots\) değerlerini tek seferde önceden hesaplamaya izin verir.

Birinci Prefix Toplamı: Değer Bloklarının Sınırları

Önce

$$S(k)=\sum_{i=1}^{k}G(i)$$

tanımını yapalım. \(k\) değeri tam \(G(k)\) kez göründüğü için, dizide \(k\)'nın bulunduğu konumlar tam olarak

$$S(k-1)+1,\ S(k-1)+2,\ \dots,\ S(k)$$

olur. Dolayısıyla

$$\boxed{G(n)=k\iff S(k-1)\lt n\le S(k).}$$

Yani \(S(k)\), \(k\) değerli bloğun bittiği konumu gösterir.

İkinci Prefix Toplamı: İndis Bloklarını Sıkıştırmak

Şimdi \(G(m)=k\) olan indisleri düşünelim. Yukarıdaki eşdeğerliğe göre bunlar

$$m=S(k-1)+1,\dots,S(k)$$

şeklindedir ve sayıları tam \(G(k)\)'dır. Bu indislerin her biri diziye \(k\) değerini katkılar; dolayısıyla tüm \(k\)-bloğunun uzunluğu \(k\cdot G(k)\) olur. Bu yüzden ağırlıklı prefix toplamı

$$P(k)=\sum_{i=1}^{k} i\,G(i)$$

tanımlanır. \(k\)'ya kadar olan tüm bloklar işlendiğinde kaplanan dizi konumları tam olarak \(1,\dots,P(k)\) olur.

\(P\)-Bloğunu Terslemek

Şimdi

$$P(k-1)\lt x\le P(k)$$

olsun. Bu durumda \(G(x)\),

$$S(k-1)+1,\dots,S(k)$$

değerlerinden oluşan blok içindedir ve bu değerlerin her biri tam \(k\) kez tekrar edilir. Eğer

$$d=x-P(k-1)$$

dersek, her \(k\) ardışık konum cevapta 1 artış yapar. Bu yüzden blok içi kayma

$$\left\lceil\frac{d}{k}\right\rceil$$

olur ve hızlı sorgu formülü elde edilir:

$$\boxed{G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.}$$

golomb_at fonksiyonunun kullandığı temel kimlik tam olarak budur.

Örnek: \(x=10\)

İlk sınırlar

$$S(1)=1,\qquad S(2)=3,\qquad S(3)=5,$$

ağırlıklı sınırlar ise

$$P(1)=1,\qquad P(2)=5,\qquad P(3)=11$$

şeklindedir. \(P(2)=5\lt 10\le 11=P(3)\) olduğundan sorgu \(k=3\) bloğundadır. O halde

$$G(10)=S(2)+\left\lceil\frac{10-P(2)}{3}\right\rceil =3+\left\lceil\frac{5}{3}\right\rceil =5.$$

Gerçekten de Golomb dizisinin onuncu terimi 5'tir.

Artan Tarama Neden Hızlı?

Problem yalnızca \(x=n^3\) türünden sorgular ister ve bu sorgular sıkı biçimde artan sırada gelir. Bu nedenle doğru blok indeksi \(k\) hiçbir zaman geriye gitmez. Kod, mevcut \(k\), \(S(k)\) ve \(P(k)\) değerlerini bir durum nesnesinde saklar; bir sonraki küp yeni bir blok sınırını geçene kadar bunları yeniden kullanır. Ön hesaplama da \(P(k)\), en büyük gerekli küp olan \((10^6-1)^3\)'e ulaştığında durur.

Kodun Çalışma Mantığı

precompute_golomb_until(max_query_position), yukarıdaki bağıntı ile golomb[n]=G(n) tablosunu üretir. Aynı anda products değişkenini de günceller; bu değişken doğrudan ağırlıklı prefix toplamı \(P(k)\)'dır. Böylece yalnızca en büyük sorguyu kapsayacak kadar Golomb değeri üretilir.

QueryState, mevcut blok indeksi \(k\), güncel sınırlar \(S(k)\) ve \(P(k)\), ayrıca önceki sınırlar \(S(k-1)\) ve \(P(k-1)\) bilgisini taşır. golomb_at(x, golomb, st) içinde önce \(P(k-1)\lt x\le P(k)\) koşulunu sağlayan tek blok bulunur; ardından

$$G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil$$

formülü uygulanır. Tavana yuvarlama saf tamsayı aritmetiğiyle

$$\left\lceil\frac{a}{b}\right\rceil=\left\lfloor\frac{a+b-1}{b}\right\rfloor$$

eşitliği üzerinden yapılır. Son olarak sum_golomb_cubes(limit, golomb), \(1\le n\lt \text{limit}\) için tüm \(G(n^3)\) değerlerini toplayarak sonucu üretir.

Karmaşıklık Analizi

\(P(K)\ge (10^6-1)^3\) olacak en büyük gerekli blok indisi \(K\) olsun. \(G(1),\dots,G(K)\) ön hesabı \(O(K)\) zaman ve \(O(K)\) bellek kullanır. Küpler üzerinde yapılan tarama amortize olarak \(O(10^6+K)\) maliyetlidir; çünkü blok göstergesi yalnızca ileri gider. Her sorgu için ayrı bir ikili arama yapılmaz ve dizi \(10^{18}\) civarı bir indekse kadar fiziksel olarak açılmaz.

Dipnotlar ve Kaynakça

  1. Problem sayfası: https://projecteuler.net/problem=341
  2. Golomb dizisi: Wikipedia - Golomb sequence
  3. Prefix toplamlar: Wikipedia - Prefix sum

Resumen del Problema

La sucesión de Golomb \(G(n)\) es la única sucesión no decreciente de enteros positivos en la que el número \(n\) aparece exactamente \(G(n)\) veces. Sus primeros términos son

$$1,2,2,3,3,4,4,4,5,5,5,6,6,6,6,\dots$$

El problema pide calcular

$$\sum_{1\le n \lt 10^6} G(n^3).$$

Los puntos de control \(G(10^3)=86\), \(G(10^6)=6137\) y \(\sum_{1\le n \lt 10^3}G(n^3)=153506976\) muestran que las posiciones consultadas crecen muy deprisa. Desplegar la sucesión directamente hasta \((10^6-1)^3\) sería inviable.

Enfoque Matemático

Estructura Auto-descriptiva

Por definición, el valor \(k\) aparece exactamente \(G(k)\) veces. La implementación genera la sucesión con la recurrencia clásica

$$G(1)=1,\qquad G(n)=1+G\bigl(n-G(G(n-1))\bigr)\quad(n\ge 2).$$

Esta forma es equivalente a la definición habitual y permite precomputar \(G(1),G(2),\dots\) una sola vez.

Primer Prefijo: Fronteras de los Valores

Definimos

$$S(k)=\sum_{i=1}^{k}G(i).$$

Como el valor \(k\) aparece \(G(k)\) veces, sus posiciones son exactamente

$$S(k-1)+1,\ S(k-1)+2,\ \dots,\ S(k).$$

Por tanto,

$$\boxed{G(n)=k\iff S(k-1)\lt n\le S(k).}$$

Segundo Prefijo: Agrupar los Índices

Ahora consideremos los índices \(m\) tales que \(G(m)=k\). Son precisamente

$$m=S(k-1)+1,\dots,S(k),$$

y hay exactamente \(G(k)\) de ellos. Cada uno aporta el valor \(k\), así que la longitud total del bloque \(k\) es \(k\cdot G(k)\). Esto motiva el prefijo ponderado

$$P(k)=\sum_{i=1}^{k} i\,G(i).$$

Tras procesar todos los bloques hasta \(k\), las posiciones ocupadas son exactamente \(1,\dots,P(k)\).

Invertir un Bloque \(P\)

Si

$$P(k-1)\lt x\le P(k),$$

entonces \(G(x)\) está dentro del bloque formado por

$$S(k-1)+1,\dots,S(k),$$

donde cada valor se repite exactamente \(k\) veces. Sea

$$d=x-P(k-1).$$

Cada grupo consecutivo de \(k\) posiciones incrementa la respuesta en 1, de modo que el desplazamiento interno es

$$\left\lceil\frac{d}{k}\right\rceil.$$

Así obtenemos

$$\boxed{G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.}$$

Ejemplo: \(x=10\)

Tenemos

$$S(1)=1,\qquad S(2)=3,\qquad S(3)=5,$$

y también

$$P(1)=1,\qquad P(2)=5,\qquad P(3)=11.$$

Como \(5=P(2)\lt 10\le P(3)=11\), estamos en el bloque \(k=3\). Entonces

$$G(10)=S(2)+\left\lceil\frac{10-P(2)}{3}\right\rceil =3+\left\lceil\frac{5}{3}\right\rceil =5.$$

Por Qué el Barrido Monótono es Eficiente

Las consultas son solo de la forma \(x=n^3\) y llegan en orden creciente. Por eso el índice correcto de bloque \(k\) nunca retrocede. El código mantiene \(k\), \(S(k)\) y \(P(k)\) en un estado reutilizable y solo los avanza cuando hace falta. La precomputación termina en cuanto \(P(k)\) alcanza el mayor cubo requerido \((10^6-1)^3\).

Cómo Funciona el Código

precompute_golomb_until(max_query_position) construye la tabla golomb[n]=G(n) usando la recurrencia anterior. A la vez, mantiene products, que coincide con el prefijo ponderado \(P(k)\). Así solo se generan los términos necesarios.

QueryState guarda el bloque actual \(k\), las fronteras \(S(k)\) y \(P(k)\), y las fronteras previas \(S(k-1)\) y \(P(k-1)\). En golomb_at(x, golomb, st) se avanza hasta el único bloque con \(P(k-1)\lt x\le P(k)\), y luego se evalúa

$$G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.$$

La función techo se implementa con aritmética entera:

$$\left\lceil\frac{a}{b}\right\rceil=\left\lfloor\frac{a+b-1}{b}\right\rfloor.$$

Por último, sum_golomb_cubes(limit, golomb) recorre \(n=1,2,\dots,\text{limit}-1\), consulta \(G(n^3)\) y acumula la suma.

Complejidad

Sea \(K\) el mayor índice necesario tal que \(P(K)\ge (10^6-1)^3\). La precomputación cuesta \(O(K)\) tiempo y \(O(K)\) memoria. El barrido de todas las consultas cuesta \(O(10^6+K)\) amortizado, porque el puntero de bloques solo avanza. No hay búsqueda binaria por consulta ni expansión directa hasta índices del orden de \(10^{18}\).

Lecturas y Referencias

  1. Problema: https://projecteuler.net/problem=341
  2. Sucesión de Golomb: Wikipedia - Sucesión de Golomb
  3. Sumas prefijo: Wikipedia - Prefix sum

问题概述

Golomb 序列 \(G(n)\) 是唯一的非递减正整数序列,并满足整数 \(n\) 在序列中恰好出现 \(G(n)\) 次。前几项为

$$1,2,2,3,3,4,4,4,5,5,5,6,6,6,6,\dots$$

题目要求计算

$$\sum_{1\le n \lt 10^6} G(n^3).$$

官方给出的检查点 \(G(10^3)=86\)、\(G(10^6)=6137\) 和 \(\sum_{1\le n \lt 10^3}G(n^3)=153506976\) 说明查询位置增长极快。若直接把序列展开到最大查询位置 \((10^6-1)^3\),无论时间还是内存都不可接受。

数学方法

自描述结构

按照定义,值 \(k\) 恰好出现 \(G(k)\) 次。代码用经典递推式生成序列:

$$G(1)=1,\qquad G(n)=1+G\bigl(n-G(G(n-1))\bigr)\quad(n\ge 2).$$

这与题目的标准定义等价,可以一次性预计算 \(G(1),G(2),\dots\)。

第一层前缀和:值块边界

定义

$$S(k)=\sum_{i=1}^{k}G(i).$$

由于值 \(k\) 一共出现 \(G(k)\) 次,所以所有等于 \(k\) 的位置恰好是

$$S(k-1)+1,\ S(k-1)+2,\ \dots,\ S(k).$$

因此

$$\boxed{G(n)=k\iff S(k-1)\lt n\le S(k).}$$

第二层加权前缀和:压缩索引块

再看满足 \(G(m)=k\) 的那些索引 \(m\)。根据上式,它们正是

$$m=S(k-1)+1,\dots,S(k),$$

共 \(G(k)\) 个。每个这样的索引都对应值 \(k\),所以整个 \(k\)-块贡献的总长度为 \(k\cdot G(k)\)。于是定义加权前缀和

$$P(k)=\sum_{i=1}^{k} i\,G(i).$$

处理完前 \(k\) 个块后,被覆盖的序列位置正好是 \(1,\dots,P(k)\)。

\(P\)-块反演公式

$$P(k-1)\lt x\le P(k),$$

则 \(G(x)\) 落在由

$$S(k-1)+1,\dots,S(k)$$

组成的那一块中,而且这些值都恰好重复 \(k\) 次。记

$$d=x-P(k-1).$$

每经过连续 \(k\) 个位置,答案就增加 1,因此块内偏移为

$$\left\lceil\frac{d}{k}\right\rceil.$$

于是得到快速查询公式

$$\boxed{G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.}$$

例子:\(x=10\)

前几个边界为

$$S(1)=1,\qquad S(2)=3,\qquad S(3)=5,$$

对应的加权边界为

$$P(1)=1,\qquad P(2)=5,\qquad P(3)=11.$$

因为 \(P(2)=5\lt 10\le 11=P(3)\),所以它位于 \(k=3\) 这一块。于是

$$G(10)=S(2)+\left\lceil\frac{10-P(2)}{3}\right\rceil =3+\left\lceil\frac{5}{3}\right\rceil =5.$$

为什么递增扫描很快

题目只查询 \(x=n^3\),而这些查询严格递增。因此正确的块编号 \(k\) 永远不会回退。代码把当前 \(k\)、\(S(k)\)、\(P(k)\) 保存在可复用状态中,只在下一个立方查询越过块边界时才继续推进。预计算也会在 \(P(k)\) 达到最大所需立方 \((10^6-1)^3\) 后停止。

代码原理

precompute_golomb_until(max_query_position) 根据上面的递推式构造数组 golomb[n]=G(n)。同时维护变量 products,它正是加权前缀和 \(P(k)\)。因此只会生成覆盖最大查询所需的那些项。

QueryState 保存当前块编号 \(k\)、当前边界 \(S(k)\) 与 \(P(k)\),以及前一块边界 \(S(k-1)\) 与 \(P(k-1)\)。在 golomb_at(x, golomb, st) 中,先通过 while 循环找到唯一满足 \(P(k-1)\lt x\le P(k)\) 的块,再计算

$$G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.$$

其中上取整通过整数恒等式实现:

$$\left\lceil\frac{a}{b}\right\rceil=\left\lfloor\frac{a+b-1}{b}\right\rfloor.$$

最后,sum_golomb_cubes(limit, golomb) 依次枚举 \(n=1,2,\dots,\text{limit}-1\),查询 \(G(n^3)\) 并累加答案。

复杂度分析

设 \(K\) 是满足 \(P(K)\ge (10^6-1)^3\) 的最大所需块编号。预计算 \(G(1),\dots,G(K)\) 需要 \(O(K)\) 时间和 \(O(K)\) 空间。所有立方查询的扫描阶段摊还复杂度为 \(O(10^6+K)\),因为块指针只会向前移动。算法既不需要对每次查询做二分搜索,也不需要真的把序列展开到 \(10^{18}\) 级别的位置。

参考资料

  1. 题目页面: https://projecteuler.net/problem=341
  2. Golomb 序列: Wikipedia - Golomb sequence
  3. 前缀和: Wikipedia - Prefix sum

Краткое описание

Последовательность Голомба \(G(n)\) - это единственная неубывающая последовательность натуральных чисел, в которой число \(n\) встречается ровно \(G(n)\) раз. Начало последовательности:

$$1,2,2,3,3,4,4,4,5,5,5,6,6,6,6,\dots$$

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

$$\sum_{1\le n \lt 10^6} G(n^3).$$

Контрольные значения \(G(10^3)=86\), \(G(10^6)=6137\) и \(\sum_{1\le n \lt 10^3}G(n^3)=153506976\) показывают, что нужные позиции растут очень быстро. Прямое построение последовательности до \((10^6-1)^3\) было бы нереалистичным.

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

Самоописывающая структура

По определению значение \(k\) встречается ровно \(G(k)\) раз. Код генерирует последовательность по классической рекуррентной формуле

$$G(1)=1,\qquad G(n)=1+G\bigl(n-G(G(n-1))\bigr)\quad(n\ge 2).$$

Она эквивалентна стандартной формулировке задачи и позволяет один раз предварительно вычислить \(G(1),G(2),\dots\).

Первая префиксная сумма: границы блоков значений

Введем

$$S(k)=\sum_{i=1}^{k}G(i).$$

Так как значение \(k\) встречается \(G(k)\) раз, его позиции - это ровно

$$S(k-1)+1,\ S(k-1)+2,\ \dots,\ S(k).$$

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

$$\boxed{G(n)=k\iff S(k-1)\lt n\le S(k).}$$

Вторая префиксная сумма: сжатие блоков индексов

Теперь посмотрим на индексы \(m\), для которых \(G(m)=k\). По предыдущему критерию это

$$m=S(k-1)+1,\dots,S(k),$$

и таких индексов ровно \(G(k)\). Каждый из них дает значение \(k\), поэтому полный \(k\)-блок имеет длину \(k\cdot G(k)\). Отсюда естественно ввести взвешенную префиксную сумму

$$P(k)=\sum_{i=1}^{k} i\,G(i).$$

После обработки всех блоков до \(k\) занятыми оказываются ровно позиции \(1,\dots,P(k)\).

Инверсия по \(P\)-блоку

Пусть

$$P(k-1)\lt x\le P(k).$$

Тогда \(G(x)\) лежит внутри блока, составленного из значений

$$S(k-1)+1,\dots,S(k),$$

и каждое из них повторяется ровно \(k\) раз. Обозначим

$$d=x-P(k-1).$$

Каждая подряд идущая группа из \(k\) позиций увеличивает ответ на 1, поэтому смещение внутри блока равно

$$\left\lceil\frac{d}{k}\right\rceil.$$

Отсюда получаем формулу быстрого запроса

$$\boxed{G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.}$$

Пример: \(x=10\)

Первые границы:

$$S(1)=1,\qquad S(2)=3,\qquad S(3)=5,$$

а взвешенные границы:

$$P(1)=1,\qquad P(2)=5,\qquad P(3)=11.$$

Так как \(P(2)=5\lt 10\le 11=P(3)\), запрос попадает в блок \(k=3\). Значит,

$$G(10)=S(2)+\left\lceil\frac{10-P(2)}{3}\right\rceil =3+\left\lceil\frac{5}{3}\right\rceil =5.$$

Почему возрастающий проход эффективен

В задаче нужны только запросы вида \(x=n^3\), и они идут строго по возрастанию. Поэтому правильный индекс блока \(k\) никогда не движется назад. Код хранит текущие \(k\), \(S(k)\) и \(P(k)\) в переиспользуемом состоянии и продвигает их только тогда, когда следующий куб пересекает новую границу блока. Предвычисление останавливается, когда \(P(k)\) достигает максимального нужного куба \((10^6-1)^3\).

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

precompute_golomb_until(max_query_position) строит таблицу golomb[n]=G(n) по рекуррентной формуле выше. Одновременно поддерживается переменная products, которая в точности равна взвешенной префиксной сумме \(P(k)\). Поэтому генерируются только действительно нужные элементы.

QueryState хранит текущий индекс блока \(k\), текущие границы \(S(k)\) и \(P(k)\), а также предыдущие границы \(S(k-1)\) и \(P(k-1)\). В golomb_at(x, golomb, st) цикл while доходит до единственного блока с \(P(k-1)\lt x\le P(k)\), после чего вычисляется

$$G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.$$

Потолок реализован чисто целочисленно:

$$\left\lceil\frac{a}{b}\right\rceil=\left\lfloor\frac{a+b-1}{b}\right\rfloor.$$

Далее sum_golomb_cubes(limit, golomb) перебирает \(n=1,2,\dots,\text{limit}-1\), запрашивает \(G(n^3)\) и суммирует ответ.

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

Пусть \(K\) - максимальный нужный индекс блока, для которого \(P(K)\ge (10^6-1)^3\). Предвычисление \(G(1),\dots,G(K)\) занимает \(O(K)\) времени и \(O(K)\) памяти. Проход по всем кубам стоит \(O(10^6+K)\) амортизированно, потому что указатель по блокам движется только вперед. Ни двоичный поиск для каждого запроса, ни материализация последовательности до индекса порядка \(10^{18}\) не требуются.

Дополнительные материалы

  1. Условие задачи: https://projecteuler.net/problem=341
  2. Последовательность Голомба: Wikipedia - Последовательность Голомба
  3. Префиксные суммы: Wikipedia - Префиксная сумма

ملخص المسألة

متتالية Golomb \(G(n)\) هي المتتالية غير المتناقصة الوحيدة من الأعداد الطبيعية التي يظهر فيها العدد \(n\) بالضبط \(G(n)\) مرة. بدايتها هي

$$1,2,2,3,3,4,4,4,5,5,5,6,6,6,6,\dots$$

والمطلوب في المسألة هو حساب

$$\sum_{1\le n \lt 10^6} G(n^3).$$

قيم التحقق \(G(10^3)=86\) و\(G(10^6)=6137\) و\(\sum_{1\le n \lt 10^3}G(n^3)=153506976\) توضّح أن مواضع الاستعلام تكبر بسرعة كبيرة. لذلك فإن توليد المتتالية مباشرة حتى الموضع الأكبر \((10^6-1)^3\) غير عملي إطلاقًا.

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

البنية ذاتية الوصف

بحسب التعريف، القيمة \(k\) تظهر بالضبط \(G(k)\) مرة. ويولّد الكود المتتالية بواسطة العلاقة التكرارية الكلاسيكية

$$G(1)=1,\qquad G(n)=1+G\bigl(n-G(G(n-1))\bigr)\quad(n\ge 2).$$

وهذه الصيغة مكافئة لصياغة Project Euler المعتادة، وتسمح بحساب \(G(1),G(2),\dots\) مسبقًا مرة واحدة.

المجموع التراكمي الأول: حدود كتل القيم

لنعرّف

$$S(k)=\sum_{i=1}^{k}G(i).$$

وبما أن القيمة \(k\) تتكرر \(G(k)\) مرة، فإن جميع المواضع التي تحمل \(k\) هي بالضبط

$$S(k-1)+1,\ S(k-1)+2,\ \dots,\ S(k).$$

ومن ثم

$$\boxed{G(n)=k\iff S(k-1)\lt n\le S(k).}$$

المجموع التراكمي الثاني: ضغط كتل الفهارس

الآن انظر إلى الفهارس \(m\) التي تحقق \(G(m)=k\). من العلاقة السابقة تكون هذه الفهارس

$$m=S(k-1)+1,\dots,S(k),$$

وعددها بالضبط \(G(k)\). وكل فهرس من هذه الفهارس يضيف القيمة \(k\)، لذلك يكون طول كتلة \(k\) كاملًا مساويًا لـ \(k\cdot G(k)\). ولهذا نعرّف المجموع التراكمي الموزون

$$P(k)=\sum_{i=1}^{k} i\,G(i).$$

وبعد إنهاء جميع الكتل حتى \(k\)، تكون المواضع المغطاة في المتتالية هي تمامًا \(1,\dots,P(k)\).

عكس كتلة \(P\)

إذا كان

$$P(k-1)\lt x\le P(k),$$

فإن \(G(x)\) يقع داخل الكتلة المكوّنة من القيم

$$S(k-1)+1,\dots,S(k),$$

مع تكرار كل قيمة منها \(k\) مرات تمامًا. ولنضع

$$d=x-P(k-1).$$

كل مجموعة متتالية من \(k\) مواضع ترفع الجواب بمقدار 1، ولذلك يكون الإزاحة داخل الكتلة

$$\left\lceil\frac{d}{k}\right\rceil.$$

فنحصل على صيغة الاستعلام السريع

$$\boxed{G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.}$$

مثال: \(x=10\)

لدينا

$$S(1)=1,\qquad S(2)=3,\qquad S(3)=5,$$

وكذلك

$$P(1)=1,\qquad P(2)=5,\qquad P(3)=11.$$

بما أن \(P(2)=5\lt 10\le 11=P(3)\)، فإن الاستعلام يقع في الكتلة \(k=3\). إذًا

$$G(10)=S(2)+\left\lceil\frac{10-P(2)}{3}\right\rceil =3+\left\lceil\frac{5}{3}\right\rceil =5.$$

لماذا المسح التصاعدي سريع؟

المسألة لا تطلب إلا استعلامات من الشكل \(x=n^3\)، وهذه تأتي بترتيب تصاعدي صارم. لذلك فإن فهرس الكتلة الصحيح \(k\) لا يعود إلى الخلف أبدًا. يحفظ الكود القيم الحالية \(k\) و\(S(k)\) و\(P(k)\) داخل حالة قابلة لإعادة الاستخدام، ولا يحرّكها إلا عندما يتجاوز المكعب التالي حد كتلة جديدة. كما تتوقف المعالجة المسبقة عندما يصل \(P(k)\) إلى أكبر مكعب مطلوب وهو \((10^6-1)^3\).

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

الدالة precompute_golomb_until(max_query_position) تبني الجدول golomb[n]=G(n) باستخدام العلاقة التكرارية السابقة. وخلال ذلك تحفظ المتغير products، وهو يساوي مباشرة المجموع التراكمي الموزون \(P(k)\). لذلك لا يتم توليد إلا العدد اللازم من حدود Golomb.

أما QueryState فيحفظ فهرس الكتلة الحالي \(k\)، والحدود الحالية \(S(k)\) و\(P(k)\)، والحدود السابقة \(S(k-1)\) و\(P(k-1)\). داخل golomb_at(x, golomb, st) يتم التقدم بحلقة while حتى الكتلة الوحيدة التي تحقق \(P(k-1)\lt x\le P(k)\)، ثم تُطبّق الصيغة

$$G(x)=S(k-1)+\left\lceil\frac{x-P(k-1)}{k}\right\rceil.$$

ويُنَفَّذ التقريب للأعلى بواسطة حساب صحيح صرف:

$$\left\lceil\frac{a}{b}\right\rceil=\left\lfloor\frac{a+b-1}{b}\right\rfloor.$$

وأخيرًا تجمع الدالة sum_golomb_cubes(limit, golomb) جميع القيم \(G(n^3)\) عندما \(1\le n\lt \text{limit}\).

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

ليكن \(K\) أكبر فهرس كتلة نحتاجه بحيث \(P(K)\ge (10^6-1)^3\). عندئذ تكلف المعالجة المسبقة \(O(K)\) زمنًا و\(O(K)\) ذاكرة. أما المرور على جميع المكعبات فيكلف \(O(10^6+K)\) على نحو مُموَّه، لأن مؤشر الكتلة يتحرك إلى الأمام فقط. لا توجد عملية بحث ثنائي لكل استعلام، ولا حاجة لتوسيع المتتالية حتى موضع يقارب \(10^{18}\).

مراجع إضافية

  1. صفحة المسألة: https://projecteuler.net/problem=341
  2. متتالية Golomb: Wikipedia - Golomb sequence
  3. المجاميع التراكمية: Wikipedia - Prefix sum