The quantity to compute is \(B(k,n)\) modulo \(M=10^9+7\). The C++, Python, and Java implementations do not enumerate matrix states directly. Instead, they rewrite the target value as a coefficient extraction problem built from factorial-based series. The difficult part is a family of sums involving sixteenth powers of binomial coefficients, and the implementations evaluate all of those sums at once by fast polynomial convolution.
For the target instance, \(k\) is as large as \(100000\), so a direct quadratic summation would be far too slow. The entire strategy is therefore organized around turning the problem into two convolutions of length about \(k\), followed by a final modular reconstruction.
Write
$$q=\left\lfloor\frac{n}{k}\right\rfloor,\qquad M=10^9+7,\qquad \mu\equiv 2^q-2 \pmod{M}.$$
The implemented formula can be derived cleanly in the following steps.
All dependence on \(n\) enters through the single scalar \(\mu\). Because \(M\) is prime and \(2\not\equiv 0\pmod{M}\), Fermat's little theorem gives
$$2^q \equiv 2^{\,q\bmod (M-1)} \pmod{M}.$$
So even when \(q\) is enormous, the implementation only needs the reduced exponent \(q\bmod(M-1)\). After this reduction, the rest of the computation depends only on \(k\) and on the already-computed value of \(\mu\).
Define the sequence
$$c_i=(i!)^{-16}\pmod{M}\qquad (0\le i\le k).$$
Now form its self-convolution:
$$d_s=\sum_{i=0}^{s} c_i\,c_{s-i} \qquad (0\le s\le k).$$
Multiplying by \((s!)^{16}\) transforms this into
$$a_s=(s!)^{16}d_s=(s!)^{16}\sum_{i=0}^{s}\frac{1}{(i!)^{16}((s-i)!)^{16}}=\sum_{i=0}^{s}\left(\frac{s!}{i!(s-i)!}\right)^{16}.$$
Therefore
$$\boxed{a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}}.$$
This identity is the key simplification: instead of evaluating each \(a_s\) separately, one self-convolution computes every value \(a_0,a_1,\dots,a_k\) in one shot.
After the quantities \(a_s\) are known, define
$$g_s=\frac{a_s}{s!},\qquad e_j=\frac{\mu^j}{j!}.$$
Convolving these two sequences gives
$$f_t=\sum_{s=0}^{t} g_s\,e_{t-s}=\sum_{s=0}^{t}\frac{a_s}{s!}\frac{\mu^{\,t-s}}{(t-s)!}.$$
The final answer is the \(t=k\) case multiplied by \(k!\):
$$B(k,n)=k!\,f_k.$$
Expanding the factorials yields the explicit summation formula
$$\boxed{B(k,n)=\sum_{s=0}^{k}\binom{k}{s}a_s\,\mu^{\,k-s}\pmod{M}.}$$
Equivalently, if we package the coefficients into an exponential generating function, then
$$\frac{B(k,n)}{k!}=\left[x^k\right]\left(\sum_{s\ge 0}a_s\frac{x^s}{s!}\right)e^{\mu x}.$$
So the whole problem becomes a coefficient lookup after a second convolution.
A naive evaluation of all \(a_s\) from the formula \(\sum_{i=0}^{s}\binom{s}{i}^{16}\) would cost
$$\sum_{s=0}^{k} O(s)=O(k^2),$$
which is too slow when \(k=100000\). Both required convolutions therefore use the Number Theoretic Transform. The target modulus \(M=10^9+7\) is not suitable for an NTT of the needed lengths, so the implementation performs each convolution under three NTT-friendly primes:
$$998244353,\qquad 1004535809,\qquad 469762049.$$
Each coefficient is first computed modulo these three primes and then reconstructed by the Chinese Remainder Theorem. The product of the three auxiliary moduli is far larger than the raw coefficient sizes that occur here, so the reconstruction is exact before the final reduction modulo \(M\).
The small checkpoint \(B(2,4)\) illustrates every layer of the formula. Here
$$q=\left\lfloor\frac{4}{2}\right\rfloor=2,\qquad \mu=2^2-2=2.$$
Next compute the first three values of \(a_s\):
$$a_0=\binom{0}{0}^{16}=1,$$
$$a_1=\binom{1}{0}^{16}+\binom{1}{1}^{16}=1+1=2,$$
$$a_2=\binom{2}{0}^{16}+\binom{2}{1}^{16}+\binom{2}{2}^{16}=1+2^{16}+1=65538.$$
Substituting into the closed form gives
$$\begin{aligned} B(2,4) &=\binom{2}{0}a_0\mu^2+\binom{2}{1}a_1\mu+\binom{2}{2}a_2\\ &=1\cdot 1\cdot 4+2\cdot 2\cdot 2+1\cdot 65538\\ &=4+8+65538\\ &=65550. \end{aligned}$$
This matches the small-value check used by the implementation.
The C++, Python, and Java implementations follow the same pipeline. First they precompute factorials and inverse factorials modulo \(M\) up to \(k\). That makes it cheap to evaluate \((i!)^{-16}\), \((s!)^{16}\), and the factorial denominators that appear in the two generating series.
Next they compute \(\mu=2^{\lfloor n/k\rfloor}-2\pmod{M}\) with fast modular exponentiation, reducing the exponent modulo \(M-1\) before the power is taken. Then they build the sequence \(c_i=(i!)^{-16}\), convolve it with itself, and rescale the result to recover every value \(a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}\).
After that, the implementation forms one series from \(a_s/s!\) and a second series from \(\mu^j/j!\), performs a second convolution, and multiplies the coefficient of degree \(k\) by \(k!\). The C++ version additionally parallelizes the three auxiliary-modulus convolutions and the CRT recombination, while the Python and Java versions execute the same mathematical steps sequentially.
Precomputing factorials, inverse factorials, and the powers needed for the exponential series costs \(O(k)\) time and \(O(k)\) memory. Each NTT-based convolution costs \(O(k\log k)\) time. There are two logical convolutions, each evaluated under a fixed set of three auxiliary primes, so the asymptotic running time remains \(O(k\log k)\) and the memory usage remains \(O(k)\). Parallel execution in the C++ implementation improves wall-clock time but does not change the asymptotic order.
Gesucht ist \(B(k,n)\) modulo \(M=10^9+7\). Die C++-, Python- und Java-Implementierungen durchsuchen keine Matrixkonfigurationen direkt, sondern formen die Zielgröße in ein Koeffizientenproblem für faktoriell gewichtete Reihen um. Der teure Teil ist eine Familie von Summen mit sechzehnten Potenzen von Binomialkoeffizienten, und genau diese Summen werden per schneller Polynomfaltung in einem Schritt berechnet.
Für die Zielparameter ist \(k\) so groß wie \(100000\). Eine naive quadratische Berechnung wäre daher unbrauchbar. Die Lösung reduziert das Problem auf zwei Faltungen von Länge ungefähr \(k\) und eine anschließende modulare Rekonstruktion.
Wir schreiben
$$q=\left\lfloor\frac{n}{k}\right\rfloor,\qquad M=10^9+7,\qquad \mu\equiv 2^q-2 \pmod{M}.$$
Die implementierte Formel lässt sich in den folgenden Schritten herleiten.
Die gesamte Abhängigkeit von \(n\) steckt in dem Skalar \(\mu\). Weil \(M\) prim ist und \(2\not\equiv 0\pmod{M}\), gilt nach dem kleinen Satz von Fermat
$$2^q \equiv 2^{\,q\bmod (M-1)} \pmod{M}.$$
Selbst für sehr große Werte von \(q\) braucht die Implementierung also nur den Rest \(q\bmod(M-1)\). Danach hängt der übrige Rechenweg nur noch von \(k\) und dem bereits bestimmten Wert \(\mu\) ab.
Definiere die Folge
$$c_i=(i!)^{-16}\pmod{M}\qquad (0\le i\le k).$$
Ihre Selbstfaltung ist
$$d_s=\sum_{i=0}^{s} c_i\,c_{s-i} \qquad (0\le s\le k).$$
Multipliziert man mit \((s!)^{16}\), erhält man
$$a_s=(s!)^{16}d_s=(s!)^{16}\sum_{i=0}^{s}\frac{1}{(i!)^{16}((s-i)!)^{16}}=\sum_{i=0}^{s}\left(\frac{s!}{i!(s-i)!}\right)^{16}.$$
Also gilt
$$\boxed{a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}}.$$
Das ist die zentrale Umformung: Statt jedes \(a_s\) separat zu berechnen, liefert eine einzige Selbstfaltung alle Werte \(a_0,a_1,\dots,a_k\) gleichzeitig.
Nachdem die Werte \(a_s\) bekannt sind, setzen wir
$$g_s=\frac{a_s}{s!},\qquad e_j=\frac{\mu^j}{j!}.$$
Die Faltung dieser beiden Folgen ist
$$f_t=\sum_{s=0}^{t} g_s\,e_{t-s}=\sum_{s=0}^{t}\frac{a_s}{s!}\frac{\mu^{\,t-s}}{(t-s)!}.$$
Die gesuchte Größe ist dann
$$B(k,n)=k!\,f_k.$$
Durch Ausmultiplizieren der Fakultäten folgt die explizite Summenformel
$$\boxed{B(k,n)=\sum_{s=0}^{k}\binom{k}{s}a_s\,\mu^{\,k-s}\pmod{M}.}$$
In Sprache der exponentiellen erzeugenden Funktionen ist das gleichbedeutend mit
$$\frac{B(k,n)}{k!}=\left[x^k\right]\left(\sum_{s\ge 0}a_s\frac{x^s}{s!}\right)e^{\mu x}.$$
Das Problem ist also nach der ersten Umformung nur noch ein weiterer Koeffizientenabruf nach einer zweiten Faltung.
Eine naive Berechnung aller \(a_s\) aus \(\sum_{i=0}^{s}\binom{s}{i}^{16}\) kostet insgesamt
$$\sum_{s=0}^{k} O(s)=O(k^2),$$
und ist bei \(k=100000\) nicht praktikabel. Deshalb verwenden beide nötigen Faltungen die Number Theoretic Transform. Der Zielmodulus \(M=10^9+7\) ist für die benötigten Längen nicht NTT-freundlich, also wird jede Faltung unter drei geeigneten Primzahlen ausgeführt:
$$998244353,\qquad 1004535809,\qquad 469762049.$$
Die Koeffizienten werden zunächst unter diesen drei Moduli berechnet und danach per chinesischem Restsatz zusammengesetzt. Das Produkt der drei Hilfsmoduli ist deutlich größer als die hier auftretenden Rohkoeffizienten, daher ist die Rekonstruktion vor der abschließenden Reduktion modulo \(M\) exakt.
Der kleine Kontrollwert \(B(2,4)\) zeigt die Formel vollständig in Aktion. Hier gilt
$$q=\left\lfloor\frac{4}{2}\right\rfloor=2,\qquad \mu=2^2-2=2.$$
Die ersten drei Werte von \(a_s\) sind
$$a_0=\binom{0}{0}^{16}=1,$$
$$a_1=\binom{1}{0}^{16}+\binom{1}{1}^{16}=1+1=2,$$
$$a_2=\binom{2}{0}^{16}+\binom{2}{1}^{16}+\binom{2}{2}^{16}=1+2^{16}+1=65538.$$
Einsetzen in die geschlossene Formel ergibt
$$\begin{aligned} B(2,4) &=\binom{2}{0}a_0\mu^2+\binom{2}{1}a_1\mu+\binom{2}{2}a_2\\ &=1\cdot 1\cdot 4+2\cdot 2\cdot 2+1\cdot 65538\\ &=4+8+65538\\ &=65550. \end{aligned}$$
Genau dieser Wert erscheint auch in der kleinen Kontrolle der Implementierung.
Die C++-, Python- und Java-Implementierungen nutzen denselben Ablauf. Zuerst werden Fakultäten und inverse Fakultäten modulo \(M\) bis \(k\) vorab berechnet. Damit lassen sich \((i!)^{-16}\), \((s!)^{16}\) und die Fakultätsnenner der beiden Reihen sehr schnell auswerten.
Anschließend wird \(\mu=2^{\lfloor n/k\rfloor}-2\pmod{M}\) per schneller modularer Potenzierung berechnet, wobei der Exponent vorher modulo \(M-1\) reduziert wird. Danach wird die Folge \(c_i=(i!)^{-16}\) aufgebaut, mit sich selbst gefaltet und passend reskaliert, um alle Werte \(a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}\) zu gewinnen.
Im letzten Teil wird eine Reihe aus \(a_s/s!\) und eine zweite Reihe aus \(\mu^j/j!\) gebildet. Nach einer weiteren Faltung wird der Koeffizient vom Grad \(k\) mit \(k!\) multipliziert. Die C++-Version parallelisiert zusätzlich die drei Hilfsfaltungen und die CRT-Rekombination; Python und Java führen dieselbe Mathematik sequentiell aus.
Die Vorberechnung von Fakultäten, inversen Fakultäten und den Potenzen für die Exponentialreihe kostet \(O(k)\) Zeit und \(O(k)\) Speicher. Jede NTT-basierte Faltung benötigt \(O(k\log k)\) Zeit. Es gibt zwei logische Faltungen, die jeweils unter einer festen Menge von drei Hilfsmoduli ausgeführt werden, daher bleibt die asymptotische Laufzeit \(O(k\log k)\) und der Speicherverbrauch \(O(k)\). Parallelisierung in C++ verbessert die Laufzeit in der Praxis, ändert aber nicht die asymptotische Ordnung.
Hesaplanması gereken büyüklük, \(M=10^9+7\) modunda \(B(k,n)\) değeridir. C++, Python ve Java uygulamaları matris durumlarını doğrudan saymaz. Bunun yerine hedef değer, faktöriyel ağırlıklı serilerden elde edilen bir katsayı çıkarımı problemine dönüştürülür. Zor kısım, binom katsayılarının on altıncı kuvvetlerini içeren toplamlar ailesidir; uygulama bu toplamların hepsini hızlı polinom konvolüsyonu ile aynı anda hesaplar.
Hedef örnekte \(k=100000\) olduğu için karesel zamanlı doğrudan toplama pratik değildir. Bu nedenle tüm yaklaşım, yaklaşık \(k\) uzunluklu iki konvolüsyon ve sonrasında yapılan modüler yeniden kurulum etrafında kuruludur.
Şöyle yazalım:
$$q=\left\lfloor\frac{n}{k}\right\rfloor,\qquad M=10^9+7,\qquad \mu\equiv 2^q-2 \pmod{M}.$$
Uygulanan formül aşağıdaki adımlarla düzenli biçimde türetilir.
\(n\)'ye bağlı tüm bilgi tek bir skaler olan \(\mu\) içinde toplanır. \(M\) asal olduğu ve \(2\not\equiv 0\pmod{M}\) olduğu için Fermat'nın küçük teoremi
$$2^q \equiv 2^{\,q\bmod (M-1)} \pmod{M}$$
sonucunu verir. Dolayısıyla \(q\) çok büyük olsa bile uygulamanın ihtiyacı olan tek şey \(q\bmod(M-1)\) değeridir. Bu indirgemeden sonra hesabın geri kalanı yalnızca \(k\) ve önceden bulunmuş \(\mu\) değeriyle yürür.
Şu diziyi tanımlayalım:
$$c_i=(i!)^{-16}\pmod{M}\qquad (0\le i\le k).$$
Bu dizinin kendisiyle konvolüsyonu
$$d_s=\sum_{i=0}^{s} c_i\,c_{s-i} \qquad (0\le s\le k)$$
olsun. Bunu \((s!)^{16}\) ile çarptığımızda
$$a_s=(s!)^{16}d_s=(s!)^{16}\sum_{i=0}^{s}\frac{1}{(i!)^{16}((s-i)!)^{16}}=\sum_{i=0}^{s}\left(\frac{s!}{i!(s-i)!}\right)^{16}$$
elde edilir. Yani
$$\boxed{a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}}.$$
Temel sadeleşme budur: her \(a_s\) değerini ayrı ayrı hesaplamak yerine tek bir öz-konvolüsyon ile \(a_0,a_1,\dots,a_k\) değerlerinin tamamı aynı anda elde edilir.
\(a_s\) değerleri bilindikten sonra
$$g_s=\frac{a_s}{s!},\qquad e_j=\frac{\mu^j}{j!}$$
tanımlarını yapalım. Bu iki dizinin konvolüsyonu
$$f_t=\sum_{s=0}^{t} g_s\,e_{t-s}=\sum_{s=0}^{t}\frac{a_s}{s!}\frac{\mu^{\,t-s}}{(t-s)!}$$
olur. Son cevap ise
$$B(k,n)=k!\,f_k$$
şeklindedir. Faktöriyeller açıldığında açık toplam formu
$$\boxed{B(k,n)=\sum_{s=0}^{k}\binom{k}{s}a_s\,\mu^{\,k-s}\pmod{M}.}$$
Aynı ifade üstel üretici fonksiyon diliyle
$$\frac{B(k,n)}{k!}=\left[x^k\right]\left(\sum_{s\ge 0}a_s\frac{x^s}{s!}\right)e^{\mu x}$$
biçiminde de yazılabilir. Böylece bütün problem ikinci bir konvolüsyon sonrası tek bir katsayının okunmasına indirgenir.
\(a_s\) değerlerini doğrudan \(\sum_{i=0}^{s}\binom{s}{i}^{16}\) formülünden hesaplamak toplamda
$$\sum_{s=0}^{k} O(s)=O(k^2)$$
zaman alır ve \(k=100000\) için çok yavaştır. Bu yüzden gereken iki konvolüsyon da Number Theoretic Transform ile yapılır. Hedef modül \(M=10^9+7\), gerekli uzunluklar için NTT dostu değildir; bu nedenle her konvolüsyon üç ayrı uygun asal altında hesaplanır:
$$998244353,\qquad 1004535809,\qquad 469762049.$$
Her katsayı önce bu üç mod altında bulunur, ardından Çin Kalan Teoremi ile yeniden kurulur. Yardımcı modüllerin çarpımı burada oluşan ham katsayı büyüklüklerinden çok daha büyük olduğundan, son olarak \(M\)'ye göre indirgemeden önce yeniden kurulum tam olarak doğrudur.
Küçük kontrol değeri \(B(2,4)\), formülün bütün katmanlarını gösterir. Bu durumda
$$q=\left\lfloor\frac{4}{2}\right\rfloor=2,\qquad \mu=2^2-2=2.$$
İlk üç \(a_s\) değeri şöyledir:
$$a_0=\binom{0}{0}^{16}=1,$$
$$a_1=\binom{1}{0}^{16}+\binom{1}{1}^{16}=1+1=2,$$
$$a_2=\binom{2}{0}^{16}+\binom{2}{1}^{16}+\binom{2}{2}^{16}=1+2^{16}+1=65538.$$
Kapalı forma yerleştirince
$$\begin{aligned} B(2,4) &=\binom{2}{0}a_0\mu^2+\binom{2}{1}a_1\mu+\binom{2}{2}a_2\\ &=1\cdot 1\cdot 4+2\cdot 2\cdot 2+1\cdot 65538\\ &=4+8+65538\\ &=65550 \end{aligned}$$
elde edilir. Bu da uygulamanın küçük doğrulama kontrolüyle aynı değerdir.
C++, Python ve Java uygulamaları aynı işlem hattını kullanır. Önce \(k\)'ya kadar olan faktöriyeller ve ters faktöriyeller \(M\) modunda önceden hesaplanır. Böylece \((i!)^{-16}\), \((s!)^{16}\) ve iki serideki faktöriyel paydalar ucuz hale gelir.
Daha sonra \(\mu=2^{\lfloor n/k\rfloor}-2\pmod{M}\) hızlı modüler üs alma ile hesaplanır; bu sırada üs önce \(M-1\) moduna göre indirgenir. Ardından \(c_i=(i!)^{-16}\) dizisi kurulur, kendisiyle konvolüsyon alınır ve uygun ölçekleme yapılarak tüm \(a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}\) değerleri çıkarılır.
Son aşamada bir seri \(a_s/s!\) katsayılarından, diğer seri ise \(\mu^j/j!\) katsayılarından oluşturulur. İkinci konvolüsyondan sonra derece \(k\) katsayısı \(k!\) ile çarpılır. C++ sürümü buna ek olarak üç yardımcı moddaki konvolüsyonları ve CRT birleştirmesini paralel yürütür; Python ve Java aynı matematiği sıralı biçimde uygular.
Faktöriyellerin, ters faktöriyellerin ve üstel seri için gereken kuvvetlerin ön hesaplaması \(O(k)\) zaman ve \(O(k)\) bellek gerektirir. Her NTT tabanlı konvolüsyon \(O(k\log k)\) zamandadır. İki mantıksal konvolüsyon vardır ve her biri sabit sayıda üç yardımcı mod altında yapılır; bu nedenle toplam asimptotik süre \(O(k\log k)\), bellek kullanımı ise \(O(k)\) olarak kalır. C++ içindeki paralelleştirme pratik çalışma süresini azaltır ama asimptotik düzeni değiştirmez.
La cantidad que se debe calcular es \(B(k,n)\) módulo \(M=10^9+7\). Las implementaciones en C++, Python y Java no enumeran directamente estados de la matriz. En su lugar, reescriben la cantidad objetivo como un problema de extracción de coeficiente construido a partir de series con factoriales. La parte costosa es una familia de sumas con potencias decimosextas de coeficientes binomiales, y esas sumas se calculan todas a la vez mediante convoluciones rápidas de polinomios.
En la instancia objetivo, \(k\) llega a \(100000\), así que una suma cuadrática directa sería inviable. Por eso toda la estrategia se apoya en dos convoluciones de longitud aproximada \(k\) y en una reconstrucción modular al final.
Escribimos
$$q=\left\lfloor\frac{n}{k}\right\rfloor,\qquad M=10^9+7,\qquad \mu\equiv 2^q-2 \pmod{M}.$$
La fórmula implementada se organiza de forma natural en los pasos siguientes.
Toda la dependencia respecto de \(n\) queda concentrada en el escalar \(\mu\). Como \(M\) es primo y \(2\not\equiv 0\pmod{M}\), el pequeño teorema de Fermat da
$$2^q \equiv 2^{\,q\bmod (M-1)} \pmod{M}.$$
Así, incluso cuando \(q\) es enorme, la implementación solo necesita el residuo \(q\bmod(M-1)\). Después de esa reducción, el resto del cálculo depende únicamente de \(k\) y del valor ya obtenido de \(\mu\).
Definimos la sucesión
$$c_i=(i!)^{-16}\pmod{M}\qquad (0\le i\le k).$$
Su autoconvolución es
$$d_s=\sum_{i=0}^{s} c_i\,c_{s-i} \qquad (0\le s\le k).$$
Multiplicando por \((s!)^{16}\) se obtiene
$$a_s=(s!)^{16}d_s=(s!)^{16}\sum_{i=0}^{s}\frac{1}{(i!)^{16}((s-i)!)^{16}}=\sum_{i=0}^{s}\left(\frac{s!}{i!(s-i)!}\right)^{16}.$$
Por tanto,
$$\boxed{a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}}.$$
Esta es la simplificación clave: en vez de evaluar cada \(a_s\) por separado, una sola autoconvolución produce simultáneamente todos los valores \(a_0,a_1,\dots,a_k\).
Una vez conocidos los \(a_s\), definimos
$$g_s=\frac{a_s}{s!},\qquad e_j=\frac{\mu^j}{j!}.$$
La convolución de estas dos sucesiones es
$$f_t=\sum_{s=0}^{t} g_s\,e_{t-s}=\sum_{s=0}^{t}\frac{a_s}{s!}\frac{\mu^{\,t-s}}{(t-s)!}.$$
La respuesta final es
$$B(k,n)=k!\,f_k.$$
Al expandir los factoriales aparece la fórmula explícita
$$\boxed{B(k,n)=\sum_{s=0}^{k}\binom{k}{s}a_s\,\mu^{\,k-s}\pmod{M}.}$$
De manera equivalente, en términos de funciones generadoras exponenciales,
$$\frac{B(k,n)}{k!}=\left[x^k\right]\left(\sum_{s\ge 0}a_s\frac{x^s}{s!}\right)e^{\mu x}.$$
Así, todo el problema se reduce a leer un solo coeficiente tras una segunda convolución.
Calcular todos los \(a_s\) directamente desde \(\sum_{i=0}^{s}\binom{s}{i}^{16}\) cuesta en total
$$\sum_{s=0}^{k} O(s)=O(k^2),$$
lo cual es demasiado lento para \(k=100000\). Por eso las dos convoluciones necesarias se realizan mediante Number Theoretic Transform. El módulo objetivo \(M=10^9+7\) no es adecuado para una NTT de las longitudes requeridas, de modo que cada convolución se calcula bajo tres primos aptos:
$$998244353,\qquad 1004535809,\qquad 469762049.$$
Cada coeficiente se obtiene primero módulo esos tres primos y luego se recompone mediante el Teorema Chino del Resto. El producto de los tres módulos auxiliares es mucho mayor que los coeficientes enteros reales que aparecen aquí, así que la reconstrucción es exacta antes de la reducción final módulo \(M\).
El valor pequeño \(B(2,4)\) muestra la fórmula completa en acción. En este caso
$$q=\left\lfloor\frac{4}{2}\right\rfloor=2,\qquad \mu=2^2-2=2.$$
Los tres primeros valores de \(a_s\) son
$$a_0=\binom{0}{0}^{16}=1,$$
$$a_1=\binom{1}{0}^{16}+\binom{1}{1}^{16}=1+1=2,$$
$$a_2=\binom{2}{0}^{16}+\binom{2}{1}^{16}+\binom{2}{2}^{16}=1+2^{16}+1=65538.$$
Al sustituir en la forma cerrada se obtiene
$$\begin{aligned} B(2,4) &=\binom{2}{0}a_0\mu^2+\binom{2}{1}a_1\mu+\binom{2}{2}a_2\\ &=1\cdot 1\cdot 4+2\cdot 2\cdot 2+1\cdot 65538\\ &=4+8+65538\\ &=65550. \end{aligned}$$
Ese resultado coincide con la comprobación pequeña usada por la implementación.
Las implementaciones en C++, Python y Java siguen la misma tubería. Primero precalculan factoriales e inversos factoriales módulo \(M\) hasta \(k\). Eso permite evaluar con rapidez \((i!)^{-16}\), \((s!)^{16}\) y los denominadores factoriales de las dos series.
Después calculan \(\mu=2^{\lfloor n/k\rfloor}-2\pmod{M}\) mediante exponenciación modular rápida, reduciendo antes el exponente módulo \(M-1\). Luego construyen la sucesión \(c_i=(i!)^{-16}\), la autoconvolucionan y reescalan el resultado para recuperar todos los valores \(a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}\).
Finalmente, la implementación forma una serie con coeficientes \(a_s/s!\) y otra con coeficientes \(\mu^j/j!\), realiza una segunda convolución y multiplica por \(k!\) el coeficiente de grado \(k\). La versión en C++ además paraleliza las tres convoluciones en módulos auxiliares y la recombinación por CRT; las versiones en Python y Java ejecutan los mismos pasos de forma secuencial.
El precálculo de factoriales, inversos factoriales y potencias necesarias para la serie exponencial cuesta \(O(k)\) tiempo y \(O(k)\) memoria. Cada convolución basada en NTT cuesta \(O(k\log k)\) tiempo. Hay dos convoluciones lógicas, cada una evaluada bajo un conjunto fijo de tres primos auxiliares, de modo que el tiempo total sigue siendo \(O(k\log k)\) y la memoria total \(O(k)\). La paralelización en C++ mejora el tiempo real de ejecución, pero no cambia el orden asintótico.
需要计算的是 \(B(k,n)\) 在模 \(M=10^9+7\) 下的值。C++、Python 和 Java 实现都没有直接去枚举矩阵状态，而是把目标量改写成一个由阶乘型级数组成的系数提取问题。真正昂贵的部分，是一整族包含二项式系数十六次幂的求和；实现的核心是用快速多项式卷积一次性把这些值全部算出来。
在目标参数下，\(k\) 达到 \(100000\)，因此任何 \(O(k^2)\) 级别的直接求和都不可行。整个方法于是被压缩成两次长度约为 \(k\) 的卷积，再加上最后一步模意义下的重建。
记
$$q=\left\lfloor\frac{n}{k}\right\rfloor,\qquad M=10^9+7,\qquad \mu\equiv 2^q-2 \pmod{M}.$$
实现中的公式可以按下面几个步骤理解。
关于 \(n\) 的全部依赖，都浓缩在单个标量 \(\mu\) 中。因为 \(M\) 是素数，而且 \(2\not\equiv 0\pmod{M}\)，由费马小定理可得
$$2^q \equiv 2^{\,q\bmod (M-1)} \pmod{M}.$$
这意味着即使 \(q\) 极大，实现也只需要 \(q\bmod(M-1)\) 这个余数，而不需要处理完整的指数。完成这一步之后，后续计算只依赖于 \(k\) 和已经得到的 \(\mu\)。
定义数列
$$c_i=(i!)^{-16}\pmod{M}\qquad (0\le i\le k).$$
把它与自身卷积，得到
$$d_s=\sum_{i=0}^{s} c_i\,c_{s-i}\qquad (0\le s\le k).$$
然后再乘上 \((s!)^{16}\)，就有
$$a_s=(s!)^{16}d_s=(s!)^{16}\sum_{i=0}^{s}\frac{1}{(i!)^{16}((s-i)!)^{16}}=\sum_{i=0}^{s}\left(\frac{s!}{i!(s-i)!}\right)^{16}.$$
也就是说
$$\boxed{a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}}.$$
这是最关键的化简。原本如果对每个 \(s\) 单独求这个和，总成本会很高；而现在一次自卷积就能同时得到 \(a_0,a_1,\dots,a_k\) 全部值。
在 \(a_s\) 已知后，定义
$$g_s=\frac{a_s}{s!},\qquad e_j=\frac{\mu^j}{j!}.$$
这两个数列的卷积为
$$f_t=\sum_{s=0}^{t} g_s\,e_{t-s}=\sum_{s=0}^{t}\frac{a_s}{s!}\frac{\mu^{\,t-s}}{(t-s)!}.$$
最终答案正是
$$B(k,n)=k!\,f_k.$$
把阶乘整理开以后，就得到显式公式
$$\boxed{B(k,n)=\sum_{s=0}^{k}\binom{k}{s}a_s\,\mu^{\,k-s}\pmod{M}.}$$
如果改用指数型生成函数的语言，它等价于
$$\frac{B(k,n)}{k!}=\left[x^k\right]\left(\sum_{s\ge 0}a_s\frac{x^s}{s!}\right)e^{\mu x}.$$
因此整个问题本质上变成了第二次卷积之后读取一个单独的系数。
如果直接根据 \(\sum_{i=0}^{s}\binom{s}{i}^{16}\) 去计算全部 \(a_s\)，总代价是
$$\sum_{s=0}^{k} O(s)=O(k^2),$$
对 \(k=100000\) 来说过于缓慢。因此两次所需卷积都用 Number Theoretic Transform 完成。目标模数 \(M=10^9+7\) 并不适合所需长度的 NTT，所以实现改在三个 NTT 友好的素数下分别计算：
$$998244353,\qquad 1004535809,\qquad 469762049.$$
每个系数先在这三个模数下求出，再用中国剩余定理合并。由于这三个辅助模数的乘积远大于这里出现的原始整数系数范围，所以在最后再模 \(M\) 之前，重建出来的系数就是精确值。
小规模检查值 \(B(2,4)\) 可以把整套公式完整展示出来。此时
$$q=\left\lfloor\frac{4}{2}\right\rfloor=2,\qquad \mu=2^2-2=2.$$
前三个 \(a_s\) 为
$$a_0=\binom{0}{0}^{16}=1,$$
$$a_1=\binom{1}{0}^{16}+\binom{1}{1}^{16}=1+1=2,$$
$$a_2=\binom{2}{0}^{16}+\binom{2}{1}^{16}+\binom{2}{2}^{16}=1+2^{16}+1=65538.$$
代回闭式可得
$$\begin{aligned} B(2,4) &=\binom{2}{0}a_0\mu^2+\binom{2}{1}a_1\mu+\binom{2}{2}a_2\\ &=1\cdot 1\cdot 4+2\cdot 2\cdot 2+1\cdot 65538\\ &=4+8+65538\\ &=65550. \end{aligned}$$
这正好与实现中的小规模校验值一致。
C++、Python 和 Java 实现遵循同一条计算流水线。首先，它们会预先计算到 \(k\) 为止的阶乘与逆阶乘模 \(M\) 的值。这样就能高效得到 \((i!)^{-16}\)、\((s!)^{16}\)，以及两组级数中出现的所有阶乘分母。
随后，程序通过快速模幂计算 \(\mu=2^{\lfloor n/k\rfloor}-2\pmod{M}\)，并且先把指数对 \(M-1\) 取模。接着构造数列 \(c_i=(i!)^{-16}\)，做一次自卷积，再经过适当的阶乘缩放，恢复出全部 \(a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}\)。
最后，程序再构造一组系数为 \(a_s/s!\) 的级数，以及另一组系数为 \(\mu^j/j!\) 的级数，做第二次卷积，并把 \(k\) 次项乘上 \(k!\)。其中 C++ 版本还会并行处理三个辅助模数下的卷积与 CRT 合并；Python 和 Java 版本则按同样的数学步骤顺序执行。
阶乘、逆阶乘以及指数级数所需幂的预计算需要 \(O(k)\) 时间和 \(O(k)\) 空间。每次基于 NTT 的卷积需要 \(O(k\log k)\) 时间。整个算法在逻辑上有两次卷积，每次都只是在固定的三个辅助素数下重复，因此总的渐近时间复杂度仍是 \(O(k\log k)\)，空间复杂度是 \(O(k)\)。C++ 中的并行化只改善实际运行时间，不改变渐近阶。
Требуется вычислить \(B(k,n)\) по модулю \(M=10^9+7\). Реализации на C++, Python и Java не перебирают состояния матрицы напрямую. Вместо этого целевая величина переписывается как задача извлечения коэффициента из рядов с факториальными весами. Самая дорогая часть состоит из семейства сумм с шестнадцатыми степенями биномиальных коэффициентов, и именно эти суммы вычисляются сразу для всех значений при помощи быстрой свертки полиномов.
В целевом случае \(k=100000\), поэтому квадратичный прямой подсчет был бы непригоден. Вся схема решения сводится к двум сверткам длины порядка \(k\) и заключительной реконструкции коэффициентов по модулю.
Запишем
$$q=\left\lfloor\frac{n}{k}\right\rfloor,\qquad M=10^9+7,\qquad \mu\equiv 2^q-2 \pmod{M}.$$
Формула, реализованная в программе, естественно раскладывается на следующие шаги.
Вся зависимость от \(n\) сосредоточена в одном скаляре \(\mu\). Так как \(M\) простое и \(2\not\equiv 0\pmod{M}\), малая теорема Ферма дает
$$2^q \equiv 2^{\,q\bmod (M-1)} \pmod{M}.$$
Значит, даже при огромном \(q\) реализации нужен только остаток \(q\bmod(M-1)\). После этого шага остальная часть вычисления зависит лишь от \(k\) и уже найденного значения \(\mu\).
Определим последовательность
$$c_i=(i!)^{-16}\pmod{M}\qquad (0\le i\le k).$$
Ее самосвертка имеет вид
$$d_s=\sum_{i=0}^{s} c_i\,c_{s-i}\qquad (0\le s\le k).$$
Если умножить результат на \((s!)^{16}\), получим
$$a_s=(s!)^{16}d_s=(s!)^{16}\sum_{i=0}^{s}\frac{1}{(i!)^{16}((s-i)!)^{16}}=\sum_{i=0}^{s}\left(\frac{s!}{i!(s-i)!}\right)^{16}.$$
Следовательно,
$$\boxed{a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}}.$$
Это ключевое упрощение: вместо отдельного вычисления каждого \(a_s\) одна самосвертка сразу дает все значения \(a_0,a_1,\dots,a_k\).
Когда значения \(a_s\) уже известны, вводим
$$g_s=\frac{a_s}{s!},\qquad e_j=\frac{\mu^j}{j!}.$$
Свертка этих последовательностей равна
$$f_t=\sum_{s=0}^{t} g_s\,e_{t-s}=\sum_{s=0}^{t}\frac{a_s}{s!}\frac{\mu^{\,t-s}}{(t-s)!}.$$
Искомая величина получается как
$$B(k,n)=k!\,f_k.$$
После раскрытия факториалов выходит явная формула
$$\boxed{B(k,n)=\sum_{s=0}^{k}\binom{k}{s}a_s\,\mu^{\,k-s}\pmod{M}.}$$
То же самое можно записать через экспоненциальную производящую функцию:
$$\frac{B(k,n)}{k!}=\left[x^k\right]\left(\sum_{s\ge 0}a_s\frac{x^s}{s!}\right)e^{\mu x}.$$
Тем самым вся задача сводится к чтению одного коэффициента после второй свертки.
Если вычислять все \(a_s\) напрямую по формуле \(\sum_{i=0}^{s}\binom{s}{i}^{16}\), суммарная стоимость составит
$$\sum_{s=0}^{k} O(s)=O(k^2),$$
что слишком медленно при \(k=100000\). Поэтому обе нужные свертки выполняются через Number Theoretic Transform. Целевой модуль \(M=10^9+7\) не подходит для NTT нужной длины, поэтому каждая свертка считается под тремя удобными простыми:
$$998244353,\qquad 1004535809,\qquad 469762049.$$
Сначала каждый коэффициент вычисляется по этим трем модулям, а затем восстанавливается по китайской теореме об остатках. Произведение трех вспомогательных модулей существенно больше реальных целых коэффициентов, возникающих в этой задаче, поэтому реконструкция точна еще до окончательного приведения по модулю \(M\).
Малое контрольное значение \(B(2,4)\) хорошо показывает всю формулу. Здесь
$$q=\left\lfloor\frac{4}{2}\right\rfloor=2,\qquad \mu=2^2-2=2.$$
Первые три значения \(a_s\) равны
$$a_0=\binom{0}{0}^{16}=1,$$
$$a_1=\binom{1}{0}^{16}+\binom{1}{1}^{16}=1+1=2,$$
$$a_2=\binom{2}{0}^{16}+\binom{2}{1}^{16}+\binom{2}{2}^{16}=1+2^{16}+1=65538.$$
Подставляя это в закрытую формулу, получаем
$$\begin{aligned} B(2,4) &=\binom{2}{0}a_0\mu^2+\binom{2}{1}a_1\mu+\binom{2}{2}a_2\\ &=1\cdot 1\cdot 4+2\cdot 2\cdot 2+1\cdot 65538\\ &=4+8+65538\\ &=65550. \end{aligned}$$
Это совпадает с малой проверкой, включенной в реализацию.
Реализации на C++, Python и Java используют один и тот же вычислительный конвейер. Сначала они заранее считают факториалы и обратные факториалы по модулю \(M\) до \(k\). Благодаря этому можно быстро получать \((i!)^{-16}\), \((s!)^{16}\) и факториальные знаменатели двух рядов.
Затем они вычисляют \(\mu=2^{\lfloor n/k\rfloor}-2\pmod{M}\) с помощью быстрого возведения в степень по модулю, предварительно уменьшая показатель по модулю \(M-1\). После этого строится последовательность \(c_i=(i!)^{-16}\), выполняется самосвертка и применяется нужное факториальное масштабирование, чтобы восстановить все значения \(a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}\).
На финальном этапе формируется один ряд с коэффициентами \(a_s/s!\) и второй ряд с коэффициентами \(\mu^j/j!\), затем выполняется вторая свертка, а коэффициент при степени \(k\) умножается на \(k!\). Версия на C++ дополнительно распараллеливает свертки по трем вспомогательным модулям и объединение через CRT; версии на Python и Java выполняют те же математические шаги последовательно.
Предварительный расчет факториалов, обратных факториалов и степеней для экспоненциального ряда занимает \(O(k)\) времени и \(O(k)\) памяти. Каждая свертка на основе NTT требует \(O(k\log k)\) времени. Логически здесь есть две свертки, каждая выполняется под фиксированным набором из трех вспомогательных простых, поэтому общая асимптотика остается \(O(k\log k)\) по времени и \(O(k)\) по памяти. Параллелизм в C++ улучшает практическое время выполнения, но не меняет асимптотический порядок.
المطلوب هو حساب \(B(k,n)\) بترديد \(M=10^9+7\). تنفيذات ++C وPython وJava لا تعد حالات المصفوفة مباشرة، بل تعيد صياغة الكمية المطلوبة على أنها مسألة استخراج معامل من سلاسل مبنية على العاملات. الجزء المكلف فعلا هو عائلة من المجاميع التي تحتوي على القوى السادسة عشرة لمعاملات ثنائية الحد، ولهذا تحسب هذه القيم كلها دفعة واحدة باستعمال التفاف كثيرات حدود سريع.
في المسألة الهدف يصل \(k\) إلى \(100000\)، ولذلك فإن أي جمع تربيعي مباشر سيكون بطيئا جدا. لهذا تنظم الفكرة كلها حول التفافين بطول يقارب \(k\)، ثم إعادة تركيب نهائية في الحساب المعياري.
نكتب
$$q=\left\lfloor\frac{n}{k}\right\rfloor,\qquad M=10^9+7,\qquad \mu\equiv 2^q-2 \pmod{M}.$$
ويمكن فهم الصيغة التي تطبقها البرامج عبر الخطوات الآتية.
كل اعتماد على \(n\) يمر عبر العدد المفرد \(\mu\). وبما أن \(M\) أولي و\(2\not\equiv 0\pmod{M}\)، فإن مبرهنة فيرما الصغرى تعطي
$$2^q \equiv 2^{\,q\bmod (M-1)} \pmod{M}.$$
إذن حتى لو كان \(q\) ضخما جدا، فالتنفيذ يحتاج فقط إلى الباقي \(q\bmod(M-1)\). بعد هذه الخطوة يصبح ما تبقى من الحساب معتمدا على \(k\) وعلى القيمة المحسوبة مسبقا لـ \(\mu\) فقط.
نعرّف المتتالية
$$c_i=(i!)^{-16}\pmod{M}\qquad (0\le i\le k).$$
ثم نأخذ التفافها مع نفسها:
$$d_s=\sum_{i=0}^{s} c_i\,c_{s-i}\qquad (0\le s\le k).$$
وعند الضرب في \((s!)^{16}\) نحصل على
$$a_s=(s!)^{16}d_s=(s!)^{16}\sum_{i=0}^{s}\frac{1}{(i!)^{16}((s-i)!)^{16}}=\sum_{i=0}^{s}\left(\frac{s!}{i!(s-i)!}\right)^{16}.$$
ومن ثم
$$\boxed{a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}}.$$
هذا هو التبسيط المركزي: بدلا من حساب كل \(a_s\) على حدة، يعطي التفاف ذاتي واحد كل القيم \(a_0,a_1,\dots,a_k\) دفعة واحدة.
بعد معرفة القيم \(a_s\)، نضع
$$g_s=\frac{a_s}{s!},\qquad e_j=\frac{\mu^j}{j!}.$$
ويكون التفاف هاتين المتتاليتين
$$f_t=\sum_{s=0}^{t} g_s\,e_{t-s}=\sum_{s=0}^{t}\frac{a_s}{s!}\frac{\mu^{\,t-s}}{(t-s)!}.$$
أما الجواب النهائي فهو
$$B(k,n)=k!\,f_k.$$
وبعد ترتيب العاملات نحصل على الصيغة الصريحة
$$\boxed{B(k,n)=\sum_{s=0}^{k}\binom{k}{s}a_s\,\mu^{\,k-s}\pmod{M}.}$$
وبلغة الدوال المولدة الأسية يمكن كتابة ذلك على الصورة
$$\frac{B(k,n)}{k!}=\left[x^k\right]\left(\sum_{s\ge 0}a_s\frac{x^s}{s!}\right)e^{\mu x}.$$
وهكذا تتحول المسألة كلها إلى قراءة معامل واحد بعد التفاف ثان.
الحساب المباشر لكل القيم \(a_s\) من الصيغة \(\sum_{i=0}^{s}\binom{s}{i}^{16}\) يكلف إجمالا
$$\sum_{s=0}^{k} O(s)=O(k^2),$$
وذلك بطيء جدا عندما \(k=100000\). لذلك يعتمد الالتفافان المطلوبان على Number Theoretic Transform. لكن الترديد الهدف \(M=10^9+7\) غير مناسب لـ NTT بالأطوال المطلوبة، ولهذا يحسب كل التفاف تحت ثلاثة أعداد أولية ملائمة:
$$998244353,\qquad 1004535809,\qquad 469762049.$$
يحسب كل معامل أولا تحت هذه الترديدات الثلاثة، ثم يعاد تركيبه بنظرية البواقي الصينية. وبما أن حاصل ضرب هذه المودات المساعدة أكبر بكثير من أحجام المعاملات الصحيحة الخام التي تظهر هنا، فإن إعادة التركيب تكون دقيقة تماما قبل الاختزال النهائي modulo \(M\).
تعطي القيمة الصغيرة \(B(2,4)\) مثالا واضحا على الصيغة كاملة. لدينا هنا
$$q=\left\lfloor\frac{4}{2}\right\rfloor=2,\qquad \mu=2^2-2=2.$$
وأول ثلاث قيم من \(a_s\) هي
$$a_0=\binom{0}{0}^{16}=1,$$
$$a_1=\binom{1}{0}^{16}+\binom{1}{1}^{16}=1+1=2,$$
$$a_2=\binom{2}{0}^{16}+\binom{2}{1}^{16}+\binom{2}{2}^{16}=1+2^{16}+1=65538.$$
وبالتعويض في الصيغة المغلقة نحصل على
$$\begin{aligned} B(2,4) &=\binom{2}{0}a_0\mu^2+\binom{2}{1}a_1\mu+\binom{2}{2}a_2\\ &=1\cdot 1\cdot 4+2\cdot 2\cdot 2+1\cdot 65538\\ &=4+8+65538\\ &=65550. \end{aligned}$$
وهذه هي بالضبط قيمة التحقق الصغيرة التي تعطيها التطبيقات.
تنفيذات ++C وPython وJava تتبع المسار الحسابي نفسه. أولا تحسب العاملات والعاملات العكسية modulo \(M\) حتى \(k\). وهذا يجعل تقييم \((i!)^{-16}\) و\((s!)^{16}\) وجميع المقامات العاملية في السلسلتين أمرا سريعا.
بعد ذلك تحسب \(\mu=2^{\lfloor n/k\rfloor}-2\pmod{M}\) باستخدام الرفع السريع modulo، مع اختزال الأس أولا modulo \(M-1\). ثم تبني المتتالية \(c_i=(i!)^{-16}\)، وتجري التفافا ذاتيا لها، ثم تعيد التحجيم بعوامل مناسبة لاسترجاع جميع القيم \(a_s=\sum_{i=0}^{s}\binom{s}{i}^{16}\).
في المرحلة الأخيرة، تبني الخوارزمية سلسلة معاملاتُها \(a_s/s!\) وسلسلة أخرى معاملاتُها \(\mu^j/j!\)، ثم تجري التفافا ثانيا وتأخذ معامل الدرجة \(k\) وتضربه في \(k!\). نسخة ++C تضيف إلى ذلك تنفيذا متوازيا لالتفافات المودات الثلاثة ولإعادة التركيب عبر CRT، بينما تنفذ نسختا Python وJava الخطوات الرياضية نفسها بشكل متسلسل.
التهيئة المسبقة للعاملات والعاملات العكسية والقوى اللازمة للسلسلة الأسية تكلف \(O(k)\) زمنا و\(O(k)\) ذاكرة. وكل التفاف مبني على NTT يكلف \(O(k\log k)\) زمنا. توجد التفافتان منطقيتان، وكل واحدة تنفذ تحت مجموعة ثابتة من ثلاثة أعداد أولية مساعدة، لذلك يبقى التعقيد الكلي \(O(k\log k)\) والذاكرة \(O(k)\). أما التوازي في نسخة ++C فيحسن الزمن الفعلي فقط ولا يغير الرتبة التقاربية.