The computation is naturally split into two stages. First we build a sequence \(L_n\) of lobster-tree counts through formal power-series algebra. Then we apply the Euler transform to \(L_n\) and obtain the target sequence \(T_n\). The requested answer is \(T_{2019}\bmod 10^9+7\), so the real task is to derive the generating functions in a form that can be truncated efficiently up to degree \(2019\).
Let \(M=10^9+7\). Write
$$L(x)=\sum_{n\ge 1} L_nx^n,\qquad T(x)=\sum_{n\ge 0} T_nx^n.$$
The C++, Python, and Java implementations all follow the same algebraic pipeline.
The base series is the ordinary generating function for integer partitions:
$$P(x)=\prod_{m\ge 1}\frac{1}{1-x^m}=\sum_{n\ge 0} p(n)x^n.$$
Its first coefficients are
$$P(x)=1+x+2x^2+3x^3+5x^4+7x^5+11x^6+15x^7+\cdots.$$
The implementations also use the shifted series
$$X(x)=xP(x)=x+x^2+2x^3+3x^4+5x^5+7x^6+\cdots.$$
This partition baseline is the raw supply of branch-multiset data from which the lobster-tree series is assembled.
Define two auxiliary series
$$A(x)=P(x)-\frac{1}{1-x},\qquad B(x)=P(x^2)-\frac{1}{1-x^2}.$$
From them the implementation constructs
$$U(x)=\frac{A(x)^2}{1-X(x)},$$
$$V(x)=\frac{B(x)\left(1+X(x)\right)}{1-x^2P(x^2)}.$$
These are the two large rational pieces that feed the lobster-tree count. The code averages the two cases, so the combined contribution later appears as \(\frac{1}{2}(U(x)+V(x))\).
The averaged rational contribution is shifted by two vertices, so it enters \(L(x)\) as
$$\frac{x^2}{2}\left(U(x)+V(x)\right).$$
Then the simpler family \(X(x)\) is added back in. Finally a correction series is subtracted. In the implementations the correction coefficients \(r_n\) satisfy
$$r_0=1,\qquad r_n=r_{n-1}+r_{n-2}-r_{n-3}\quad (n\ge 1),$$
with the missing negative indices treated as \(0\). Therefore the correction generating function is
$$R(x)=\sum_{n\ge 0} r_nx^n=\frac{1}{1-x-x^2+x^3}=\frac{1}{(1-x)(1-x^2)}.$$
So the lobster-tree generating function used by the implementations is
$$L(x)=X(x)+\frac{x^2}{2}\left(U(x)+V(x)\right)-x^3R(x).$$
Substituting the definitions of \(U\), \(V\), and \(R\) gives the closed form
$$L(x)=xP(x)+\frac{x^2}{2}\left(\frac{\left(P(x)-\frac{1}{1-x}\right)^2}{1-xP(x)}+\frac{\left(P(x^2)-\frac{1}{1-x^2}\right)\left(1+xP(x)\right)}{1-x^2P(x^2)}\right)-\frac{x^3}{(1-x)(1-x^2)}.$$
Expanding the rational pieces shows
$$U(x)=x^4+5x^5+18x^6+53x^7+\cdots,$$
$$V(x)=x^4+x^5+4x^6+5x^7+\cdots,$$
$$\frac{x^3}{(1-x)(1-x^2)}=x^3+x^4+2x^5+2x^6+3x^7+3x^8+\cdots.$$
Combining these with \(X(x)\) gives
$$L(x)=x+x^2+x^3+2x^4+3x^5+6x^6+11x^7+23x^8+47x^9+105x^{10}+\cdots.$$
Hence
$$L_1,\dots,L_{10}=1,1,1,2,3,6,11,23,47,105,$$
which matches the checked prefix in the implementations.
The target series is the Euler transform of the lobster counts:
$$T(x)=\prod_{d\ge 1}(1-x^d)^{-L_d}.$$
Define the divisor sum
$$c_m=\sum_{d\mid m} d\,L_d.$$
Taking a logarithmic derivative gives
$$\log T(x)=\sum_{d\ge 1} L_d\sum_{j\ge 1}\frac{x^{dj}}{j},$$
so
$$\frac{xT'(x)}{T(x)}=\sum_{m\ge 1} c_mx^m.$$
Comparing coefficients with \(T(x)=\sum_{n\ge 0}T_nx^n\) yields
$$T_0=1,\qquad nT_n=\sum_{k=1}^{n} c_kT_{n-k}.$$
Since the modulus is prime, division by \(n\) is performed with modular inverses:
$$T_n=\frac{1}{n}\sum_{k=1}^{n} c_kT_{n-k}\pmod{M}.$$
Using the first lobster counts \(L_1=1\), \(L_2=1\), \(L_3=1\), \(L_4=2\), we obtain
$$c_1=1,\qquad c_2=1+2=3,\qquad c_3=1+3=4,\qquad c_4=1+2+8=11.$$
Now compute the first transformed values:
$$T_1=\frac{c_1T_0}{1}=1,$$
$$T_2=\frac{c_1T_1+c_2T_0}{2}=\frac{1+3}{2}=2,$$
$$T_3=\frac{c_1T_2+c_2T_1+c_3T_0}{3}=\frac{2+3+4}{3}=3,$$
$$T_4=\frac{c_1T_3+c_2T_2+c_3T_1+c_4T_0}{4}=\frac{3+6+4+11}{4}=6.$$
Continuing in the same way gives \(T_7=37\), \(T_{10}=328\), and \(T_{20}=1416269\), exactly the values used as checkpoints by the implementation.
The C++, Python, and Java implementations precompute modular inverses \(1^{-1},2^{-1},\dots,n^{-1}\) because the formulas divide by \(2\) and by \(n\). They then build the partition series \(P(x)\) with the standard dynamic program for partition numbers. Next they create the two rational series above, invert the denominators coefficient-by-coefficient as formal power series, and multiply series with truncated convolutions so only terms up to degree \(n\) are kept.
After assembling \(L_1,\dots,L_n\), the implementation computes every divisor sum \(c_m\) by adding the contribution of each divisor to all of its multiples. Finally it evaluates the Euler-transform recurrence from \(T_0=1\) up to the requested index. The C++ version optionally parallelizes the largest convolution loops, but all three languages implement the same mathematics and return the same coefficients.
Let \(n\) be the requested index. Building the partition series costs \(O(n^2)\) time and \(O(n)\) memory. Each truncated convolution and each formal power-series inversion is also \(O(n^2)\), and only a constant number of such arrays are stored at once. The divisor-sum stage costs \(O(n\log n)\), while the final Euler-transform recurrence is \(O(n^2)\). Therefore the overall running time is \(O(n^2)\) and the memory usage is \(O(n)\).
Die Rechnung zerfällt natürlich in zwei Stufen. Zuerst wird eine Folge \(L_n\) von Lobster-Baum-Zahlen mit formalen Potenzreihen aufgebaut. Danach wird auf \(L_n\) die Euler-Transformation angewandt, wodurch die gesuchte Folge \(T_n\) entsteht. Gesucht ist also \(T_{2019}\bmod 10^9+7\); mathematisch bedeutet das, eine brauchbare Darstellung von \(L(x)\) herzuleiten und daraus \(T(x)\) effizient bis Grad \(2019\) zu gewinnen.
Sei \(M=10^9+7\). Wir schreiben
$$L(x)=\sum_{n\ge 1} L_nx^n,\qquad T(x)=\sum_{n\ge 0} T_nx^n.$$
Die C++-, Python- und Java-Implementierungen benutzen exakt dieselbe algebraische Kette.
Die Basis ist die gewöhnliche erzeugende Funktion der Partitionen:
$$P(x)=\prod_{m\ge 1}\frac{1}{1-x^m}=\sum_{n\ge 0} p(n)x^n.$$
Die ersten Koeffizienten sind
$$P(x)=1+x+2x^2+3x^3+5x^4+7x^5+11x^6+15x^7+\cdots.$$
Außerdem wird die verschobene Reihe
$$X(x)=xP(x)=x+x^2+2x^3+3x^4+5x^5+7x^6+\cdots$$
fortlaufend verwendet. Diese Partitionsbasis liefert die Rohdaten der Verzweigungs-Multimengen, aus denen die Lobster-Baum-Reihe zusammengesetzt wird.
Definiere die Hilfsreihen
$$A(x)=P(x)-\frac{1}{1-x},\qquad B(x)=P(x^2)-\frac{1}{1-x^2}.$$
Daraus konstruiert die Implementierung
$$U(x)=\frac{A(x)^2}{1-X(x)},$$
$$V(x)=\frac{B(x)\left(1+X(x)\right)}{1-x^2P(x^2)}.$$
Das sind die beiden großen rationalen Beiträge zur Lobster-Zählung. Anschließend werden beide Fälle gemittelt; deshalb erscheint später \(\frac{1}{2}(U(x)+V(x))\).
Der gemittelte rationale Beitrag wird um zwei Knoten verschoben und geht daher als
$$\frac{x^2}{2}\left(U(x)+V(x)\right)$$
in \(L(x)\) ein. Danach wird die einfachere Familie \(X(x)\) addiert. Schließlich wird eine Korrekturreihenfolge subtrahiert. In den Implementierungen erfüllen die Korrekturkoeffizienten \(r_n\)
$$r_0=1,\qquad r_n=r_{n-1}+r_{n-2}-r_{n-3}\quad (n\ge 1),$$
wobei fehlende negative Indizes als \(0\) behandelt werden. Damit ist die zugehörige erzeugende Funktion
$$R(x)=\sum_{n\ge 0} r_nx^n=\frac{1}{1-x-x^2+x^3}=\frac{1}{(1-x)(1-x^2)}.$$
Die von den Implementierungen benutzte Lobster-Reihe lautet also
$$L(x)=X(x)+\frac{x^2}{2}\left(U(x)+V(x)\right)-x^3R(x).$$
Nach Einsetzen von \(U\), \(V\) und \(R\) ergibt sich die geschlossene Form
$$L(x)=xP(x)+\frac{x^2}{2}\left(\frac{\left(P(x)-\frac{1}{1-x}\right)^2}{1-xP(x)}+\frac{\left(P(x^2)-\frac{1}{1-x^2}\right)\left(1+xP(x)\right)}{1-x^2P(x^2)}\right)-\frac{x^3}{(1-x)(1-x^2)}.$$
Die Entwicklung der rationalen Teile beginnt mit
$$U(x)=x^4+5x^5+18x^6+53x^7+\cdots,$$
$$V(x)=x^4+x^5+4x^6+5x^7+\cdots,$$
$$\frac{x^3}{(1-x)(1-x^2)}=x^3+x^4+2x^5+2x^6+3x^7+3x^8+\cdots.$$
Zusammen mit \(X(x)\) folgt
$$L(x)=x+x^2+x^3+2x^4+3x^5+6x^6+11x^7+23x^8+47x^9+105x^{10}+\cdots.$$
Also gilt
$$L_1,\dots,L_{10}=1,1,1,2,3,6,11,23,47,105,$$
genau wie im geprüften Präfix der Implementierungen.
Die Zielreihe ist die Euler-Transformation der Lobster-Zahlen:
$$T(x)=\prod_{d\ge 1}(1-x^d)^{-L_d}.$$
Definiere die Teilersumme
$$c_m=\sum_{d\mid m} d\,L_d.$$
Die logarithmische Ableitung liefert
$$\log T(x)=\sum_{d\ge 1} L_d\sum_{j\ge 1}\frac{x^{dj}}{j},$$
also
$$\frac{xT'(x)}{T(x)}=\sum_{m\ge 1} c_mx^m.$$
Koeffizientenvergleich mit \(T(x)=\sum_{n\ge 0}T_nx^n\) ergibt
$$T_0=1,\qquad nT_n=\sum_{k=1}^{n} c_kT_{n-k}.$$
Weil der Modulus prim ist, wird die Division durch \(n\) mittels modularer Inversen ausgeführt:
$$T_n=\frac{1}{n}\sum_{k=1}^{n} c_kT_{n-k}\pmod{M}.$$
Aus \(L_1=1\), \(L_2=1\), \(L_3=1\), \(L_4=2\) folgt
$$c_1=1,\qquad c_2=1+2=3,\qquad c_3=1+3=4,\qquad c_4=1+2+8=11.$$
Damit erhält man die ersten transformierten Werte:
$$T_1=\frac{c_1T_0}{1}=1,$$
$$T_2=\frac{c_1T_1+c_2T_0}{2}=\frac{1+3}{2}=2,$$
$$T_3=\frac{c_1T_2+c_2T_1+c_3T_0}{3}=\frac{2+3+4}{3}=3,$$
$$T_4=\frac{c_1T_3+c_2T_2+c_3T_1+c_4T_0}{4}=\frac{3+6+4+11}{4}=6.$$
Führt man das Verfahren fort, so erhält man \(T_7=37\), \(T_{10}=328\) und \(T_{20}=1416269\), exakt die Kontrollwerte der Implementierung.
Die C++-, Python- und Java-Implementierungen berechnen zuerst die modularen Inversen \(1^{-1},2^{-1},\dots,n^{-1}\), weil in den Formeln durch \(2\) und durch \(n\) geteilt wird. Danach wird die Partitionsreihe \(P(x)\) per Standard-Dynamik für Partitionszahlen aufgebaut. Anschließend werden die beiden rationalen Reihen erzeugt, ihre Nenner koeffizientenweise als formale Potenzreihen invertiert und die nötigen Produkte per gekürzter Faltung nur bis Grad \(n\) ausmultipliziert.
Wenn \(L_1,\dots,L_n\) feststehen, berechnet die Implementierung jede Teilersumme \(c_m\), indem der Beitrag jedes Teilers zu allen seinen Vielfachen addiert wird. Danach wird die Euler-Rekurrenz von \(T_0=1\) bis zum gewünschten Index ausgewertet. Die C++-Version kann die größten Faltungsschleifen parallelisieren; mathematisch und numerisch sind jedoch alle drei Sprachfassungen identisch.
Sei \(n\) der verlangte Index. Der Aufbau der Partitionsreihe kostet \(O(n^2)\) Zeit und \(O(n)\) Speicher. Jede gekürzte Faltung und jede Inversion einer formalen Potenzreihe ist ebenfalls \(O(n^2)\), und es wird nur eine konstante Anzahl solcher Arrays gehalten. Die Teilersummenphase benötigt \(O(n\log n)\), die abschließende Euler-Rekurrenz wiederum \(O(n^2)\). Insgesamt ergibt sich also \(O(n^2)\) Laufzeit bei \(O(n)\) Speicherverbrauch.
Hesap iki doğal aşamaya ayrılıyor. Önce formel kuvvet serileri kullanılarak lobster-ağaç sayıları \(L_n\) üretiliyor. Sonra bu diziye Euler dönüşümü uygulanarak hedef dizi \(T_n\) elde ediliyor. İstenen cevap \(T_{2019}\bmod 10^9+7\) olduğundan, asıl matematiksel iş \(L(x)\) için uygun bir kapalı ifade kurmak ve bunu derece \(2019\)'a kadar verimli biçimde kırpmaktır.
\(M=10^9+7\) olsun. Şöyle yazalım:
$$L(x)=\sum_{n\ge 1} L_nx^n,\qquad T(x)=\sum_{n\ge 0} T_nx^n.$$
C++, Python ve Java uygulamaları aynı cebirsel hattı izler.
Temel seri, tamsayı partisyonlarının olağan üreteç fonksiyonudur:
$$P(x)=\prod_{m\ge 1}\frac{1}{1-x^m}=\sum_{n\ge 0} p(n)x^n.$$
İlk katsayılar
$$P(x)=1+x+2x^2+3x^3+5x^4+7x^5+11x^6+15x^7+\cdots$$
şeklindedir. Ayrıca kaydırılmış seri
$$X(x)=xP(x)=x+x^2+2x^3+3x^4+5x^5+7x^6+\cdots$$
sürekli kullanılır. Bu partisyon tabanı, lobster-ağaç serisini kurarken gereken dal-çokluğu bilgisinin ham kaynağıdır.
Şu yardımcı serileri tanımlayalım:
$$A(x)=P(x)-\frac{1}{1-x},\qquad B(x)=P(x^2)-\frac{1}{1-x^2}.$$
Uygulama bunlardan
$$U(x)=\frac{A(x)^2}{1-X(x)},$$
$$V(x)=\frac{B(x)\left(1+X(x)\right)}{1-x^2P(x^2)}$$
serilerini kurar. Bunlar lobster sayımına giren iki büyük rasyonel bileşendir. Sonraki aşamada bu iki durumun ortalaması alındığı için katkı \(\frac{1}{2}(U(x)+V(x))\) olarak görünür.
Ortalaması alınan rasyonel katkı iki tepe kaydırılarak
$$\frac{x^2}{2}\left(U(x)+V(x)\right)$$
şeklinde \(L(x)\)'e girer. Ardından daha basit aile \(X(x)\) geri eklenir. Son olarak bir düzeltme serisi çıkarılır. Uygulamalarda bu düzeltmenin katsayıları \(r_n\)
$$r_0=1,\qquad r_n=r_{n-1}+r_{n-2}-r_{n-3}\quad (n\ge 1)$$
reküransını sağlar; negatif indeksler \(0\) kabul edilir. Buna karşılık gelen üreteç fonksiyonu
$$R(x)=\sum_{n\ge 0} r_nx^n=\frac{1}{1-x-x^2+x^3}=\frac{1}{(1-x)(1-x^2)}$$
olur. Dolayısıyla uygulamanın kullandığı lobster-ağaç serisi
$$L(x)=X(x)+\frac{x^2}{2}\left(U(x)+V(x)\right)-x^3R(x)$$
şeklindedir. \(U\), \(V\) ve \(R\) yerine yazınca kapalı form
$$L(x)=xP(x)+\frac{x^2}{2}\left(\frac{\left(P(x)-\frac{1}{1-x}\right)^2}{1-xP(x)}+\frac{\left(P(x^2)-\frac{1}{1-x^2}\right)\left(1+xP(x)\right)}{1-x^2P(x^2)}\right)-\frac{x^3}{(1-x)(1-x^2)}$$
olur.
Rasyonel parçaların ilk açılımları şöyledir:
$$U(x)=x^4+5x^5+18x^6+53x^7+\cdots,$$
$$V(x)=x^4+x^5+4x^6+5x^7+\cdots,$$
$$\frac{x^3}{(1-x)(1-x^2)}=x^3+x^4+2x^5+2x^6+3x^7+3x^8+\cdots.$$
Bunları \(X(x)\) ile birleştirince
$$L(x)=x+x^2+x^3+2x^4+3x^5+6x^6+11x^7+23x^8+47x^9+105x^{10}+\cdots$$
elde edilir. Yani
$$L_1,\dots,L_{10}=1,1,1,2,3,6,11,23,47,105$$
olur; bu da uygulamaların kontrol ettiği başlangıç dizisiyle aynıdır.
Hedef seri, lobster katsayılarının Euler dönüşümüdür:
$$T(x)=\prod_{d\ge 1}(1-x^d)^{-L_d}.$$
Şu bölen toplamını tanımlayalım:
$$c_m=\sum_{d\mid m} d\,L_d.$$
Logaritmik türev alınırsa
$$\log T(x)=\sum_{d\ge 1} L_d\sum_{j\ge 1}\frac{x^{dj}}{j},$$
dolayısıyla
$$\frac{xT'(x)}{T(x)}=\sum_{m\ge 1} c_mx^m$$
olur. \(T(x)=\sum_{n\ge 0}T_nx^n\) ile katsayı karşılaştırması yapınca
$$T_0=1,\qquad nT_n=\sum_{k=1}^{n} c_kT_{n-k}$$
elde edilir. Modül asal olduğu için bölme, modüler terslerle yapılır:
$$T_n=\frac{1}{n}\sum_{k=1}^{n} c_kT_{n-k}\pmod{M}.$$
\(L_1=1\), \(L_2=1\), \(L_3=1\), \(L_4=2\) için
$$c_1=1,\qquad c_2=1+2=3,\qquad c_3=1+3=4,\qquad c_4=1+2+8=11$$
olur. İlk dönüşmüş değerler de
$$T_1=\frac{c_1T_0}{1}=1,$$
$$T_2=\frac{c_1T_1+c_2T_0}{2}=\frac{1+3}{2}=2,$$
$$T_3=\frac{c_1T_2+c_2T_1+c_3T_0}{3}=\frac{2+3+4}{3}=3,$$
$$T_4=\frac{c_1T_3+c_2T_2+c_3T_1+c_4T_0}{4}=\frac{3+6+4+11}{4}=6$$
şeklinde bulunur. Aynı yöntem sürdürülünce \(T_7=37\), \(T_{10}=328\) ve \(T_{20}=1416269\) elde edilir; bunlar uygulamanın kullandığı kontrol değerleridir.
C++, Python ve Java uygulamaları önce \(1^{-1},2^{-1},\dots,n^{-1}\) modüler terslerini hesaplar; çünkü formüllerde hem \(2\)'ye hem de \(n\)'ye bölme vardır. Sonra partisyon sayıları için standart dinamik programla \(P(x)\) serisini oluştururlar. Ardından yukarıdaki iki rasyonel seri kurulur, paydalar formel kuvvet serisi olarak katsayı katsayı terslenir ve gerekli çarpımlar sadece derece \(n\)'e kadar tutulacak biçimde kırpılmış konvolüsyonlarla hesaplanır.
\(L_1,\dots,L_n\) hazır olduğunda, her \(c_m\) değeri bir bölenin katkısını tüm katlarına ekleyerek bulunur. Son adımda Euler dönüşümü reküransı \(T_0=1\)'den başlayıp istenen indekse kadar yürütülür. C++ sürümü en büyük konvolüsyon döngülerini isteğe bağlı paralelleştirebilir; fakat üç dilin yaptığı matematik ve ürettiği katsayılar aynıdır.
\(n\) istenen indeks olsun. Partisyon serisini kurmak \(O(n^2)\) zaman ve \(O(n)\) bellek gerektirir. Her kırpılmış konvolüsyon ve her formel seri tersleme işlemi de \(O(n^2)\) düzeyindedir; aynı anda yalnızca sabit sayıda dizi tutulur. Bölen toplamı aşaması \(O(n\log n)\), son Euler dönüşümü reküransı ise \(O(n^2)\) sürer. Bu yüzden toplam karmaşıklık \(O(n^2)\) zaman ve \(O(n)\) bellek olur.
El cálculo se divide de forma natural en dos etapas. Primero se construye una sucesión \(L_n\) de conteos de árboles lobster mediante álgebra de series formales. Después se aplica la transformada de Euler a \(L_n\) para obtener la sucesión objetivo \(T_n\). La respuesta pedida es \(T_{2019}\bmod 10^9+7\), así que el trabajo matemático consiste en derivar una expresión manejable para \(L(x)\) y truncarla eficientemente hasta grado \(2019\).
Sea \(M=10^9+7\). Escribimos
$$L(x)=\sum_{n\ge 1} L_nx^n,\qquad T(x)=\sum_{n\ge 0} T_nx^n.$$
Las implementaciones en C++, Python y Java siguen exactamente la misma cadena algebraica.
La serie base es la función generadora ordinaria de las particiones enteras:
$$P(x)=\prod_{m\ge 1}\frac{1}{1-x^m}=\sum_{n\ge 0} p(n)x^n.$$
Sus primeros coeficientes son
$$P(x)=1+x+2x^2+3x^3+5x^4+7x^5+11x^6+15x^7+\cdots.$$
También aparece constantemente la serie desplazada
$$X(x)=xP(x)=x+x^2+2x^3+3x^4+5x^5+7x^6+\cdots.$$
Esta base de particiones aporta la materia prima combinatoria con la que se arma la serie de árboles lobster.
Definimos las series auxiliares
$$A(x)=P(x)-\frac{1}{1-x},\qquad B(x)=P(x^2)-\frac{1}{1-x^2}.$$
A partir de ellas, la implementación forma
$$U(x)=\frac{A(x)^2}{1-X(x)},$$
$$V(x)=\frac{B(x)\left(1+X(x)\right)}{1-x^2P(x^2)}.$$
Estas son las dos grandes piezas racionales que alimentan el conteo lobster. Después se promedian, de modo que su aporte conjunto aparece como \(\frac{1}{2}(U(x)+V(x))\).
La contribución racional promediada se desplaza dos vértices, así que entra en \(L(x)\) como
$$\frac{x^2}{2}\left(U(x)+V(x)\right).$$
Luego se vuelve a sumar la familia más simple \(X(x)\). Finalmente se resta una serie de corrección. En las implementaciones, sus coeficientes \(r_n\) verifican
$$r_0=1,\qquad r_n=r_{n-1}+r_{n-2}-r_{n-3}\quad (n\ge 1),$$
tomando los índices negativos faltantes como \(0\). Por tanto, la función generadora de corrección es
$$R(x)=\sum_{n\ge 0} r_nx^n=\frac{1}{1-x-x^2+x^3}=\frac{1}{(1-x)(1-x^2)}.$$
Así, la serie de árboles lobster usada por las implementaciones es
$$L(x)=X(x)+\frac{x^2}{2}\left(U(x)+V(x)\right)-x^3R(x).$$
Sustituyendo \(U\), \(V\) y \(R\), obtenemos la forma cerrada
$$L(x)=xP(x)+\frac{x^2}{2}\left(\frac{\left(P(x)-\frac{1}{1-x}\right)^2}{1-xP(x)}+\frac{\left(P(x^2)-\frac{1}{1-x^2}\right)\left(1+xP(x)\right)}{1-x^2P(x^2)}\right)-\frac{x^3}{(1-x)(1-x^2)}.$$
Las expansiones iniciales de las piezas racionales son
$$U(x)=x^4+5x^5+18x^6+53x^7+\cdots,$$
$$V(x)=x^4+x^5+4x^6+5x^7+\cdots,$$
$$\frac{x^3}{(1-x)(1-x^2)}=x^3+x^4+2x^5+2x^6+3x^7+3x^8+\cdots.$$
Al combinarlas con \(X(x)\), resulta
$$L(x)=x+x^2+x^3+2x^4+3x^5+6x^6+11x^7+23x^8+47x^9+105x^{10}+\cdots.$$
Por tanto,
$$L_1,\dots,L_{10}=1,1,1,2,3,6,11,23,47,105,$$
que coincide con el prefijo verificado por las implementaciones.
La serie objetivo es la transformada de Euler de los conteos lobster:
$$T(x)=\prod_{d\ge 1}(1-x^d)^{-L_d}.$$
Definimos la suma por divisores
$$c_m=\sum_{d\mid m} d\,L_d.$$
Tomando derivada logarítmica se obtiene
$$\log T(x)=\sum_{d\ge 1} L_d\sum_{j\ge 1}\frac{x^{dj}}{j},$$
y por tanto
$$\frac{xT'(x)}{T(x)}=\sum_{m\ge 1} c_mx^m.$$
Comparando coeficientes con \(T(x)=\sum_{n\ge 0}T_nx^n\), llegamos a
$$T_0=1,\qquad nT_n=\sum_{k=1}^{n} c_kT_{n-k}.$$
Como el módulo es primo, la división por \(n\) se hace con inversos modulares:
$$T_n=\frac{1}{n}\sum_{k=1}^{n} c_kT_{n-k}\pmod{M}.$$
Con \(L_1=1\), \(L_2=1\), \(L_3=1\), \(L_4=2\), se tiene
$$c_1=1,\qquad c_2=1+2=3,\qquad c_3=1+3=4,\qquad c_4=1+2+8=11.$$
Los primeros valores transformados son entonces
$$T_1=\frac{c_1T_0}{1}=1,$$
$$T_2=\frac{c_1T_1+c_2T_0}{2}=\frac{1+3}{2}=2,$$
$$T_3=\frac{c_1T_2+c_2T_1+c_3T_0}{3}=\frac{2+3+4}{3}=3,$$
$$T_4=\frac{c_1T_3+c_2T_2+c_3T_1+c_4T_0}{4}=\frac{3+6+4+11}{4}=6.$$
Si se continúa igual, se obtienen \(T_7=37\), \(T_{10}=328\) y \(T_{20}=1416269\), exactamente los valores de control usados por la implementación.
Las implementaciones en C++, Python y Java precalculan primero los inversos modulares \(1^{-1},2^{-1},\dots,n^{-1}\), porque las fórmulas dividen por \(2\) y por \(n\). Luego construyen la serie de particiones \(P(x)\) mediante el programa dinámico estándar para números de partición. Después forman las dos series racionales anteriores, invierten sus denominadores coeficiente a coeficiente como series formales y multiplican las series con convoluciones truncadas, conservando solo los términos hasta grado \(n\).
Una vez ensamblados \(L_1,\dots,L_n\), la implementación calcula cada suma \(c_m\) añadiendo la contribución de cada divisor a todos sus múltiplos. Por último evalúa la recurrencia de la transformada de Euler desde \(T_0=1\) hasta el índice pedido. La versión en C++ puede paralelizar las convoluciones más grandes, pero la matemática y los coeficientes devueltos son los mismos en los tres lenguajes.
Sea \(n\) el índice solicitado. Construir la serie de particiones cuesta \(O(n^2)\) tiempo y \(O(n)\) memoria. Cada convolución truncada y cada inversión de una serie formal también cuesta \(O(n^2)\), y solo se mantiene un número constante de arreglos de este tipo. La fase de sumas por divisores cuesta \(O(n\log n)\), mientras que la recurrencia final de la transformada de Euler vuelve a ser \(O(n^2)\). En consecuencia, el coste total es \(O(n^2)\) tiempo y \(O(n)\) memoria.
这个问题的计算流程天然分成两段。第一段先用形式幂级数构造一列 lobster 树计数 \(L_n\)。第二段再对 \(L_n\) 做 Euler 变换,得到目标序列 \(T_n\)。最终要求的是 \(T_{2019}\bmod 10^9+7\),所以真正要解决的是:怎样把 \(L(x)\) 写成适合截断到 \(2019\) 次的生成函数表达式,并据此恢复 \(T(x)\)。
记 \(M=10^9+7\)。设
$$L(x)=\sum_{n\ge 1}L_nx^n,\qquad T(x)=\sum_{n\ge 0}T_nx^n.$$
C++、Python 和 Java 三个实现采用的是同一条代数路线。
基础级数是整数拆分的普通生成函数:
$$P(x)=\prod_{m\ge 1}\frac{1}{1-x^m}=\sum_{n\ge 0}p(n)x^n.$$
它的前几项是
$$P(x)=1+x+2x^2+3x^3+5x^4+7x^5+11x^6+15x^7+\cdots.$$
实现里还会反复使用它的移位版本
$$X(x)=xP(x)=x+x^2+2x^3+3x^4+5x^5+7x^6+\cdots.$$
可以把这一步理解为:先把所有“分支多重集”编码进一个最基本的级数里,后面的 lobster 计数都是在这个基础上组合出来的。
先定义两个辅助级数
$$A(x)=P(x)-\frac{1}{1-x},\qquad B(x)=P(x^2)-\frac{1}{1-x^2}.$$
实现接着构造
$$U(x)=\frac{A(x)^2}{1-X(x)},$$
$$V(x)=\frac{B(x)\left(1+X(x)\right)}{1-x^2P(x^2)}.$$
这两个有理级数就是 lobster 计数中的两大主体贡献。随后程序把两种情况取平均,因此后面会出现 \(\frac{1}{2}(U(x)+V(x))\)。
取平均后的有理部分会整体向后平移两个顶点,所以它对 \(L(x)\) 的贡献是
$$\frac{x^2}{2}\left(U(x)+V(x)\right).$$
然后再把较简单的一族 \(X(x)\) 加回来。最后还要减去一个修正级数。实现中这个修正序列 \(r_n\) 满足
$$r_0=1,\qquad r_n=r_{n-1}+r_{n-2}-r_{n-3}\quad (n\ge 1),$$
缺失的负下标视为 \(0\)。因此它的生成函数是
$$R(x)=\sum_{n\ge 0}r_nx^n=\frac{1}{1-x-x^2+x^3}=\frac{1}{(1-x)(1-x^2)}.$$
于是,实现真正使用的 lobster 树生成函数就是
$$L(x)=X(x)+\frac{x^2}{2}\left(U(x)+V(x)\right)-x^3R(x).$$
把 \(U\)、\(V\)、\(R\) 全部代回去,可以写成单个闭式:
$$L(x)=xP(x)+\frac{x^2}{2}\left(\frac{\left(P(x)-\frac{1}{1-x}\right)^2}{1-xP(x)}+\frac{\left(P(x^2)-\frac{1}{1-x^2}\right)\left(1+xP(x)\right)}{1-x^2P(x^2)}\right)-\frac{x^3}{(1-x)(1-x^2)}.$$
把上面的有理部分展开,前几项分别是
$$U(x)=x^4+5x^5+18x^6+53x^7+\cdots,$$
$$V(x)=x^4+x^5+4x^6+5x^7+\cdots,$$
$$\frac{x^3}{(1-x)(1-x^2)}=x^3+x^4+2x^5+2x^6+3x^7+3x^8+\cdots.$$
再与 \(X(x)\) 合并,就得到
$$L(x)=x+x^2+x^3+2x^4+3x^5+6x^6+11x^7+23x^8+47x^9+105x^{10}+\cdots.$$
因此
$$L_1,\dots,L_{10}=1,1,1,2,3,6,11,23,47,105,$$
这与实现里校验的前缀完全一致。
目标级数是 lobster 计数的 Euler 变换:
$$T(x)=\prod_{d\ge 1}(1-x^d)^{-L_d}.$$
定义按约数求和的量
$$c_m=\sum_{d\mid m} d\,L_d.$$
对 \(T(x)\) 取对数导数,有
$$\log T(x)=\sum_{d\ge 1}L_d\sum_{j\ge 1}\frac{x^{dj}}{j},$$
所以
$$\frac{xT'(x)}{T(x)}=\sum_{m\ge 1}c_mx^m.$$
再与 \(T(x)=\sum_{n\ge 0}T_nx^n\) 对比系数,就得到
$$T_0=1,\qquad nT_n=\sum_{k=1}^{n}c_kT_{n-k}.$$
由于模数 \(M\) 是素数,所以除以 \(n\) 可以通过模逆来完成:
$$T_n=\frac{1}{n}\sum_{k=1}^{n}c_kT_{n-k}\pmod{M}.$$
已知 \(L_1=1\)、\(L_2=1\)、\(L_3=1\)、\(L_4=2\),于是
$$c_1=1,\qquad c_2=1+2=3,\qquad c_3=1+3=4,\qquad c_4=1+2+8=11.$$
然后递推得到
$$T_1=\frac{c_1T_0}{1}=1,$$
$$T_2=\frac{c_1T_1+c_2T_0}{2}=\frac{1+3}{2}=2,$$
$$T_3=\frac{c_1T_2+c_2T_1+c_3T_0}{3}=\frac{2+3+4}{3}=3,$$
$$T_4=\frac{c_1T_3+c_2T_2+c_3T_1+c_4T_0}{4}=\frac{3+6+4+11}{4}=6.$$
继续往后推,就会得到 \(T_7=37\)、\(T_{10}=328\)、\(T_{20}=1416269\),这正是实现中使用的检查值。
C++、Python 和 Java 实现首先预处理 \(1^{-1},2^{-1},\dots,n^{-1}\) 的模逆,因为公式里需要除以 \(2\) 和 \(n\)。接着用标准的整数拆分动态规划构造 \(P(x)\)。之后程序建立上面的两个有理级数,将分母按形式幂级数逐项求逆,并通过截断卷积只保留到 \(n\) 次项为止。
当 \(L_1,\dots,L_n\) 组装完成后,实现会把每个约数的贡献加到它所有倍数上,从而得到全部 \(c_m\)。最后从 \(T_0=1\) 开始计算 Euler 变换递推式,一直到目标下标。C++ 版本可以选择并行化最大的卷积循环,但三种语言执行的数学步骤和输出系数完全相同。
设目标下标为 \(n\)。构造拆分级数需要 \(O(n^2)\) 时间和 \(O(n)\) 空间。每一次截断卷积以及每一次形式幂级数求逆也都是 \(O(n^2)\),而同时维护的数组个数只是常数个。约数求和阶段是 \(O(n\log n)\),最后的 Euler 变换递推仍然是 \(O(n^2)\)。因此总时间复杂度为 \(O(n^2)\),总空间复杂度为 \(O(n)\)。
Вычисление естественно распадается на два шага. Сначала строится последовательность \(L_n\), считающая lobster-деревья, с помощью алгебры формальных степенных рядов. Затем к \(L_n\) применяется преобразование Эйлера, и получается целевая последовательность \(T_n\). Требуется число \(T_{2019}\bmod 10^9+7\), поэтому главная математическая задача состоит в том, чтобы вывести удобную формулу для \(L(x)\) и эффективно обрезать ее до степени \(2019\).
Пусть \(M=10^9+7\). Обозначим
$$L(x)=\sum_{n\ge 1}L_nx^n,\qquad T(x)=\sum_{n\ge 0}T_nx^n.$$
Реализации на C++, Python и Java используют один и тот же алгебраический конвейер.
Исходный ряд - обычная производящая функция для целочисленных разбиений:
$$P(x)=\prod_{m\ge 1}\frac{1}{1-x^m}=\sum_{n\ge 0}p(n)x^n.$$
Первые коэффициенты равны
$$P(x)=1+x+2x^2+3x^3+5x^4+7x^5+11x^6+15x^7+\cdots.$$
Также постоянно используется сдвинутый ряд
$$X(x)=xP(x)=x+x^2+2x^3+3x^4+5x^5+7x^6+\cdots.$$
Этот базовый ряд играет роль сырья: из него собираются все ветвящиеся конфигурации, из которых затем строится серия lobster-деревьев.
Введем вспомогательные ряды
$$A(x)=P(x)-\frac{1}{1-x},\qquad B(x)=P(x^2)-\frac{1}{1-x^2}.$$
Из них реализация строит
$$U(x)=\frac{A(x)^2}{1-X(x)},$$
$$V(x)=\frac{B(x)\left(1+X(x)\right)}{1-x^2P(x^2)}.$$
Это две большие рациональные составляющие счета lobster-объектов. Затем программа усредняет оба случая, поэтому дальше появляется вклад \(\frac{1}{2}(U(x)+V(x))\).
Усредненный рациональный вклад сдвигается на две вершины и входит в \(L(x)\) как
$$\frac{x^2}{2}\left(U(x)+V(x)\right).$$
После этого добавляется более простая семья \(X(x)\). Наконец, вычитается корректирующий ряд. В реализациях его коэффициенты \(r_n\) удовлетворяют рекуррентному соотношению
$$r_0=1,\qquad r_n=r_{n-1}+r_{n-2}-r_{n-3}\quad (n\ge 1),$$
где отсутствующие отрицательные индексы считаются нулевыми. Значит, его производящая функция равна
$$R(x)=\sum_{n\ge 0}r_nx^n=\frac{1}{1-x-x^2+x^3}=\frac{1}{(1-x)(1-x^2)}.$$
Итак, используемый в решении ряд lobster-деревьев имеет вид
$$L(x)=X(x)+\frac{x^2}{2}\left(U(x)+V(x)\right)-x^3R(x).$$
Подставляя определения \(U\), \(V\) и \(R\), получаем замкнутую формулу
$$L(x)=xP(x)+\frac{x^2}{2}\left(\frac{\left(P(x)-\frac{1}{1-x}\right)^2}{1-xP(x)}+\frac{\left(P(x^2)-\frac{1}{1-x^2}\right)\left(1+xP(x)\right)}{1-x^2P(x^2)}\right)-\frac{x^3}{(1-x)(1-x^2)}.$$
Разложение рациональных частей начинается так:
$$U(x)=x^4+5x^5+18x^6+53x^7+\cdots,$$
$$V(x)=x^4+x^5+4x^6+5x^7+\cdots,$$
$$\frac{x^3}{(1-x)(1-x^2)}=x^3+x^4+2x^5+2x^6+3x^7+3x^8+\cdots.$$
Вместе с \(X(x)\) это дает
$$L(x)=x+x^2+x^3+2x^4+3x^5+6x^6+11x^7+23x^8+47x^9+105x^{10}+\cdots.$$
Следовательно,
$$L_1,\dots,L_{10}=1,1,1,2,3,6,11,23,47,105,$$
и это в точности совпадает с проверяемым префиксом в реализациях.
Целевой ряд есть преобразование Эйлера от последовательности \(L_n\):
$$T(x)=\prod_{d\ge 1}(1-x^d)^{-L_d}.$$
Введем сумму по делителям
$$c_m=\sum_{d\mid m} d\,L_d.$$
Логарифмическая производная дает
$$\log T(x)=\sum_{d\ge 1}L_d\sum_{j\ge 1}\frac{x^{dj}}{j},$$
откуда
$$\frac{xT'(x)}{T(x)}=\sum_{m\ge 1}c_mx^m.$$
Сравнение коэффициентов с \(T(x)=\sum_{n\ge 0}T_nx^n\) приводит к формуле
$$T_0=1,\qquad nT_n=\sum_{k=1}^{n} c_kT_{n-k}.$$
Так как модуль \(M\) прост, деление на \(n\) выполняется через модульные обратные:
$$T_n=\frac{1}{n}\sum_{k=1}^{n} c_kT_{n-k}\pmod{M}.$$
Если взять \(L_1=1\), \(L_2=1\), \(L_3=1\), \(L_4=2\), то получаем
$$c_1=1,\qquad c_2=1+2=3,\qquad c_3=1+3=4,\qquad c_4=1+2+8=11.$$
Тогда первые значения преобразованной последовательности равны
$$T_1=\frac{c_1T_0}{1}=1,$$
$$T_2=\frac{c_1T_1+c_2T_0}{2}=\frac{1+3}{2}=2,$$
$$T_3=\frac{c_1T_2+c_2T_1+c_3T_0}{3}=\frac{2+3+4}{3}=3,$$
$$T_4=\frac{c_1T_3+c_2T_2+c_3T_1+c_4T_0}{4}=\frac{3+6+4+11}{4}=6.$$
Продолжая тот же процесс, получаем \(T_7=37\), \(T_{10}=328\) и \(T_{20}=1416269\), то есть именно те контрольные значения, которые использует реализация.
Реализации на C++, Python и Java сначала предвычисляют модульные обратные \(1^{-1},2^{-1},\dots,n^{-1}\), потому что в формулах есть деление на \(2\) и на \(n\). Затем стандартной динамикой для разбиений строится ряд \(P(x)\). После этого формируются два указанных рациональных ряда, их знаменатели инвертируются покоэффициентно как формальные степенные ряды, а произведения считаются усеченными свертками только до степени \(n\).
Когда коэффициенты \(L_1,\dots,L_n\) уже собраны, программа вычисляет каждую сумму \(c_m\), добавляя вклад каждого делителя ко всем его кратным. На последнем этапе рекуррентная формула преобразования Эйлера вычисляется от \(T_0=1\) до нужного индекса. Версия на C++ при желании распараллеливает самые крупные свертки, но математически все три реализации полностью совпадают.
Пусть \(n\) - нужный индекс. Построение ряда разбиений требует \(O(n^2)\) времени и \(O(n)\) памяти. Каждая усеченная свертка и каждая инверсия формального ряда также стоят \(O(n^2)\), причем одновременно хранится лишь константное число массивов. Этап сумм по делителям работает за \(O(n\log n)\), а финальная рекурсия преобразования Эйлера снова занимает \(O(n^2)\). Следовательно, полная сложность равна \(O(n^2)\) по времени и \(O(n)\) по памяти.
ينقسم الحساب طبيعيًا إلى مرحلتين. في المرحلة الأولى نبني المتتالية \(L_n\) التي تعد أشجار lobster باستخدام جبر السلاسل الشكلية. وفي المرحلة الثانية نطبق تحويل أويلر على \(L_n\) فنحصل على المتتالية الهدف \(T_n\). المطلوب في النهاية هو \(T_{2019}\bmod 10^9+7\)، ولذلك فجوهر الحل الرياضي هو إيجاد صيغة مناسبة للدالة المولدة \(L(x)\) ثم قصها بكفاءة حتى الدرجة \(2019\).
لنكتب \(M=10^9+7\)، ثم نعرف
$$L(x)=\sum_{n\ge 1}L_nx^n,\qquad T(x)=\sum_{n\ge 0}T_nx^n.$$
تتبع تطبيقات C++ وPython وJava نفس السلسلة الجبرية تمامًا.
الأساس هو دالة التوليد العادية لتقسيمات الأعداد الصحيحة:
$$P(x)=\prod_{m\ge 1}\frac{1}{1-x^m}=\sum_{n\ge 0}p(n)x^n.$$
وأول معاملات هذه السلسلة هي
$$P(x)=1+x+2x^2+3x^3+5x^4+7x^5+11x^6+15x^7+\cdots.$$
كما تستعمل الشيفرة كثيرًا السلسلة المزاحة
$$X(x)=xP(x)=x+x^2+2x^3+3x^4+5x^5+7x^6+\cdots.$$
يمكن النظر إلى هذه الخطوة على أنها توفير المخزون الخام لجميع تراكيب الفروع التي ستدخل لاحقًا في بناء سلسلة أشجار lobster.
نعرف السلسلتين المساعدتين
$$A(x)=P(x)-\frac{1}{1-x},\qquad B(x)=P(x^2)-\frac{1}{1-x^2}.$$
ثم تبني التطبيقات
$$U(x)=\frac{A(x)^2}{1-X(x)},$$
$$V(x)=\frac{B(x)\left(1+X(x)\right)}{1-x^2P(x^2)}.$$
هذان هما المكونان النسبيان الكبيران في عد أشجار lobster. وبعد ذلك يؤخذ متوسط الحالتين، ولهذا يظهر لاحقًا الحد \(\frac{1}{2}(U(x)+V(x))\).
بعد أخذ المتوسط، يُزاح هذا المكون النسبي بمقدار رأسين، لذا يساهم في \(L(x)\) على الصورة
$$\frac{x^2}{2}\left(U(x)+V(x)\right).$$
بعدها تضاف العائلة الأبسط \(X(x)\) من جديد، ثم يُطرح حد تصحيحي. في التطبيقات تحقق معاملات هذا الحد \(r_n\) العلاقة
$$r_0=1,\qquad r_n=r_{n-1}+r_{n-2}-r_{n-3}\quad (n\ge 1),$$
مع اعتبار الحدود ذات الفهارس السالبة المفقودة مساوية للصفر. لذلك تكون دالته المولدة
$$R(x)=\sum_{n\ge 0}r_nx^n=\frac{1}{1-x-x^2+x^3}=\frac{1}{(1-x)(1-x^2)}.$$
إذًا الدالة المولدة لأشجار lobster كما تنفذها التطبيقات هي
$$L(x)=X(x)+\frac{x^2}{2}\left(U(x)+V(x)\right)-x^3R(x).$$
وبالتعويض عن \(U\) و\(V\) و\(R\) نحصل على الصيغة المغلقة
$$L(x)=xP(x)+\frac{x^2}{2}\left(\frac{\left(P(x)-\frac{1}{1-x}\right)^2}{1-xP(x)}+\frac{\left(P(x^2)-\frac{1}{1-x^2}\right)\left(1+xP(x)\right)}{1-x^2P(x^2)}\right)-\frac{x^3}{(1-x)(1-x^2)}.$$
عند توسيع الجزأين النسبيين نحصل أولًا على
$$U(x)=x^4+5x^5+18x^6+53x^7+\cdots,$$
$$V(x)=x^4+x^5+4x^6+5x^7+\cdots,$$
$$\frac{x^3}{(1-x)(1-x^2)}=x^3+x^4+2x^5+2x^6+3x^7+3x^8+\cdots.$$
وبدمج هذه الحدود مع \(X(x)\) نحصل على
$$L(x)=x+x^2+x^3+2x^4+3x^5+6x^6+11x^7+23x^8+47x^9+105x^{10}+\cdots.$$
ومن ثم
$$L_1,\dots,L_{10}=1,1,1,2,3,6,11,23,47,105,$$
وهو تمامًا نفس المقطع الابتدائي الذي تتحقق منه التطبيقات.
السلسلة الهدف هي تحويل أويلر لمعاملات أشجار lobster:
$$T(x)=\prod_{d\ge 1}(1-x^d)^{-L_d}.$$
ونعرف مجموع القواسم
$$c_m=\sum_{d\mid m} d\,L_d.$$
بأخذ الاشتقاق اللوغاريتمي نجد أن
$$\log T(x)=\sum_{d\ge 1}L_d\sum_{j\ge 1}\frac{x^{dj}}{j},$$
ومن ثم
$$\frac{xT'(x)}{T(x)}=\sum_{m\ge 1}c_mx^m.$$
وبمقارنة المعاملات مع \(T(x)=\sum_{n\ge 0}T_nx^n\) نحصل على
$$T_0=1,\qquad nT_n=\sum_{k=1}^{n}c_kT_{n-k}.$$
وبما أن \(M\) عدد أولي، فإن القسمة على \(n\) تنفذ باستخدام المعكوسات الضربية:
$$T_n=\frac{1}{n}\sum_{k=1}^{n}c_kT_{n-k}\pmod{M}.$$
إذا أخذنا \(L_1=1\) و\(L_2=1\) و\(L_3=1\) و\(L_4=2\)، فإن
$$c_1=1,\qquad c_2=1+2=3,\qquad c_3=1+3=4,\qquad c_4=1+2+8=11.$$
وعندها تصبح القيم الأولى بعد التحويل
$$T_1=\frac{c_1T_0}{1}=1,$$
$$T_2=\frac{c_1T_1+c_2T_0}{2}=\frac{1+3}{2}=2,$$
$$T_3=\frac{c_1T_2+c_2T_1+c_3T_0}{3}=\frac{2+3+4}{3}=3,$$
$$T_4=\frac{c_1T_3+c_2T_2+c_3T_1+c_4T_0}{4}=\frac{3+6+4+11}{4}=6.$$
وبمتابعة العملية نفسها نحصل على \(T_7=37\) و\(T_{10}=328\) و\(T_{20}=1416269\)، وهي بالضبط قيم التحقق الموجودة في التنفيذ.
تبدأ تطبيقات C++ وPython وJava بحساب المعكوسات الضربية \(1^{-1},2^{-1},\dots,n^{-1}\) لأن الصيغ تحتوي على قسمة على \(2\) وعلى \(n\). بعد ذلك تُبنى السلسلة \(P(x)\) بواسطة البرمجة الديناميكية القياسية لأعداد التقسيمات. ثم تُنشأ السلسلتان النسبيتان السابقتان، وتُقلب المقامات حدًا بعد حد كسلاسل شكلية، وتُحسب الضروبات بواسطة التفافات مقصوصة تحتفظ فقط بالحدود حتى الدرجة \(n\).
وبعد تجميع \(L_1,\dots,L_n\)، تحسب التطبيقات كل قيمة \(c_m\) بإضافة مساهمة كل قاسم إلى جميع مضاعفاته. وأخيرًا تُقيَّم علاقة تحويل أويلر بدءًا من \(T_0=1\) وحتى الفهرس المطلوب. ويمكن لنسخة C++ أن توازي أكبر حلقات الالتفاف، لكن الخطوات الرياضية والمعاملات النهائية متطابقة في اللغات الثلاث.
إذا كان \(n\) هو الفهرس المطلوب، فإن بناء سلسلة التقسيمات يكلف \(O(n^2)\) زمنًا و\(O(n)\) ذاكرة. وكل التفاف مقصوص وكل قلب لسلسلة شكلية يكلف أيضًا \(O(n^2)\)، مع الاحتفاظ بعدد ثابت فقط من المصفوفات في الذاكرة. أما مرحلة جمع القواسم فتكلف \(O(n\log n)\)، في حين أن علاقة تحويل أويلر النهائية تعود إلى \(O(n^2)\). لذلك فالتعقيد الكلي هو \(O(n^2)\) زمنًا و\(O(n)\) ذاكرة.