The problem defines a sequence \(D(n)\) associated with counting amoeba configurations in a three-dimensional grid. The computational goal used by the implementations is to evaluate \(D(10000)\) modulo \(10^9\), while several smaller values are checked exactly before the large modular run.
A direct enumeration of all relevant 3D shapes would grow far too quickly. The implemented method therefore replaces explicit shape generation with a layered dynamic program whose states encode admissible boundary profiles and whose scalar recurrences extract the final sequence value.
For exposition, denote the two triangular state families by \(A_{n,k}(m)\) and \(B_{n,k}(m)\), and denote the scalar sequences by \(P(m)\) and \(Q(m)\). The desired quantity is
$$D(n)=Q(n-1).$$
The implementations never materialize the underlying amoebas one by one. Instead, they count how many partial configurations can end in each admissible profile and then close the whole system with two one-dimensional recurrences.
The first structural quantity is the triangular threshold
$$T(n)=\frac{(n+1)(n+2)}{2}.$$
This is the smallest size parameter at which the level indexed by \(n\) can contribute. Consequently, for every \(n\ge 1\),
$$m<T(n)\implies A_{n,k}(m)=B_{n,k}(m)=0.$$
So when the outer loop is processing a fixed \(m\), only levels with \(T(n)\le m\) are relevant. Because \(T(n)\sim n^2/2\), the largest active level is only \(O(\sqrt m)\).
Negative sizes are impossible, so every quantity evaluated at \(m<0\) is defined to be \(0\). At level \(n=0\), the two state families do not store a full row of profile data. Instead they reduce to the same scalar stream:
$$A_{0,0}(m)=B_{0,0}(m)=P(m),$$
while
$$A_{0,k}(m)=B_{0,k}(m)=0\qquad (k\ne 0).$$
This identification is essential because the higher-level recurrences refer to \(n-1\), so the whole system must know what happens when the index reaches zero. The initial scalar condition is
$$Q(0)=1,$$
and all other missing values are supplied automatically by the zero rules for negative shifts.
For \(1\le k\le n\), define the folded index and the edge multiplier by
$$r(n,k)=\begin{cases} k,&k<n,\\ k-1,&k=n, \end{cases} \qquad \lambda(n,k)=\begin{cases} 1,&k<n,\\ 2,&k=n. \end{cases}$$
Then for every active layer \(m\ge T(n)\), the two profile families satisfy
$$\begin{aligned} A_{n,k}(m)=&\,A_{n,k}(m-n-2)+B_{n+1,1}(m-n-3)+A_{n+1,k+1}(m-n-3)\\ &+A_{n-1,r(n,k)}(m-n-1)+B_{n,1}(m-n-2)+A_{n,r(n,k)+1}(m-n-2), \end{aligned}$$
$$\begin{aligned} B_{n,k}(m)=&\,B_{n,k}(m-n-2)+B_{n+1,k+1}(m-n-3)+\lambda(n,k)\,A_{n+1,1}(m-n-3)\\ &+B_{n-1,r(n,k)}(m-n-1)+B_{n,r(n,k)+1}(m-n-2)+\lambda(n,k)\,A_{n,1}(m-n-2). \end{aligned}$$
These are exactly the transitions evaluated by the implementations. The special case \(k=n\) is the only point where the coefficient \(2\) appears, so the edge of the triangular index set has a slightly different combinatorial weight from the interior.
Once all \(A\) and \(B\) states at size \(m\) have been filled, the implementations update the scalar layer by
$$P(m)=Q(m-1)+4P(m-2)+2A_{1,1}(m-3)+B_{1,1}(m-3),$$
and, for \(m\ge 1\),
$$Q(m)=3Q(m-1)+3P(m-2).$$
The second recurrence produces the sequence whose shifted values are the final answers. Therefore
$$D(n)=Q(n-1).$$
From the dynamic-programming point of view, \(P(m)\) is the bridge between the boundary level \(n=0\) and the genuinely two-parameter profile states with \(n\ge 1\).
Every right-hand side only looks backward by offsets such as \(m-n-1\), \(m-n-2\), \(m-n-3\), \(m-1\), \(m-2\), and \(m-3\). No transition needs the full history of all earlier layers.
If
$$n_{\max}(M)=\max\{n:T(n)\le M\},$$
then all look-backs for the final target size \(M\) fit inside a circular buffer of length \(n_{\max}(M)+8\). That is why the program can keep the quadratic-time recurrence while avoiding quadratic storage in \(m\).
The initial layers can be computed by hand directly from the recurrence.
At \(m=0\), no triangular state is active, so
$$P(0)=0,\qquad Q(0)=1.$$
At \(m=1\), still no \(n\ge 1\) state survives because \(T(1)=3\). Hence
$$P(1)=Q(0)=1,\qquad Q(1)=3Q(0)=3.$$
At \(m=2\), the same argument gives
$$P(2)=Q(1)+4P(0)=3,\qquad Q(2)=3Q(1)+3P(0)=9.$$
At \(m=3\), the first triangular level appears because \(T(1)=3\). For the only index pair \((n,k)=(1,1)\), we have \(r(1,1)=0\) and \(\lambda(1,1)=2\). Using the boundary rule \(A_{0,0}(1)=B_{0,0}(1)=P(1)=1\), every invalid or negatively shifted term vanishes, so
$$A_{1,1}(3)=1,\qquad B_{1,1}(3)=1.$$
The scalar update then yields
$$P(3)=Q(2)+4P(1)=9+4=13,$$
$$Q(3)=3Q(2)+3P(1)=27+3=30.$$
Therefore
$$D(1)=1,\qquad D(2)=3,\qquad D(3)=9,\qquad D(4)=30.$$
This small example shows exactly how the first nonzero triangular state enters when the threshold is met.
The C++, Python, and Java implementations all evaluate the same recurrence. The compiled solvers allocate two triangular state tables for admissible profile pairs \((n,k)\), together with two one-dimensional arrays for the scalar layer. The outer loop advances the size parameter \(m\) from \(0\) to the target, clears the current column of the circular history buffer, fills every active profile state, and finally updates the two scalar sequences.
Before the dynamic program starts, the code computes the largest active level \(n_{\max}\) from the triangular threshold. That value determines both how many profile rows are needed and how large the circular buffer must be. Because the recurrence reaches one row above the current level, the allocation includes a small safety margin beyond \(n_{\max}\).
The C++ implementation is used in two arithmetic modes: exact integer evaluation for small verification checkpoints and modular evaluation for the large final query. The Java implementation keeps the modular variant, which is enough for the published answer. The Python implementation is a thin launcher that builds and runs the compiled solver and then extracts the final numeric output, so the three language versions share the same mathematical core.
The embedded checkpoints are
$$D(10)=44499,\qquad D(20)=9204559704,\qquad D(32)=22037102049132222,$$
and
$$D(100)\equiv 780166455\pmod{10^9}.$$
These values verify that the recurrence has been wired correctly before the computation of \(D(10000)\bmod 10^9\).
Let \(M\) be the largest size parameter reached by the program, so the final query uses \(M=n-1\). The active range satisfies \(T(n_{\max})\le M<T(n_{\max}+1)\), hence \(n_{\max}=\Theta(\sqrt M)\).
For each \(m\), the implementation visits every admissible pair \((n,k)\) with \(1\le n\le n_{\max}\) and \(1\le k\le n\). The number of profile updates per layer is therefore
$$\sum_{n=1}^{n_{\max}} n=\frac{n_{\max}(n_{\max}+1)}{2}=O(n_{\max}^2),$$
which makes the total running time
$$O(M\,n_{\max}^2)=O(M^2).$$
For memory, the circular buffer length is \(O(n_{\max})\), and the total number of stored profile entries is
$$O\!\left(n_{\max}\sum_{n=1}^{n_{\max}} n\right)=O(n_{\max}^3)=O(M^{3/2}).$$
The two scalar sequences add only \(O(n_{\max})\) more cells, so they do not change the leading asymptotic term.
Das Problem definiert eine Folge \(D(n)\), die mit der Zählung von Amoeben-Konfigurationen in einem dreidimensionalen Gitter verbunden ist. Die in den Implementierungen berechnete Zielgröße ist \(D(10000)\) modulo \(10^9\); zusätzlich werden einige kleinere Werte vorab exakt verifiziert.
Eine direkte Enumeration aller relevanten 3D-Formen wächst viel zu schnell. Deshalb ersetzt die Lösung die explizite Formerzeugung durch eine geschichtete dynamische Programmierung, deren Zustände zulässige Randprofile kodieren und deren skalare Rekurrenzen am Ende den gesuchten Folgenwert liefern.
Zur Darstellung nennen wir die beiden dreieckigen Zustandsfamilien \(A_{n,k}(m)\) und \(B_{n,k}(m)\), und die skalaren Folgen \(P(m)\) und \(Q(m)\). Die gesuchte Größe ist dann
$$D(n)=Q(n-1).$$
Die Implementierungen konstruieren die zugrunde liegenden Amoeben nicht einzeln. Stattdessen zählen sie, wie viele partielle Konfigurationen in jedem zulässigen Profil enden können, und schließen das System anschließend mit zwei eindimensionalen Rekurrenzen.
Die erste strukturelle Größe ist die Dreiecksschwelle
$$T(n)=\frac{(n+1)(n+2)}{2}.$$
Das ist der kleinste Größenparameter, bei dem die durch \(n\) indizierte Ebene überhaupt beitragen kann. Für jedes \(n\ge 1\) gilt daher
$$m<T(n)\implies A_{n,k}(m)=B_{n,k}(m)=0.$$
Wenn die äußere Schleife also einen festen Wert \(m\) verarbeitet, sind nur Ebenen mit \(T(n)\le m\) relevant. Da \(T(n)\sim n^2/2\) gilt, ist das größte aktive \(n\) nur von Ordnung \(O(\sqrt m)\).
Negative Größen sind unmöglich, also wird jede Größe für \(m<0\) als \(0\) definiert. Auf der Ebene \(n=0\) speichern die beiden Zustandsfamilien keine volle Profilzeile, sondern fallen auf denselben skalaren Strom zurück:
$$A_{0,0}(m)=B_{0,0}(m)=P(m),$$
während
$$A_{0,k}(m)=B_{0,k}(m)=0\qquad (k\ne 0).$$
Diese Identifikation ist entscheidend, weil die Rekurrenzen der höheren Ebenen auf \(n-1\) verweisen und das System daher wissen muss, was beim Erreichen der Null-Ebene passiert. Die anfängliche skalare Bedingung lautet
$$Q(0)=1,$$
und alle übrigen fehlenden Werte entstehen automatisch aus den Nullregeln für negative Verschiebungen.
Für \(1\le k\le n\) definieren wir den gefalteten Index und den Randmultiplikator durch
$$r(n,k)=\begin{cases} k,&k<n,\\ k-1,&k=n, \end{cases} \qquad \lambda(n,k)=\begin{cases} 1,&k<n,\\ 2,&k=n. \end{cases}$$
Dann gelten für jede aktive Ebene \(m\ge T(n)\) die beiden Profilrekurrenzen
$$\begin{aligned} A_{n,k}(m)=&\,A_{n,k}(m-n-2)+B_{n+1,1}(m-n-3)+A_{n+1,k+1}(m-n-3)\\ &+A_{n-1,r(n,k)}(m-n-1)+B_{n,1}(m-n-2)+A_{n,r(n,k)+1}(m-n-2), \end{aligned}$$
$$\begin{aligned} B_{n,k}(m)=&\,B_{n,k}(m-n-2)+B_{n+1,k+1}(m-n-3)+\lambda(n,k)\,A_{n+1,1}(m-n-3)\\ &+B_{n-1,r(n,k)}(m-n-1)+B_{n,r(n,k)+1}(m-n-2)+\lambda(n,k)\,A_{n,1}(m-n-2). \end{aligned}$$
Genau diese Übergänge werden von den Implementierungen ausgewertet. Der Spezialfall \(k=n\) ist die einzige Stelle, an der der Koeffizient \(2\) auftaucht; der Rand des dreieckigen Indexbereichs erhält also ein anderes kombinatorisches Gewicht als das Innere.
Sobald alle \(A\)- und \(B\)-Zustände für die Größe \(m\) berechnet sind, aktualisieren die Implementierungen die skalare Ebene mittels
$$P(m)=Q(m-1)+4P(m-2)+2A_{1,1}(m-3)+B_{1,1}(m-3),$$
und für \(m\ge 1\) gilt zusätzlich
$$Q(m)=3Q(m-1)+3P(m-2).$$
Die zweite Rekurrenz erzeugt die Folge, deren verschobene Werte die endgültigen Antworten liefern. Also ist
$$D(n)=Q(n-1).$$
Aus Sicht der dynamischen Programmierung ist \(P(m)\) die Brücke zwischen der Randebene \(n=0\) und den echten Zweiparameter-Zuständen mit \(n\ge 1\).
Jede rechte Seite schaut nur auf Rücksprünge wie \(m-n-1\), \(m-n-2\), \(m-n-3\), \(m-1\), \(m-2\) und \(m-3\). Kein Übergang benötigt die vollständige Historie aller früheren Ebenen.
Ist
$$n_{\max}(M)=\max\{n:T(n)\le M\},$$
dann passen für die Zielgröße \(M\) alle Rücksprünge in einen Ringpuffer der Länge \(n_{\max}(M)+8\). Deshalb kann das Programm die quadratische Laufzeit der Rekurrenz beibehalten, ohne zugleich quadratischen Speicher in \(m\) zu benötigen.
Die Anfangsschichten lassen sich direkt aus der Rekurrenz von Hand berechnen.
Bei \(m=0\) ist noch kein dreieckiger Zustand aktiv, also
$$P(0)=0,\qquad Q(0)=1.$$
Bei \(m=1\) überlebt ebenfalls kein Zustand mit \(n\ge 1\), denn \(T(1)=3\). Daher
$$P(1)=Q(0)=1,\qquad Q(1)=3Q(0)=3.$$
Bei \(m=2\) ergibt dasselbe Argument
$$P(2)=Q(1)+4P(0)=3,\qquad Q(2)=3Q(1)+3P(0)=9.$$
Bei \(m=3\) erscheint die erste dreieckige Ebene, weil \(T(1)=3\). Für das einzige Indexpaar \((n,k)=(1,1)\) gilt \(r(1,1)=0\) und \(\lambda(1,1)=2\). Mit der Randbedingung \(A_{0,0}(1)=B_{0,0}(1)=P(1)=1\) verschwinden alle unzulässigen oder negativ verschobenen Terme, sodass
$$A_{1,1}(3)=1,\qquad B_{1,1}(3)=1.$$
Das skalare Update liefert dann
$$P(3)=Q(2)+4P(1)=9+4=13,$$
$$Q(3)=3Q(2)+3P(1)=27+3=30.$$
Folglich gilt
$$D(1)=1,\qquad D(2)=3,\qquad D(3)=9,\qquad D(4)=30.$$
Dieses kleine Beispiel zeigt genau, wie der erste nichtverschwindende dreieckige Zustand in dem Moment eintritt, in dem die Schwelle erreicht wird.
Die C++-, Python- und Java-Implementierungen werten dieselbe Rekurrenz aus. Die kompilierten Löser reservieren zwei dreieckige Zustandstabellen für zulässige Profilpaare \((n,k)\) sowie zwei eindimensionale Arrays für die skalare Ebene. Die äußere Schleife erhöht den Größenparameter \(m\) von \(0\) bis zum Zielwert, leert die aktuelle Spalte des Ringpuffers, füllt jeden aktiven Profilzustand und aktualisiert erst danach die beiden skalaren Folgen.
Vor dem Start der dynamischen Programmierung bestimmt der Code aus der Dreiecksschwelle die größte aktive Ebene \(n_{\max}\). Dieser Wert legt sowohl die benötigte Anzahl von Profilzeilen als auch die Länge des Ringpuffers fest. Weil die Rekurrenz eine Zeile oberhalb der aktuellen Ebene anspricht, enthält die Allokation noch einen kleinen Sicherheitsrand über \(n_{\max}\).
Die C++-Implementierung wird in zwei Rechenmodi verwendet: exakte Ganzzahlauswertung für kleine Kontrollwerte und modulare Auswertung für die große Schlussanfrage. Die Java-Implementierung behält nur die modulare Variante bei, was für die veröffentlichte Antwort ausreicht. Die Python-Implementierung ist ein schlanker Starter, der den kompilierten Löser baut, ausführt und den numerischen Endwert extrahiert; dadurch teilen alle drei Sprachfassungen denselben mathematischen Kern.
Die eingebetteten Kontrollwerte sind
$$D(10)=44499,\qquad D(20)=9204559704,\qquad D(32)=22037102049132222,$$
und
$$D(100)\equiv 780166455\pmod{10^9}.$$
Damit wird verifiziert, dass die Rekurrenz korrekt verdrahtet ist, bevor \(D(10000)\bmod 10^9\) berechnet wird.
Sei \(M\) der größte vom Programm berechnete Größenparameter, also für die Endanfrage \(M=n-1\). Der aktive Bereich erfüllt \(T(n_{\max})\le M<T(n_{\max}+1)\), also gilt \(n_{\max}=\Theta(\sqrt M)\).
Für jedes \(m\) besucht die Implementierung jedes zulässige Paar \((n,k)\) mit \(1\le n\le n_{\max}\) und \(1\le k\le n\). Die Anzahl der Profil-Updates pro Ebene ist daher
$$\sum_{n=1}^{n_{\max}} n=\frac{n_{\max}(n_{\max}+1)}{2}=O(n_{\max}^2),$$
woraus sich für die Gesamtlaufzeit
$$O(M\,n_{\max}^2)=O(M^2)$$
ergibt. Für den Speicher ist die Ringpufferlänge \(O(n_{\max})\), und die Gesamtzahl der gespeicherten Profileinträge beträgt
$$O\!\left(n_{\max}\sum_{n=1}^{n_{\max}} n\right)=O(n_{\max}^3)=O(M^{3/2}).$$
Die beiden skalaren Folgen benötigen nur \(O(n_{\max})\) zusätzliche Zellen und ändern damit den führenden asymptotischen Term nicht.
Problem, üç boyutlu bir ızgaradaki amip konfigürasyonlarını sayan \(D(n)\) dizisini tanımlar. Uygulamalarda hedeflenen hesaplama \(D(10000)\) değerinin \(10^9\) moduna göre bulunmasıdır; büyük modlu çalıştırmadan önce daha küçük bazı değerler de tam olarak doğrulanır.
İlgili 3B şekillerin doğrudan tek tek üretilmesi çok hızlı büyür. Bu yüzden çözüm, açık şekil üretimini bırakıp geçerli sınır profillerini tutan katmanlı bir dinamik program kurar; son dizi değeri de iki skaler bağıntı üzerinden elde edilir.
Anlatım için iki üçgensel durum ailesini \(A_{n,k}(m)\) ve \(B_{n,k}(m)\), skaler dizileri de \(P(m)\) ve \(Q(m)\) ile gösterelim. Aranan nicelik
$$D(n)=Q(n-1)$$
şeklindedir. Uygulamalar amipleri tek tek kurmaz; bunun yerine her geçerli profilde kaç kısmi konfigürasyonun bitebileceğini sayar ve sistemi iki tek boyutlu rekürans ile kapatır.
İlk yapısal büyüklük şu üçgensel eştir:
$$T(n)=\frac{(n+1)(n+2)}{2}.$$
Bu, \(n\) ile indekslenen katmanın katkı verebileceği en küçük boyut parametresidir. Dolayısıyla her \(n\ge 1\) için
$$m<T(n)\implies A_{n,k}(m)=B_{n,k}(m)=0$$
olur. Yani dış döngü sabit bir \(m\) değerini işlerken yalnızca \(T(n)\le m\) olan katmanlar önemlidir. \(T(n)\sim n^2/2\) olduğundan en büyük etkin katman ancak \(O(\sqrt m)\) mertebesindedir.
Negatif boyutlar imkansızdır; bu yüzden \(m<0\) için değerlendirilen her nicelik \(0\) kabul edilir. \(n=0\) seviyesinde iki durum ailesi tam bir profil satırı tutmaz; ikisi de aynı skaler akışa indirgenir:
$$A_{0,0}(m)=B_{0,0}(m)=P(m),$$
ve
$$A_{0,k}(m)=B_{0,k}(m)=0\qquad (k\ne 0).$$
Bu özdeşleştirme kritiktir; çünkü üst seviyelerin reküransları \(n-1\) katmanına geri döner ve sistemin sıfır seviyesine ulaşıldığında ne olacağını bilmesi gerekir. Başlangıç skaler koşulu
$$Q(0)=1$$
olup geri kalan eksik değerler negatif kaymalar için sıfır kurallarından otomatik gelir.
\(1\le k\le n\) için katlanmış indeks ve kenar katsayısını şöyle tanımlayalım:
$$r(n,k)=\begin{cases} k,&k<n,\\ k-1,&k=n, \end{cases} \qquad \lambda(n,k)=\begin{cases} 1,&k<n,\\ 2,&k=n. \end{cases}$$
O zaman her etkin katman \(m\ge T(n)\) için iki profil ailesi şu bağıntıları sağlar:
$$\begin{aligned} A_{n,k}(m)=&\,A_{n,k}(m-n-2)+B_{n+1,1}(m-n-3)+A_{n+1,k+1}(m-n-3)\\ &+A_{n-1,r(n,k)}(m-n-1)+B_{n,1}(m-n-2)+A_{n,r(n,k)+1}(m-n-2), \end{aligned}$$
$$\begin{aligned} B_{n,k}(m)=&\,B_{n,k}(m-n-2)+B_{n+1,k+1}(m-n-3)+\lambda(n,k)\,A_{n+1,1}(m-n-3)\\ &+B_{n-1,r(n,k)}(m-n-1)+B_{n,r(n,k)+1}(m-n-2)+\lambda(n,k)\,A_{n,1}(m-n-2). \end{aligned}$$
Uygulamaların hesapladığı geçişler tam olarak bunlardır. \(k=n\) özel durumu, katsayının \(2\) olduğu tek yerdir; yani üçgensel indeks kümesinin kenarı, iç bölgeden biraz farklı bir kombinatorik ağırlık taşır.
Belirli bir \(m\) için tüm \(A\) ve \(B\) durumları hesaplandıktan sonra uygulamalar skaler katmanı
$$P(m)=Q(m-1)+4P(m-2)+2A_{1,1}(m-3)+B_{1,1}(m-3)$$
ve \(m\ge 1\) için
$$Q(m)=3Q(m-1)+3P(m-2)$$
bağıntılarıyla günceller. Nihai cevaplar, ikinci dizinin kaydırılmış değerlerinden gelir; dolayısıyla
$$D(n)=Q(n-1).$$
Dinamik programlama açısından bakıldığında \(P(m)\), sınır seviyesi \(n=0\) ile \(n\ge 1\) olan gerçek iki parametreli profil durumları arasındaki köprüdür.
Sağ taraftaki her terim yalnızca \(m-n-1\), \(m-n-2\), \(m-n-3\), \(m-1\), \(m-2\) ve \(m-3\) gibi geri kaymalara bakar. Geçmişin tamamını saklamak gerekmez.
Eğer
$$n_{\max}(M)=\max\{n:T(n)\le M\}$$
ise, hedef boyut \(M\) için gereken tüm geri bakışlar uzunluğu \(n_{\max}(M)+8\) olan bir döngüsel tampon içine sığar. Böylece program, \(m\) boyunca kuadratik depolama kullanmadan kuadratik zamanlı reküransı çalıştırabilir.
Başlangıç katmanları doğrudan reküranstan elle hesaplanabilir.
\(m=0\) iken henüz üçgensel hiçbir durum etkin değildir; dolayısıyla
$$P(0)=0,\qquad Q(0)=1.$$
\(m=1\) için de \(T(1)=3\) olduğundan \(n\ge 1\) durumları yoktur. Bu yüzden
$$P(1)=Q(0)=1,\qquad Q(1)=3Q(0)=3.$$
\(m=2\) için aynı gerekçe
$$P(2)=Q(1)+4P(0)=3,\qquad Q(2)=3Q(1)+3P(0)=9$$
sonucunu verir.
\(m=3\) olduğunda ilk üçgensel katman açılır, çünkü \(T(1)=3\). Tek indeks çifti \((n,k)=(1,1)\) için \(r(1,1)=0\) ve \(\lambda(1,1)=2\) olur. Sınır kuralı \(A_{0,0}(1)=B_{0,0}(1)=P(1)=1\) kullanıldığında geçersiz veya negatif kaymalı tüm terimler sıfırlanır ve
$$A_{1,1}(3)=1,\qquad B_{1,1}(3)=1$$
elde edilir. Sonra skaler güncelleme
$$P(3)=Q(2)+4P(1)=9+4=13,$$
$$Q(3)=3Q(2)+3P(1)=27+3=30$$
olur. Böylece
$$D(1)=1,\qquad D(2)=3,\qquad D(3)=9,\qquad D(4)=30.$$
Bu küçük örnek, eşiğin sağlandığı anda ilk üçgensel durumun nasıl devreye girdiğini açık biçimde gösterir.
C++, Python ve Java uygulamalarının hepsi aynı reküransı değerlendirir. Derlenen çözücüler, geçerli \((n,k)\) profil çiftleri için iki üçgensel durum tablosu ve skaler katman için iki tek boyutlu dizi ayırır. Dış döngü, boyut parametresi \(m\)'yi \(0\)'dan hedefe kadar ilerletir; önce döngüsel geçmiş tamponunun mevcut sütununu temizler, sonra tüm etkin profil durumlarını doldurur ve en sonda iki skaler diziyi günceller.
Dinamik program başlamadan önce kod, üçgensel eşikten en büyük etkin katman \(n_{\max}\) değerini çıkarır. Bu değer hem kaç profil satırına ihtiyaç olduğunu hem de döngüsel tamponun uzunluğunu belirler. Rekürans mevcut katmanın bir üstünü de kullandığı için ayırma işlemi \(n_{\max}\)'in biraz üstüne güvenlik payı ekler.
C++ uygulaması iki aritmetik kipte kullanılır: küçük kontrol değerleri için tam sayı hesabı ve büyük nihai sorgu için modlu hesap. Java uygulaması yalnızca modlu sürümü tutar; yayımlanan cevap için bu yeterlidir. Python uygulaması ise derlenen çözücüyü oluşturup çalıştıran ve son sayısal çıktıyı ayıklayan ince bir sarmaldır; böylece üç dil sürümü de aynı matematiksel çekirdeği paylaşır.
Uygulamalara gömülü kontrol değerleri
$$D(10)=44499,\qquad D(20)=9204559704,\qquad D(32)=22037102049132222,$$
ve
$$D(100)\equiv 780166455\pmod{10^9}$$
şeklindedir. Bunlar, \(D(10000)\bmod 10^9\) hesaplanmadan önce reküransın doğru kurulduğunu doğrular.
\(M\) programın ulaştığı en büyük boyut parametresi olsun; nihai sorguda \(M=n-1\) kullanılır. Etkin aralık \(T(n_{\max})\le M<T(n_{\max}+1)\) koşulunu sağlar; dolayısıyla \(n_{\max}=\Theta(\sqrt M)\).
Her \(m\) için uygulama, \(1\le n\le n_{\max}\) ve \(1\le k\le n\) olan her geçerli \((n,k)\) çiftini ziyaret eder. Katman başına profil güncelleme sayısı
$$\sum_{n=1}^{n_{\max}} n=\frac{n_{\max}(n_{\max}+1)}{2}=O(n_{\max}^2)$$
olduğundan toplam çalışma süresi
$$O(M\,n_{\max}^2)=O(M^2)$$
olur. Bellekte ise döngüsel tampon uzunluğu \(O(n_{\max})\), saklanan toplam profil girdisi sayısı da
$$O\!\left(n_{\max}\sum_{n=1}^{n_{\max}} n\right)=O(n_{\max}^3)=O(M^{3/2})$$
mertebesindedir. İki skaler dizi yalnızca \(O(n_{\max})\) ek hücre kullandığı için baskın terimi değiştirmez.
El problema define una sucesión \(D(n)\) asociada al conteo de configuraciones de amebas en una cuadrícula tridimensional. El objetivo computacional usado por las implementaciones es evaluar \(D(10000)\) módulo \(10^9\), verificando antes varios valores más pequeños de manera exacta.
Enumerar de forma directa todas las formas 3D relevantes crece demasiado rápido. Por eso la solución sustituye la generación explícita de figuras por una programación dinámica por capas, donde los estados codifican perfiles de frontera admisibles y las recurrencias escalares extraen al final el valor buscado de la sucesión.
Para la exposición, llamemos \(A_{n,k}(m)\) y \(B_{n,k}(m)\) a las dos familias triangulares de estados, y \(P(m)\), \(Q(m)\) a las sucesiones escalares. La cantidad deseada es
$$D(n)=Q(n-1).$$
Las implementaciones no construyen una por una las amebas subyacentes. En su lugar cuentan cuántas configuraciones parciales pueden terminar en cada perfil permitido y cierran después el sistema con dos recurrencias unidimensionales.
La primera magnitud estructural es el umbral triangular
$$T(n)=\frac{(n+1)(n+2)}{2}.$$
Este es el menor parámetro de tamaño en el que la capa indexada por \(n\) puede contribuir. En consecuencia, para todo \(n\ge 1\),
$$m<T(n)\implies A_{n,k}(m)=B_{n,k}(m)=0.$$
Así, cuando el bucle exterior procesa un valor fijo de \(m\), solo importan las capas con \(T(n)\le m\). Como \(T(n)\sim n^2/2\), el mayor nivel activo es solo \(O(\sqrt m)\).
Los tamaños negativos son imposibles, de modo que toda cantidad evaluada en \(m<0\) se define como \(0\). En el nivel \(n=0\), las dos familias de estados no guardan una fila completa de perfiles, sino que ambas se reducen a la misma sucesión escalar:
$$A_{0,0}(m)=B_{0,0}(m)=P(m),$$
mientras que
$$A_{0,k}(m)=B_{0,k}(m)=0\qquad (k\ne 0).$$
Esta identificación es esencial porque las recurrencias de niveles superiores retroceden hasta \(n-1\), así que el sistema necesita saber qué ocurre al llegar al nivel cero. La condición escalar inicial es
$$Q(0)=1,$$
y todos los demás valores ausentes aparecen automáticamente gracias a las reglas de anulación para desplazamientos negativos.
Para \(1\le k\le n\), definimos el índice plegado y el multiplicador de borde mediante
$$r(n,k)=\begin{cases} k,&k<n,\\ k-1,&k=n, \end{cases} \qquad \lambda(n,k)=\begin{cases} 1,&k<n,\\ 2,&k=n. \end{cases}$$
Entonces, para toda capa activa \(m\ge T(n)\), las dos familias de perfiles cumplen
$$\begin{aligned} A_{n,k}(m)=&\,A_{n,k}(m-n-2)+B_{n+1,1}(m-n-3)+A_{n+1,k+1}(m-n-3)\\ &+A_{n-1,r(n,k)}(m-n-1)+B_{n,1}(m-n-2)+A_{n,r(n,k)+1}(m-n-2), \end{aligned}$$
$$\begin{aligned} B_{n,k}(m)=&\,B_{n,k}(m-n-2)+B_{n+1,k+1}(m-n-3)+\lambda(n,k)\,A_{n+1,1}(m-n-3)\\ &+B_{n-1,r(n,k)}(m-n-1)+B_{n,r(n,k)+1}(m-n-2)+\lambda(n,k)\,A_{n,1}(m-n-2). \end{aligned}$$
Estas son exactamente las transiciones que evalúan las implementaciones. El caso especial \(k=n\) es el único punto donde aparece el coeficiente \(2\), así que el borde del conjunto triangular de índices tiene un peso combinatorio ligeramente distinto del interior.
Una vez calculados todos los estados \(A\) y \(B\) para un tamaño \(m\), las implementaciones actualizan la capa escalar mediante
$$P(m)=Q(m-1)+4P(m-2)+2A_{1,1}(m-3)+B_{1,1}(m-3),$$
y, para \(m\ge 1\),
$$Q(m)=3Q(m-1)+3P(m-2).$$
La segunda recurrencia produce la sucesión cuyas versiones desplazadas son las respuestas finales. Por tanto,
$$D(n)=Q(n-1).$$
Desde el punto de vista de la programación dinámica, \(P(m)\) es el puente entre el nivel de borde \(n=0\) y los estados de perfil genuinamente bidimensionales con \(n\ge 1\).
Cada lado derecho solo mira hacia atrás a desplazamientos como \(m-n-1\), \(m-n-2\), \(m-n-3\), \(m-1\), \(m-2\) y \(m-3\). Ninguna transición necesita la historia completa de todas las capas anteriores.
Si
$$n_{\max}(M)=\max\{n:T(n)\le M\},$$
entonces todos los retrocesos necesarios para el tamaño final \(M\) caben en un búfer circular de longitud \(n_{\max}(M)+8\). Por eso el programa puede mantener la recurrencia cuadrática en tiempo sin necesitar también memoria cuadrática en \(m\).
Las capas iniciales pueden calcularse a mano directamente a partir de la recurrencia.
En \(m=0\) todavía no hay estados triangulares activos, así que
$$P(0)=0,\qquad Q(0)=1.$$
En \(m=1\) sigue sin sobrevivir ningún estado con \(n\ge 1\) porque \(T(1)=3\). Luego
$$P(1)=Q(0)=1,\qquad Q(1)=3Q(0)=3.$$
En \(m=2\), el mismo razonamiento da
$$P(2)=Q(1)+4P(0)=3,\qquad Q(2)=3Q(1)+3P(0)=9.$$
En \(m=3\) aparece el primer nivel triangular porque \(T(1)=3\). Para el único par \((n,k)=(1,1)\), tenemos \(r(1,1)=0\) y \(\lambda(1,1)=2\). Usando la condición de borde \(A_{0,0}(1)=B_{0,0}(1)=P(1)=1\), todos los términos inválidos o con desplazamiento negativo desaparecen, de modo que
$$A_{1,1}(3)=1,\qquad B_{1,1}(3)=1.$$
La actualización escalar produce entonces
$$P(3)=Q(2)+4P(1)=9+4=13,$$
$$Q(3)=3Q(2)+3P(1)=27+3=30.$$
Por lo tanto,
$$D(1)=1,\qquad D(2)=3,\qquad D(3)=9,\qquad D(4)=30.$$
Este ejemplo pequeño muestra exactamente cómo entra el primer estado triangular no nulo cuando se alcanza el umbral.
Las implementaciones en C++, Python y Java evalúan la misma recurrencia. Los solucionadores compilados reservan dos tablas triangulares de estados para los pares de perfil admisibles \((n,k)\), junto con dos arreglos unidimensionales para la capa escalar. El bucle exterior hace avanzar el parámetro de tamaño \(m\) desde \(0\) hasta el objetivo, limpia la columna actual del búfer circular, rellena todos los estados de perfil activos y solo después actualiza las dos sucesiones escalares.
Antes de comenzar la programación dinámica, el código calcula el mayor nivel activo \(n_{\max}\) a partir del umbral triangular. Ese valor determina tanto cuántas filas de perfil se necesitan como el tamaño del búfer circular. Como la recurrencia consulta una fila por encima del nivel actual, la reserva incluye un pequeño margen de seguridad por encima de \(n_{\max}\).
La implementación en C++ se usa en dos modos aritméticos: evaluación exacta con enteros para puntos de control pequeños y evaluación modular para la gran consulta final. La implementación en Java conserva la variante modular, suficiente para la respuesta publicada. La implementación en Python es un lanzador ligero que construye y ejecuta el solucionador compilado y extrae la salida numérica final; así, las tres versiones comparten el mismo núcleo matemático.
Los puntos de control incorporados son
$$D(10)=44499,\qquad D(20)=9204559704,\qquad D(32)=22037102049132222,$$
y
$$D(100)\equiv 780166455\pmod{10^9}.$$
Estos valores verifican que la recurrencia está conectada correctamente antes de calcular \(D(10000)\bmod 10^9\).
Sea \(M\) el mayor parámetro de tamaño alcanzado por el programa; en la consulta final se usa \(M=n-1\). El rango activo cumple \(T(n_{\max})\le M<T(n_{\max}+1)\), así que \(n_{\max}=\Theta(\sqrt M)\).
Para cada \(m\), la implementación visita cada par admisible \((n,k)\) con \(1\le n\le n_{\max}\) y \(1\le k\le n\). El número de actualizaciones de perfil por capa es
$$\sum_{n=1}^{n_{\max}} n=\frac{n_{\max}(n_{\max}+1)}{2}=O(n_{\max}^2),$$
de donde el tiempo total resulta
$$O(M\,n_{\max}^2)=O(M^2).$$
En memoria, la longitud del búfer circular es \(O(n_{\max})\), y el número total de entradas de perfil almacenadas es
$$O\!\left(n_{\max}\sum_{n=1}^{n_{\max}} n\right)=O(n_{\max}^3)=O(M^{3/2}).$$
Las dos sucesiones escalares añaden solo \(O(n_{\max})\) celdas extra, así que no cambian el término asintótico dominante.
这道题定义了一个与三维网格中“amoeba”计数相关的序列 \(D(n)\)。从实现可以看出，实际计算目标是求 \(D(10000)\bmod 10^9\)，并在进入大规模模运算之前先对若干较小值做精确校验。
如果直接枚举所有相关的三维形状，状态数会膨胀得非常快，因此不可行。实现采用的做法是把显式的形状生成改写为分层动态规划：状态只记录“可行边界轮廓”的计数，再通过两个一维递推把这些局部状态汇总成最终的 \(D(n)\)。
为了叙述方便，把两族三角形状态记为 \(A_{n,k}(m)\) 和 \(B_{n,k}(m)\)，把两个标量序列记为 \(P(m)\) 与 \(Q(m)\)。最终答案写成
$$D(n)=Q(n-1).$$
整个算法从头到尾都没有把具体 amoeba 逐个构造出来，而是在问：大小参数为 \(m\) 时，有多少部分构型能落在某个允许的轮廓状态上。等所有轮廓状态都算完以后，再由标量递推收束成真正需要的序列值。
首先定义三角阈值
$$T(n)=\frac{(n+1)(n+2)}{2}.$$
这表示带有层指标 \(n\) 的状态，最早要到多大的 \(m\) 才有可能非零。因此对所有 \(n\ge 1\)，都有
$$m<T(n)\implies A_{n,k}(m)=B_{n,k}(m)=0.$$
这条规则非常重要，因为它说明当外层循环固定在某个 \(m\) 时，只有满足 \(T(n)\le m\) 的层才需要处理。由于 \(T(n)\sim n^2/2\)，所以真正活跃的最大层数只有 \(O(\sqrt m)\) 级别，而不是线性级别。
负的大小参数没有意义，所以凡是出现在 \(m<0\) 处的量都定义为 \(0\)。对于 \(n=0\) 这一层，两族状态并不保存一整行轮廓，而是统一退化为同一个标量流：
$$A_{0,0}(m)=B_{0,0}(m)=P(m),$$
并且
$$A_{0,k}(m)=B_{0,k}(m)=0\qquad (k\ne 0).$$
这个边界识别是整个系统能闭合的关键，因为高层递推里会访问 \(n-1\) 的状态；当索引降到零时，程序必须知道要接到什么对象上。初始化方面，唯一显式给出的标量初值是
$$Q(0)=1,$$
其余缺失值都由“负下标即零”的规则自动补齐。
对于 \(1\le k\le n\)，定义折叠后的索引和边界倍数
$$r(n,k)=\begin{cases} k,&k<n,\\ k-1,&k=n, \end{cases} \qquad \lambda(n,k)=\begin{cases} 1,&k<n,\\ 2,&k=n. \end{cases}$$
那么当 \(m\ge T(n)\) 时，两族状态满足
$$\begin{aligned} A_{n,k}(m)=&\,A_{n,k}(m-n-2)+B_{n+1,1}(m-n-3)+A_{n+1,k+1}(m-n-3)\\ &+A_{n-1,r(n,k)}(m-n-1)+B_{n,1}(m-n-2)+A_{n,r(n,k)+1}(m-n-2), \end{aligned}$$
$$\begin{aligned} B_{n,k}(m)=&\,B_{n,k}(m-n-2)+B_{n+1,k+1}(m-n-3)+\lambda(n,k)\,A_{n+1,1}(m-n-3)\\ &+B_{n-1,r(n,k)}(m-n-1)+B_{n,r(n,k)+1}(m-n-2)+\lambda(n,k)\,A_{n,1}(m-n-2). \end{aligned}$$
这两条式子就是实现逐项计算的核心。只有在边界位置 \(k=n\) 时，系数才从 \(1\) 变成 \(2\)，说明三角索引区域的边缘与内部在组合计数上并不完全相同。
在给定的 \(m\) 上，所有 \(A\) 与 \(B\) 状态都填完之后，实现继续计算标量层：
$$P(m)=Q(m-1)+4P(m-2)+2A_{1,1}(m-3)+B_{1,1}(m-3),$$
并且当 \(m\ge 1\) 时，
$$Q(m)=3Q(m-1)+3P(m-2).$$
第二条递推直接生成目标序列，只是答案要向前平移一位读取，因此
$$D(n)=Q(n-1).$$
从动态规划的结构来看，\(P(m)\) 扮演了“边界基层”和“真正二维状态表”之间的桥梁角色：一方面它接收来自 \(A_{1,1}\)、\(B_{1,1}\) 的贡献，另一方面又作为 \(n=0\) 的公共底层返回给高层递推使用。
观察所有右端项，可以看到它们只会访问有限的历史位置，例如 \(m-n-1\)、\(m-n-2\)、\(m-n-3\)、\(m-1\)、\(m-2\)、\(m-3\)。这意味着程序根本不需要保存从 \(0\) 到当前 \(m\) 的完整历史。
若定义
$$n_{\max}(M)=\max\{n:T(n)\le M\},$$
那么对于最终目标大小 \(M\)，实现中的所有回看距离都可以装进长度为 \(n_{\max}(M)+8\) 的环形缓冲区中。也就是说，时间复杂度虽然仍然是二次级别，但在大小参数这一维上并不需要二次内存。
这套递推的前几层完全可以手算出来。
当 \(m=0\) 时，还没有任何三角层激活，因此
$$P(0)=0,\qquad Q(0)=1.$$
当 \(m=1\) 时，由于 \(T(1)=3\)，所有 \(n\ge 1\) 的状态依然为零，所以
$$P(1)=Q(0)=1,\qquad Q(1)=3Q(0)=3.$$
当 \(m=2\) 时，情况仍然一样，于是
$$P(2)=Q(1)+4P(0)=3,\qquad Q(2)=3Q(1)+3P(0)=9.$$
当 \(m=3\) 时，第一层终于出现，因为 \(T(1)=3\)。此时唯一可能的指标对是 \((n,k)=(1,1)\)，并且 \(r(1,1)=0\)、\(\lambda(1,1)=2\)。再结合边界条件 \(A_{0,0}(1)=B_{0,0}(1)=P(1)=1\)，所有非法项和负位移项都会消失，因此
$$A_{1,1}(3)=1,\qquad B_{1,1}(3)=1.$$
接着标量层更新为
$$P(3)=Q(2)+4P(1)=9+4=13,$$
$$Q(3)=3Q(2)+3P(1)=27+3=30.$$
于是得到
$$D(1)=1,\qquad D(2)=3,\qquad D(3)=9,\qquad D(4)=30.$$
这个小例子清楚展示了：只有当三角阈值真正达到以后，第一批非零的轮廓状态才会进入系统。
C++、Python 和 Java 三个实现使用的是完全相同的递推体系。编译型实现会为所有允许的 \((n,k)\) 分配两张三角形状态表，同时再分配两个一维数组保存标量层。外层循环让大小参数 \(m\) 从 \(0\) 递增到目标值；每走一步，先清空环形历史缓冲区中当前这一列，再填充所有活跃轮廓状态，最后更新两条标量序列。
在动态规划开始之前，代码会根据三角阈值求出最大活跃层 \(n_{\max}\)。这个值既决定需要多少行轮廓状态，也决定环形缓冲区要开多长。由于递推会访问当前层之上的一层，实际分配时还会在 \(n_{\max}\) 之上留出一小段安全余量。
C++ 实现支持两种算术模式：对小型校验点使用精确整数运算，对最终大规模问题使用模 \(10^9\) 运算。Java 实现保留的是模运算版本，这已经足够得到发布答案。Python 实现则是一个轻量级启动器，负责构建并运行编译后的求解器，再从输出中提取最终数字，因此三种语言在数学核心上完全一致。
实现内嵌的校验值为
$$D(10)=44499,\qquad D(20)=9204559704,\qquad D(32)=22037102049132222,$$
以及
$$D(100)\equiv 780166455\pmod{10^9}.$$
这些结果说明，在真正计算 \(D(10000)\bmod 10^9\) 之前，递推结构已经通过了较强的一致性检查。
设程序需要计算到的最大大小参数为 \(M\)，那么最终查询对应的是 \(M=n-1\)。活跃层范围满足 \(T(n_{\max})\le M<T(n_{\max}+1)\)，因此 \(n_{\max}=\Theta(\sqrt M)\)。
对于每个 \(m\)，实现都要访问所有满足 \(1\le n\le n_{\max}\)、\(1\le k\le n\) 的状态对 \((n,k)\)。每一层的状态更新次数为
$$\sum_{n=1}^{n_{\max}} n=\frac{n_{\max}(n_{\max}+1)}{2}=O(n_{\max}^2),$$
所以总时间复杂度是
$$O(M\,n_{\max}^2)=O(M^2).$$
空间方面，环形缓冲区的长度是 \(O(n_{\max})\)，而总共保存的轮廓状态项数为
$$O\!\left(n_{\max}\sum_{n=1}^{n_{\max}} n\right)=O(n_{\max}^3)=O(M^{3/2}).$$
两条标量序列只额外增加 \(O(n_{\max})\) 个单元，因此不会改变主导阶。
В задаче определяется последовательность \(D(n)\), связанная с подсчетом конфигураций амебы в трехмерной решетке. Из реализаций видно, что вычислительная цель состоит в нахождении \(D(10000)\) по модулю \(10^9\), причем несколько меньших значений предварительно проверяются точно.
Прямой перебор всех нужных трехмерных форм растет слишком быстро. Поэтому решение заменяет явную генерацию фигур послойной динамикой, в которой состояния кодируют допустимые граничные профили, а две скалярные рекурсии в конце извлекают нужное значение последовательности.
Для изложения обозначим две треугольные семьи состояний через \(A_{n,k}(m)\) и \(B_{n,k}(m)\), а две скалярные последовательности через \(P(m)\) и \(Q(m)\). Тогда искомая величина равна
$$D(n)=Q(n-1).$$
Реализации не строят амебы по одной. Вместо этого они считают, сколько частичных конфигураций может закончиться в каждом допустимом профиле, а затем замыкают систему двумя одномерными рекурсиями.
Первой структурной величиной является треугольный порог
$$T(n)=\frac{(n+1)(n+2)}{2}.$$
Это минимальный параметр размера, при котором слой с индексом \(n\) может дать ненулевой вклад. Следовательно, для любого \(n\ge 1\)
$$m<T(n)\implies A_{n,k}(m)=B_{n,k}(m)=0.$$
Значит, когда внешний цикл обрабатывает фиксированное \(m\), нужно рассматривать только те уровни, для которых \(T(n)\le m\). Поскольку \(T(n)\sim n^2/2\), максимально активный уровень имеет лишь порядок \(O(\sqrt m)\).
Отрицательные размеры невозможны, поэтому любая величина при \(m<0\) определяется как \(0\). На уровне \(n=0\) две семьи состояний не хранят полную строку профилей, а обе сводятся к одному и тому же скалярному потоку:
$$A_{0,0}(m)=B_{0,0}(m)=P(m),$$
а
$$A_{0,k}(m)=B_{0,k}(m)=0\qquad (k\ne 0).$$
Это отождествление принципиально важно, потому что рекурсии верхних уровней обращаются к \(n-1\), и система должна знать, что происходит при достижении нулевого уровня. Начальное скалярное условие имеет вид
$$Q(0)=1,$$
а все остальные отсутствующие значения автоматически задаются правилом нуля для отрицательных сдвигов.
Для \(1\le k\le n\) введем свернутый индекс и краевой множитель:
$$r(n,k)=\begin{cases} k,&k<n,\\ k-1,&k=n, \end{cases} \qquad \lambda(n,k)=\begin{cases} 1,&k<n,\\ 2,&k=n. \end{cases}$$
Тогда для любого активного слоя \(m\ge T(n)\) выполняются соотношения
$$\begin{aligned} A_{n,k}(m)=&\,A_{n,k}(m-n-2)+B_{n+1,1}(m-n-3)+A_{n+1,k+1}(m-n-3)\\ &+A_{n-1,r(n,k)}(m-n-1)+B_{n,1}(m-n-2)+A_{n,r(n,k)+1}(m-n-2), \end{aligned}$$
$$\begin{aligned} B_{n,k}(m)=&\,B_{n,k}(m-n-2)+B_{n+1,k+1}(m-n-3)+\lambda(n,k)\,A_{n+1,1}(m-n-3)\\ &+B_{n-1,r(n,k)}(m-n-1)+B_{n,r(n,k)+1}(m-n-2)+\lambda(n,k)\,A_{n,1}(m-n-2). \end{aligned}$$
Именно эти переходы и вычисляют реализации. Особый случай \(k=n\) является единственным местом, где возникает коэффициент \(2\); следовательно, край треугольной области индексов имеет немного иной комбинаторный вес, чем ее внутренность.
После того как для данного \(m\) заполнены все состояния \(A\) и \(B\), реализации обновляют скалярный слой по формулам
$$P(m)=Q(m-1)+4P(m-2)+2A_{1,1}(m-3)+B_{1,1}(m-3),$$
и, при \(m\ge 1\),
$$Q(m)=3Q(m-1)+3P(m-2).$$
Вторая рекурсия порождает последовательность, сдвинутые значения которой и являются окончательными ответами. Поэтому
$$D(n)=Q(n-1).$$
С точки зрения динамического программирования, \(P(m)\) играет роль моста между граничным уровнем \(n=0\) и настоящими двумерно индексированными состояниями профиля при \(n\ge 1\).
Каждая правая часть смотрит назад только на сдвиги вида \(m-n-1\), \(m-n-2\), \(m-n-3\), \(m-1\), \(m-2\) и \(m-3\). Полная история всех предыдущих слоев не нужна.
Если
$$n_{\max}(M)=\max\{n:T(n)\le M\},$$
то для целевого размера \(M\) все необходимые обращения назад помещаются в кольцевой буфер длины \(n_{\max}(M)+8\). Именно поэтому программа может выполнять квадратичную по времени рекурсию, не требуя при этом квадратичной памяти по параметру \(m\).
Начальные слои можно вычислить вручную прямо из рекурсии.
При \(m=0\) еще ни один треугольный слой не активен, поэтому
$$P(0)=0,\qquad Q(0)=1.$$
При \(m=1\) по-прежнему нет состояний с \(n\ge 1\), так как \(T(1)=3\). Следовательно,
$$P(1)=Q(0)=1,\qquad Q(1)=3Q(0)=3.$$
При \(m=2\) тот же аргумент дает
$$P(2)=Q(1)+4P(0)=3,\qquad Q(2)=3Q(1)+3P(0)=9.$$
При \(m=3\) появляется первый треугольный уровень, поскольку \(T(1)=3\). Для единственной пары \((n,k)=(1,1)\) имеем \(r(1,1)=0\) и \(\lambda(1,1)=2\). Используя граничное правило \(A_{0,0}(1)=B_{0,0}(1)=P(1)=1\), все недопустимые и отрицательно сдвинутые члены обнуляются, так что
$$A_{1,1}(3)=1,\qquad B_{1,1}(3)=1.$$
После этого скалярное обновление дает
$$P(3)=Q(2)+4P(1)=9+4=13,$$
$$Q(3)=3Q(2)+3P(1)=27+3=30.$$
Значит,
$$D(1)=1,\qquad D(2)=3,\qquad D(3)=9,\qquad D(4)=30.$$
Этот небольшой расчет наглядно показывает, как первый ненулевой треугольный профиль появляется ровно в тот момент, когда порог становится достижим.
Реализации на C++, Python и Java вычисляют одну и ту же рекурсию. Компилируемые решатели выделяют две треугольные таблицы состояний для допустимых пар профиля \((n,k)\) и два одномерных массива для скалярного слоя. Внешний цикл увеличивает параметр размера \(m\) от \(0\) до целевого значения, очищает текущий столбец кольцевого буфера, заполняет все активные состояния профиля и лишь затем обновляет две скалярные последовательности.
До начала динамики код по треугольному порогу вычисляет максимальный активный уровень \(n_{\max}\). Это значение определяет и число необходимых строк профиля, и длину кольцевого буфера. Поскольку рекурсия обращается на одну строку выше текущего уровня, при выделении памяти добавляется небольшой запас сверх \(n_{\max}\).
Реализация на C++ используется в двух режимах арифметики: точные целые числа для небольших контрольных значений и модульная арифметика для большого итогового запроса. Реализация на Java оставляет только модульный вариант, которого достаточно для опубликованного ответа. Реализация на Python представляет собой тонкий запускатель, который собирает и выполняет компилируемый решатель, а затем извлекает итоговый числовой результат, поэтому все три версии имеют одно и то же математическое ядро.
Встроенные контрольные значения таковы:
$$D(10)=44499,\qquad D(20)=9204559704,\qquad D(32)=22037102049132222,$$
а также
$$D(100)\equiv 780166455\pmod{10^9}.$$
Они подтверждают, что рекурсия собрана корректно, прежде чем вычисляется \(D(10000)\bmod 10^9\).
Пусть \(M\) обозначает наибольший параметр размера, до которого доходит программа; для финального запроса это \(M=n-1\). Активный диапазон задается неравенством \(T(n_{\max})\le M<T(n_{\max}+1)\), откуда \(n_{\max}=\Theta(\sqrt M)\).
Для каждого \(m\) реализация просматривает каждую допустимую пару \((n,k)\) с \(1\le n\le n_{\max}\) и \(1\le k\le n\). Число обновлений профиля на одном слое равно
$$\sum_{n=1}^{n_{\max}} n=\frac{n_{\max}(n_{\max}+1)}{2}=O(n_{\max}^2),$$
поэтому полное время работы составляет
$$O(M\,n_{\max}^2)=O(M^2).$$
По памяти длина кольцевого буфера равна \(O(n_{\max})\), а общее число сохраненных элементов профиля равно
$$O\!\left(n_{\max}\sum_{n=1}^{n_{\max}} n\right)=O(n_{\max}^3)=O(M^{3/2}).$$
Две скалярные последовательности добавляют лишь \(O(n_{\max})\) ячеек и не меняют старший асимптотический член.
تعرّف المسألة متتالية \(D(n)\) مرتبطة بعدّ تكوينات الأميبا داخل شبكة ثلاثية الأبعاد. ومن الواضح من التطبيقات أن الهدف الحسابي هو إيجاد \(D(10000)\) بترديد \(10^9\)، مع التحقق أولاً من عدة قيم أصغر بطريقة دقيقة قبل التشغيل الكبير بالمودولو.
التعداد المباشر لكل الأشكال ثلاثية الأبعاد ذات الصلة ينمو بسرعة هائلة، ولذلك لا يصلح عملياً. لهذا تستبدل الخوارزمية بناء الأشكال الصريح ببرمجة ديناميكية طبقية، حيث تمثل الحالات ملفات حدية مسموحاً بها، ثم تُغلق المنظومة بعلاقتين عدديتين لاستخراج القيمة النهائية للمتتالية.
لأغراض الشرح، نرمز إلى عائلتي الحالات المثلثيتين بـ \(A_{n,k}(m)\) و \(B_{n,k}(m)\)، وإلى المتتاليتين العدديتين بـ \(P(m)\) و \(Q(m)\). عندئذ تكون الكمية المطلوبة هي
$$D(n)=Q(n-1).$$
لا تقوم التطبيقات ببناء الأميبات واحدة واحدة. بل إنها تعدّ عدد التكوينات الجزئية التي يمكن أن تنتهي في كل ملف مسموح، ثم تستخدم علاقتين أحاديتَي البعد لإغلاق النظام بالكامل.
أول كمية بنيوية هي العتبة المثلثية
$$T(n)=\frac{(n+1)(n+2)}{2}.$$
وهذه هي أصغر قيمة لمعامل الحجم تسمح للطبقة ذات الفهرس \(n\) بأن تسهم. لذلك، لكل \(n\ge 1\)، لدينا
$$m<T(n)\implies A_{n,k}(m)=B_{n,k}(m)=0.$$
أي أنه عندما تعالج الحلقة الخارجية قيمة ثابتة من \(m\)، فإن الطبقات المهمة فقط هي تلك التي تحقق \(T(n)\le m\). وبما أن \(T(n)\sim n^2/2\)، فإن أكبر طبقة فعالة لها رتبة \(O(\sqrt m)\) فقط.
الأحجام السالبة غير ممكنة، ولذلك تُعرَّف كل كمية عند \(m<0\) بأنها \(0\). وعند المستوى \(n=0\)، لا تخزن عائلتا الحالات صفاً كاملاً من الملفات، بل تختزلان إلى التيار القياسي نفسه:
$$A_{0,0}(m)=B_{0,0}(m)=P(m),$$
بينما
$$A_{0,k}(m)=B_{0,k}(m)=0\qquad (k\ne 0).$$
هذا التطابق أساسي لأن علاقات المستويات الأعلى ترجع إلى \(n-1\)، وبالتالي يجب على المنظومة أن تعرف ماذا يحدث عندما يصل الفهرس إلى الصفر. والشرط الابتدائي القياسي هو
$$Q(0)=1,$$
أما القيم المفقودة الأخرى فتنتج تلقائياً من قاعدة الانعدام عند الإزاحات السالبة.
لكل \(1\le k\le n\)، نعرّف الفهرس المطوي ومعامل الحافة كما يلي:
$$r(n,k)=\begin{cases} k,&k<n,\\ k-1,&k=n, \end{cases} \qquad \lambda(n,k)=\begin{cases} 1,&k<n,\\ 2,&k=n. \end{cases}$$
عندئذ، لكل طبقة فعالة تحقق \(m\ge T(n)\)، تلبي عائلتا الملفات العلاقتين
$$\begin{aligned} A_{n,k}(m)=&\,A_{n,k}(m-n-2)+B_{n+1,1}(m-n-3)+A_{n+1,k+1}(m-n-3)\\ &+A_{n-1,r(n,k)}(m-n-1)+B_{n,1}(m-n-2)+A_{n,r(n,k)+1}(m-n-2), \end{aligned}$$
$$\begin{aligned} B_{n,k}(m)=&\,B_{n,k}(m-n-2)+B_{n+1,k+1}(m-n-3)+\lambda(n,k)\,A_{n+1,1}(m-n-3)\\ &+B_{n-1,r(n,k)}(m-n-1)+B_{n,r(n,k)+1}(m-n-2)+\lambda(n,k)\,A_{n,1}(m-n-2). \end{aligned}$$
وهاتان هما بالضبط العلاقتان اللتان تنفذهما التطبيقات. ولا يظهر المعامل \(2\) إلا في الحالة الحدية \(k=n\)، مما يعني أن طرف المجال المثلثي للفهرسة يحمل وزناً تركيبياً مختلفاً قليلاً عن داخله.
بعد إكمال جميع حالات \(A\) و \(B\) عند الحجم \(m\)، تحدّث التطبيقات الطبقة القياسية وفقاً للعلاقة
$$P(m)=Q(m-1)+4P(m-2)+2A_{1,1}(m-3)+B_{1,1}(m-3),$$
وكذلك، عندما \(m\ge 1\)،
$$Q(m)=3Q(m-1)+3P(m-2).$$
العلاقة الثانية هي التي تولد المتتالية التي تُستخرج منها الإجابة بعد إزاحة بمقدار واحد، ولذلك
$$D(n)=Q(n-1).$$
ومن منظور البرمجة الديناميكية، تلعب \(P(m)\) دور الجسر بين طبقة الحافة \(n=0\) وبين حالات الملفات الحقيقية ذات المؤشرين عندما \(n\ge 1\).
كل طرف أيمن ينظر إلى إزاحات سابقة محدودة فقط، مثل \(m-n-1\) و \(m-n-2\) و \(m-n-3\) و \(m-1\) و \(m-2\) و \(m-3\). ولذلك لا توجد حاجة إلى الاحتفاظ بتاريخ كامل لجميع الطبقات السابقة.
إذا كان
$$n_{\max}(M)=\max\{n:T(n)\le M\},$$
فإن كل الرجوعات المطلوبة للحجم النهائي \(M\) تدخل داخل مخزن دائري طوله \(n_{\max}(M)+8\). ولهذا يمكن للبرنامج أن يحافظ على زمن تربيعي من دون أن يحتاج أيضاً إلى ذاكرة تربيعية في البعد \(m\).
يمكن حساب الطبقات الأولى يدوياً مباشرة من العلاقات السابقة.
عند \(m=0\)، لا تكون أي طبقة مثلثية فعالة بعد، ومن ثم
$$P(0)=0,\qquad Q(0)=1.$$
وعند \(m=1\)، ما زالت جميع الحالات ذات \(n\ge 1\) منعدمة لأن \(T(1)=3\)، وبالتالي
$$P(1)=Q(0)=1,\qquad Q(1)=3Q(0)=3.$$
وعند \(m=2\)، يعطينا المنطق نفسه
$$P(2)=Q(1)+4P(0)=3,\qquad Q(2)=3Q(1)+3P(0)=9.$$
أما عند \(m=3\)، فتظهر أول طبقة مثلثية لأن \(T(1)=3\). وللزوج الوحيد \((n,k)=(1,1)\) يكون لدينا \(r(1,1)=0\) و \(\lambda(1,1)=2\). وباستخدام شرط الحافة \(A_{0,0}(1)=B_{0,0}(1)=P(1)=1\)، تختفي جميع الحدود غير الصالحة أو ذات الإزاحة السالبة، فنحصل على
$$A_{1,1}(3)=1,\qquad B_{1,1}(3)=1.$$
ثم يصبح التحديث القياسي
$$P(3)=Q(2)+4P(1)=9+4=13,$$
$$Q(3)=3Q(2)+3P(1)=27+3=30.$$
ومن ثم
$$D(1)=1,\qquad D(2)=3,\qquad D(3)=9,\qquad D(4)=30.$$
يبين هذا المثال الصغير بدقة كيف تدخل أول حالة مثلثية غير صفرية فور تحقق العتبة اللازمة.
تقوم تطبيقات C++ و Python و Java كلها بتقييم العلاقة نفسها. وتخصص الحلول المترجمة جدولين مثلثيين للحالات من أجل الأزواج المسموح بها \((n,k)\)، بالإضافة إلى مصفوفتين أحاديتَي البعد للطبقة القياسية. تزيد الحلقة الخارجية معامل الحجم \(m\) من \(0\) حتى الهدف، وتمسح العمود الحالي من المخزن الدائري، ثم تملأ كل حالات الملفات الفعالة، وبعد ذلك فقط تحدّث المتتاليتين القياسيتين.
قبل بدء البرمجة الديناميكية، يحسب الكود أكبر طبقة فعالة \(n_{\max}\) انطلاقاً من العتبة المثلثية. وهذا يحدد عدد صفوف الملفات المطلوبة وكذلك طول المخزن الدائري. وبما أن العلاقة تصل إلى صف أعلى من الصف الحالي، فإن الحجز الفعلي يتضمن هامش أمان صغيراً فوق \(n_{\max}\).
تستخدم نسخة C++ نمطين حسابيين: حساباً صحيحاً دقيقاً لنقاط تحقق صغيرة، وحساباً بترديد \(10^9\) للاستعلام النهائي الكبير. وتحتفظ نسخة Java بالصيغة المعيارية فقط، وهي كافية للإجابة المنشورة. أما نسخة Python فهي مشغل خفيف يبني الحل المترجم ويشغله ثم يستخرج الناتج العددي النهائي، ولذلك تشترك اللغات الثلاث في الجوهر الرياضي نفسه.
قيم التحقق المضمنة هي
$$D(10)=44499,\qquad D(20)=9204559704,\qquad D(32)=22037102049132222,$$
وكذلك
$$D(100)\equiv 780166455\pmod{10^9}.$$
وتؤكد هذه القيم أن العلاقات صيغت بشكل صحيح قبل حساب \(D(10000)\bmod 10^9\).
ليكن \(M\) أكبر معامل حجم تصل إليه الخوارزمية؛ وفي الاستعلام النهائي يكون \(M=n-1\). ويحقق المجال الفعال العلاقة \(T(n_{\max})\le M<T(n_{\max}+1)\)، ومن ثم \(n_{\max}=\Theta(\sqrt M)\).
لكل قيمة \(m\)، تزور الخوارزمية كل زوج مسموح \((n,k)\) بحيث \(1\le n\le n_{\max}\) و \(1\le k\le n\). وعدد تحديثات الملفات في كل طبقة يساوي
$$\sum_{n=1}^{n_{\max}} n=\frac{n_{\max}(n_{\max}+1)}{2}=O(n_{\max}^2),$$
ولذلك يكون الزمن الكلي
$$O(M\,n_{\max}^2)=O(M^2).$$
أما من جهة الذاكرة، فطول المخزن الدائري هو \(O(n_{\max})\)، وعدد خانات الملفات المخزنة يساوي
$$O\!\left(n_{\max}\sum_{n=1}^{n_{\max}} n\right)=O(n_{\max}^3)=O(M^{3/2}).$$
ولا تضيف المتتاليتان القياسيتان سوى \(O(n_{\max})\) خانة إضافية، ولذلك لا يتغير الحد الأسيمبتوطي الغالب.