There are \(k\) bacterial types. Type \(i\) does not have just one reproduction rule; it has \(m\) equally likely rules extracted from the deterministic pseudo-random sequence
$$r_0=306,\qquad r_{t+1}\equiv r_t^2 \pmod{10007}.$$
The block \(r_{im},r_{im+1},\dots,r_{im+m-1}\) belongs to type \(i\). Starting from a single bacterium of type \(0\), we repeatedly apply these rules independently to all descendants and ask for the probability that the entire population eventually disappears.
This is a multitype branching process. The key simplification is that each type’s \(m\) rules can be compressed into five residue counts, after which the extinction problem becomes a fixed-point system on \([0,1]^k\).
For each type \(i\) and each residue \(q\in\{0,1,2,3,4\}\), define
$$c_{i,q}=\#\left\{j\in\{0,\dots,m-1\}: r_{im+j}\equiv q \pmod{5}\right\}.$$
These counts satisfy
$$\sum_{q=0}^4 c_{i,q}=m.$$
Once the counts are known, the individual rule values are irrelevant. Every bacterium of type \(i\) simply chooses uniformly from a multiset containing \(c_{i,0}\) rules of kind \(0\), \(c_{i,1}\) rules of kind \(1\), and so on.
The five residue classes encode five possible outcomes for one bacterium of type \(i\).
Residue \(0\): no offspring are produced.
Residue \(1\): two offspring of type \(i\) are produced.
Residue \(2\): one offspring of type \(2i \bmod k\) is produced.
Residue \(3\): three offspring of type \((i^2+1) \bmod k\) are produced.
Residue \(4\): one offspring of type \(i\) and one offspring of type \((i+1) \bmod k\) are produced.
If \(s_j\) denotes the extinction contribution of a single child of type \(j\), these five cases contribute the monomials
$$1,\qquad s_i^2,\qquad s_{2i \bmod k},\qquad s_{(i^2+1) \bmod k}^3,\qquad s_i s_{(i+1) \bmod k}.$$
The squares, cubes, and products come from independence: all descendants created by the chosen rule must die out.
Let \(\mathbf{s}=(s_0,\dots,s_{k-1})\). Averaging over the \(m\) equally likely rules of type \(i\) gives the offspring generating polynomial
$$f_i(\mathbf{s})=\frac{c_{i,0}+c_{i,1}s_i^2+c_{i,2}s_{2i \bmod k}+c_{i,3}s_{(i^2+1) \bmod k}^3+c_{i,4}s_i s_{(i+1) \bmod k}}{m}.$$
Now let \(p_i\) be the probability that the process eventually becomes extinct when it starts from one bacterium of type \(i\). Standard branching-process reasoning gives
$$p_i=f_i(\mathbf{p})=\frac{c_{i,0}+c_{i,1}p_i^2+c_{i,2}p_{2i \bmod k}+c_{i,3}p_{(i^2+1) \bmod k}^3+c_{i,4}p_i p_{(i+1) \bmod k}}{m}.$$
So the desired answer is the first coordinate of the extinction vector \(\mathbf{p}\).
Define the map
$$T(\mathbf{x})=\left(f_0(\mathbf{x}),f_1(\mathbf{x}),\dots,f_{k-1}(\mathbf{x})\right).$$
Because all coefficients are non-negative and each row is averaged by \(m\), the map sends \([0,1]^k\) into itself and is monotone in every coordinate. Therefore, if we start from
$$\mathbf{p}^{(0)}=(0,\dots,0),\qquad \mathbf{p}^{(t+1)}=T\left(\mathbf{p}^{(t)}\right),$$
the iterates form an increasing sequence bounded above by \((1,\dots,1)\). They converge to the minimal fixed point of the system, which is exactly the extinction vector for the branching process.
The first four pseudo-random values are
$$306,\qquad 3573,\qquad 7404,\qquad 870,$$
so the residues modulo \(5\) are
$$1,\qquad 3,\qquad 4,\qquad 0.$$
Thus type \(0\) receives one residue-\(1\) rule and one residue-\(3\) rule, while type \(1\) receives one residue-\(4\) rule and one residue-\(0\) rule. The extinction equations become
$$p_0=\frac{p_0^2+p_1^3}{2},\qquad p_1=\frac{1+p_0p_1}{2}.$$
Starting from \((0,0)\), the first synchronous iterates are
$$\begin{aligned} \mathbf{p}^{(1)}&=\left(0,\frac12\right),\\ \mathbf{p}^{(2)}&=\left(\frac1{16},\frac12\right),\\ \mathbf{p}^{(3)}&=\left(\frac{33}{512},\frac{33}{64}\right). \end{aligned}$$
Continuing this process gives
$$p_0\approx 0.07243802,$$
which matches the small verification case used by the implementation.
The C++, Python, and Java implementations all follow the same structure. First they generate the full pseudo-random stream of length \(km\), then they split it into \(k\) consecutive blocks of size \(m\) and tally the five residue counts for each type.
Next they keep two probability vectors. During one sweep over all types, every new coordinate is computed only from the previous sweep’s vector, so the update is synchronous rather than in-place. After the sweep, the vectors are swapped.
Each sweep also records the largest absolute change between old and new coordinates. The C++ implementation uses extended precision and stops once this change is below \(10^{-20}\); the Python and Java implementations use double precision and stop at \(10^{-15}\). All three impose a hard cap of \(10^6\) iterations.
After convergence, the implementation prints the extinction probability of type \(0\) rounded to eight decimal places.
Generating the pseudo-random stream and building the five-count histogram for every type costs \(O(km)\) time. If the fixed-point solver needs \(I\) sweeps, the nonlinear iteration costs \(O(Ik)\) time, so the total running time is
$$O(km+Ik).$$
As written, the implementations store the entire stream of length \(km\), the \(5k\) histogram counts, and two probability vectors of length \(k\). Therefore the memory usage is \(O(km)\).
Es gibt \(k\) Bakterientypen. Jeder Typ \(i\) besitzt nicht nur eine Reproduktionsregel, sondern \(m\) gleichwahrscheinliche Regeln, die aus der deterministischen pseudozufälligen Folge
$$r_0=306,\qquad r_{t+1}\equiv r_t^2 \pmod{10007}$$
gewonnen werden. Der Block \(r_{im},r_{im+1},\dots,r_{im+m-1}\) gehört zu Typ \(i\). Ausgehend von einem einzelnen Bakterium vom Typ \(0\) wird der Prozess unabhängig auf alle Nachkommen angewendet; gesucht ist die Wahrscheinlichkeit, dass die gesamte Population irgendwann ausstirbt.
Das Problem ist ein multityper Verzweigungsprozess. Die zentrale Vereinfachung besteht darin, die \(m\) Regeln jedes Typs zu fünf Häufigkeiten zusammenzufassen. Danach bleibt ein Fixpunktsystem auf \([0,1]^k\).
Für jeden Typ \(i\) und jeden Rest \(q\in\{0,1,2,3,4\}\) definieren wir
$$c_{i,q}=\#\left\{j\in\{0,\dots,m-1\}: r_{im+j}\equiv q \pmod{5}\right\}.$$
Diese Zahlen erfüllen
$$\sum_{q=0}^4 c_{i,q}=m.$$
Sobald diese Häufigkeiten bekannt sind, spielen die einzelnen Regelwerte keine eigene Rolle mehr. Jedes Bakterium vom Typ \(i\) wählt nur noch gleichverteilt aus einem Multiset mit \(c_{i,0}\) Regeln der Art \(0\), \(c_{i,1}\) Regeln der Art \(1\) usw.
Die fünf Restklassen kodieren fünf mögliche Ergebnisse für ein einzelnes Bakterium vom Typ \(i\).
Rest \(0\): Es entstehen keine Nachkommen.
Rest \(1\): Es entstehen zwei Nachkommen vom Typ \(i\).
Rest \(2\): Es entsteht ein Nachkomme vom Typ \(2i \bmod k\).
Rest \(3\): Es entstehen drei Nachkommen vom Typ \((i^2+1) \bmod k\).
Rest \(4\): Es entstehen ein Nachkomme vom Typ \(i\) und ein Nachkomme vom Typ \((i+1) \bmod k\).
Bezeichnet \(s_j\) den Aussterbensbeitrag eines einzelnen Kindes vom Typ \(j\), dann liefern diese fünf Fälle die Monome
$$1,\qquad s_i^2,\qquad s_{2i \bmod k},\qquad s_{(i^2+1) \bmod k}^3,\qquad s_i s_{(i+1) \bmod k}.$$
Quadrate, Kuben und Produkte entstehen aus der Unabhängigkeit der Nachkommenlinien: Alle durch die gewählte Regel erzeugten Teilprozesse müssen aussterben.
Mit \(\mathbf{s}=(s_0,\dots,s_{k-1})\) ergibt das Mittel über die \(m\) gleichwahrscheinlichen Regeln von Typ \(i\) das erzeugende Polynom
$$f_i(\mathbf{s})=\frac{c_{i,0}+c_{i,1}s_i^2+c_{i,2}s_{2i \bmod k}+c_{i,3}s_{(i^2+1) \bmod k}^3+c_{i,4}s_i s_{(i+1) \bmod k}}{m}.$$
Sei nun \(p_i\) die Wahrscheinlichkeit, dass ein von einem einzelnen Bakterium des Typs \(i\) gestarteter Prozess ausstirbt. Dann gilt nach Standardargumenten für Verzweigungsprozesse
$$p_i=f_i(\mathbf{p})=\frac{c_{i,0}+c_{i,1}p_i^2+c_{i,2}p_{2i \bmod k}+c_{i,3}p_{(i^2+1) \bmod k}^3+c_{i,4}p_i p_{(i+1) \bmod k}}{m}.$$
Gesucht ist also die erste Koordinate des Aussterbensvektors \(\mathbf{p}\).
Definiere
$$T(\mathbf{x})=\left(f_0(\mathbf{x}),f_1(\mathbf{x}),\dots,f_{k-1}(\mathbf{x})\right).$$
Da alle Koeffizienten nichtnegativ sind und jede Zeile durch \(m\) gemittelt wird, bildet diese Abbildung \([0,1]^k\) in sich ab und ist koordinatenweise monoton. Starten wir daher mit
$$\mathbf{p}^{(0)}=(0,\dots,0),\qquad \mathbf{p}^{(t+1)}=T\left(\mathbf{p}^{(t)}\right),$$
so entsteht eine wachsende Folge, die nach oben durch \((1,\dots,1)\) beschränkt ist. Sie konvergiert gegen den kleinsten Fixpunkt des Systems, und genau dieser ist der gesuchte Aussterbensvektor.
Die ersten vier pseudozufälligen Werte sind
$$306,\qquad 3573,\qquad 7404,\qquad 870,$$
also lauten die Reste modulo \(5\)
$$1,\qquad 3,\qquad 4,\qquad 0.$$
Damit erhält Typ \(0\) genau eine Regel vom Rest \(1\) und eine vom Rest \(3\), während Typ \(1\) genau eine Regel vom Rest \(4\) und eine vom Rest \(0\) erhält. Das Fixpunktsystem lautet dann
$$p_0=\frac{p_0^2+p_1^3}{2},\qquad p_1=\frac{1+p_0p_1}{2}.$$
Beginnend bei \((0,0)\) ergeben die ersten synchronen Iterationen
$$\begin{aligned} \mathbf{p}^{(1)}&=\left(0,\frac12\right),\\ \mathbf{p}^{(2)}&=\left(\frac1{16},\frac12\right),\\ \mathbf{p}^{(3)}&=\left(\frac{33}{512},\frac{33}{64}\right). \end{aligned}$$
Die Folge konvergiert weiter zu
$$p_0\approx 0.07243802,$$
genau dem kleinen Prüfwert, den die Implementierung verwendet.
Die C++-, Python- und Java-Implementierungen haben denselben Aufbau. Zuerst erzeugen sie die vollständige pseudozufällige Folge der Länge \(km\). Danach zerlegen sie diese in \(k\) aufeinanderfolgende Blöcke der Länge \(m\) und zählen für jeden Typ die fünf Restklassen.
Anschließend werden zwei Wahrscheinlichkeitsvektoren geführt. Während eines Durchlaufs über alle Typen wird jede neue Koordinate nur aus dem Vektor des vorherigen Durchlaufs berechnet; das Update geschieht also synchron und nicht direkt im selben Array.
Nach jedem Durchlauf wird die größte absolute Koordinatenänderung gemessen. Die C++-Implementierung benutzt erweiterte Präzision und stoppt unter \(10^{-20}\); die Python- und Java-Implementierungen arbeiten mit doppelter Präzision und stoppen unter \(10^{-15}\). Zusätzlich gilt in allen drei Fällen eine harte Obergrenze von \(10^6\) Iterationen.
Nach der Konvergenz wird die Aussterbenswahrscheinlichkeit für Typ \(0\) mit acht Nachkommastellen ausgegeben.
Das Erzeugen der pseudozufälligen Folge und das Aufbauen der Fünfer-Histogramme kostet \(O(km)\) Zeit. Benötigt der Fixpunktlöser \(I\) Durchläufe, so kostet die nichtlineare Iteration \(O(Ik)\) Zeit. Insgesamt ergibt sich damit
$$O(km+Ik).$$
So wie der Code geschrieben ist, speichert er die gesamte Folge der Länge \(km\), die \(5k\) Histogrammwerte und zwei Wahrscheinlichkeitsvektoren der Länge \(k\). Der Speicherverbrauch liegt daher bei \(O(km)\).
\(k\) farklı bakteri tipi vardır. Her tip \(i\) için tek bir çoğalma kuralı değil, deterministik pseudo-random diziden türetilen \(m\) tane eş olasılıklı kural bulunur:
$$r_0=306,\qquad r_{t+1}\equiv r_t^2 \pmod{10007}.$$
\(r_{im},r_{im+1},\dots,r_{im+m-1}\) bloğu tip \(i\)'ye aittir. Süreç, tip \(0\) olan tek bir bakteriyle başlar; kurallar tüm alt soylar üzerinde bağımsız biçimde tekrar tekrar uygulanır. İstenen şey, popülasyonun sonunda tamamen yok olma olasılığıdır.
Bu yapı çok tipli bir dallanma sürecidir. Ana sadeleştirme, her tipin \(m\) kuralını beş artık sınıfı sayısına indirgemektir. Böylece problem \([0,1]^k\) üzerinde bir sabit nokta sistemine dönüşür.
Her tip \(i\) ve her \(q\in\{0,1,2,3,4\}\) için
$$c_{i,q}=\#\left\{j\in\{0,\dots,m-1\}: r_{im+j}\equiv q \pmod{5}\right\}$$
tanımını yapalım. Bu sayılar
$$\sum_{q=0}^4 c_{i,q}=m$$
eşitliğini sağlar. Bu aşamadan sonra tek tek kural değerlerine ihtiyaç kalmaz; her tip yalnızca beş olası davranışın kaç kez görüldüğü ile temsil edilir.
Beş artık sınıfı, tip \(i\) olan tek bir bakterinin beş olası sonucunu temsil eder.
Artık \(0\): hiç yavru oluşmaz.
Artık \(1\): tip \(i\) olan iki yavru oluşur.
Artık \(2\): tip \(2i \bmod k\) olan bir yavru oluşur.
Artık \(3\): tip \((i^2+1) \bmod k\) olan üç yavru oluşur.
Artık \(4\): tip \(i\) olan bir yavru ve tip \((i+1) \bmod k\) olan bir yavru oluşur.
\(s_j\) bir tip \(j\) yavrusunun yok olma katkısını göstersin. O zaman bu beş durum şu monomlarla temsil edilir:
$$1,\qquad s_i^2,\qquad s_{2i \bmod k},\qquad s_{(i^2+1) \bmod k}^3,\qquad s_i s_{(i+1) \bmod k}.$$
Kare, küp ve çarpım terimleri bağımsızlıktan gelir: seçilen kuralın ürettiği tüm alt soyların birlikte yok olması gerekir.
\(\mathbf{s}=(s_0,\dots,s_{k-1})\) olsun. Tip \(i\) için \(m\) eş olasılıklı kuralın ortalaması, yavru üreten polinomu verir:
$$f_i(\mathbf{s})=\frac{c_{i,0}+c_{i,1}s_i^2+c_{i,2}s_{2i \bmod k}+c_{i,3}s_{(i^2+1) \bmod k}^3+c_{i,4}s_i s_{(i+1) \bmod k}}{m}.$$
Şimdi \(p_i\), tek bir tip \(i\) bakterisi ile başlayan sürecin sonunda yok olma olasılığı olsun. Dallanma süreci teorisine göre
$$p_i=f_i(\mathbf{p})=\frac{c_{i,0}+c_{i,1}p_i^2+c_{i,2}p_{2i \bmod k}+c_{i,3}p_{(i^2+1) \bmod k}^3+c_{i,4}p_i p_{(i+1) \bmod k}}{m}.$$
Dolayısıyla aranan cevap, yok olma vektörünün ilk bileşenidir.
Şu dönüşümü tanımlayalım:
$$T(\mathbf{x})=\left(f_0(\mathbf{x}),f_1(\mathbf{x}),\dots,f_{k-1}(\mathbf{x})\right).$$
Tüm katsayılar negatif olmadığı ve her satır \(m\)'ye bölünmüş bir ortalama olduğu için bu dönüşüm \([0,1]^k\) kümesini kendi içine gönderir ve koordinat bazında monotondur. Bu yüzden
$$\mathbf{p}^{(0)}=(0,\dots,0),\qquad \mathbf{p}^{(t+1)}=T\left(\mathbf{p}^{(t)}\right)$$
ile başladığımızda artan ve \((1,\dots,1)\) tarafından yukarıdan sınırlanan bir dizi elde ederiz. Limit, sistemin en küçük sabit noktasıdır ve tam olarak fiziksel yok olma olasılıklarını verir.
İlk dört pseudo-random değer
$$306,\qquad 3573,\qquad 7404,\qquad 870$$
olduğundan, bunların \(5\)'e göre kalanları
$$1,\qquad 3,\qquad 4,\qquad 0$$
şeklindedir. Böylece tip \(0\), bir adet artık-\(1\) ve bir adet artık-\(3\) kuralı alır; tip \(1\) ise bir adet artık-\(4\) ve bir adet artık-\(0\) kuralı alır. Denklem sistemi
$$p_0=\frac{p_0^2+p_1^3}{2},\qquad p_1=\frac{1+p_0p_1}{2}$$
haline gelir. \((0,0)\) noktasından başlayan ilk senkron iterasyonlar
$$\begin{aligned} \mathbf{p}^{(1)}&=\left(0,\frac12\right),\\ \mathbf{p}^{(2)}&=\left(\frac1{16},\frac12\right),\\ \mathbf{p}^{(3)}&=\left(\frac{33}{512},\frac{33}{64}\right) \end{aligned}$$
olur. Süreç devam ettiğinde
$$p_0\approx 0.07243802$$
değerine yakınsar; bu da implementasyondaki küçük doğrulama örneğiyle aynıdır.
C++, Python ve Java implementasyonları aynı akışı izler. Önce uzunluğu \(km\) olan pseudo-random dizi tamamen üretilir. Sonra bu dizi, uzunluğu \(m\) olan \(k\) bloğa ayrılır ve her tip için beş artık sınıfının sayıları hesaplanır.
Ardından iki olasılık vektörü tutulur. Tiplere yapılan bir tarama sırasında her yeni koordinat yalnızca bir önceki iterasyonun vektöründen hesaplanır; yani güncelleme yerinde değil, senkron biçimde yapılır. Tur bitince iki vektör yer değiştirir.
Her turdan sonra eski ve yeni değerler arasındaki en büyük mutlak fark ölçülür. C++ implementasyonu genişletilmiş hassasiyet kullanır ve bu fark \(10^{-20}\)'nin altına indiğinde durur; Python ve Java implementasyonları çift hassasiyet kullandığı için eşik \(10^{-15}\)'tir. Üçünde de ayrıca \(10^6\) iterasyonluk sert bir üst sınır vardır.
Yakınsama sağlandıktan sonra tip \(0\)'ın yok olma olasılığı sekiz ondalık basamakla yazdırılır.
Pseudo-random diziyi üretmek ve her tip için beşli histogramı kurmak \(O(km)\) zaman alır. Sabit nokta çözücüsü \(I\) tur gerektirirse iterasyon kısmı \(O(Ik)\) olur. Toplam süre bu nedenle
$$O(km+Ik)$$
şeklindedir. Kodun mevcut hali, uzunluğu \(km\) olan dizinin tamamını, \(5k\) sayımları ve uzunluğu \(k\) olan iki olasılık vektörünü sakladığı için bellek kullanımı \(O(km)\)'dir.
Hay \(k\) tipos de bacterias. Cada tipo \(i\) no tiene una sola regla de reproducción, sino \(m\) reglas equiprobables extraídas de la sucesión pseudoaleatoria determinista
$$r_0=306,\qquad r_{t+1}\equiv r_t^2 \pmod{10007}.$$
El bloque \(r_{im},r_{im+1},\dots,r_{im+m-1}\) pertenece al tipo \(i\). Partiendo de una sola bacteria del tipo \(0\), las reglas se aplican de manera independiente a toda la descendencia, y debemos calcular la probabilidad de que la población completa termine extinguiéndose.
Esto es un proceso de ramificación multitype. La simplificación clave consiste en resumir las \(m\) reglas de cada tipo mediante cinco conteos de residuos; después de eso, el problema se convierte en un sistema de punto fijo sobre \([0,1]^k\).
Para cada tipo \(i\) y cada residuo \(q\in\{0,1,2,3,4\}\), definimos
$$c_{i,q}=\#\left\{j\in\{0,\dots,m-1\}: r_{im+j}\equiv q \pmod{5}\right\}.$$
Estos conteos satisfacen
$$\sum_{q=0}^4 c_{i,q}=m.$$
Una vez conocidos, ya no importa qué valores concretos produjo la secuencia. Cada bacteria del tipo \(i\) solo necesita saber cuántas veces aparece cada una de las cinco clases de comportamiento.
Las cinco clases residuales codifican cinco posibles resultados para una bacteria del tipo \(i\).
Residuo \(0\): no se produce descendencia.
Residuo \(1\): se producen dos descendientes del tipo \(i\).
Residuo \(2\): se produce un descendiente del tipo \(2i \bmod k\).
Residuo \(3\): se producen tres descendientes del tipo \((i^2+1) \bmod k\).
Residuo \(4\): se producen un descendiente del tipo \(i\) y uno del tipo \((i+1) \bmod k\).
Si \(s_j\) representa la contribución a la extinción de un solo hijo del tipo \(j\), estos cinco casos aportan los monomios
$$1,\qquad s_i^2,\qquad s_{2i \bmod k},\qquad s_{(i^2+1) \bmod k}^3,\qquad s_i s_{(i+1) \bmod k}.$$
Los cuadrados, cubos y productos aparecen porque las líneas de descendencia son independientes: todos los descendientes creados por la regla elegida deben extinguirse.
Sea \(\mathbf{s}=(s_0,\dots,s_{k-1})\). Al promediar las \(m\) reglas equiprobables del tipo \(i\), obtenemos el polinomio generador
$$f_i(\mathbf{s})=\frac{c_{i,0}+c_{i,1}s_i^2+c_{i,2}s_{2i \bmod k}+c_{i,3}s_{(i^2+1) \bmod k}^3+c_{i,4}s_i s_{(i+1) \bmod k}}{m}.$$
Si \(p_i\) es la probabilidad de extinción al comenzar con una sola bacteria del tipo \(i\), la teoría estándar de procesos de ramificación da
$$p_i=f_i(\mathbf{p})=\frac{c_{i,0}+c_{i,1}p_i^2+c_{i,2}p_{2i \bmod k}+c_{i,3}p_{(i^2+1) \bmod k}^3+c_{i,4}p_i p_{(i+1) \bmod k}}{m}.$$
La respuesta pedida es, por tanto, la primera coordenada del vector de extinción \(\mathbf{p}\).
Definimos la aplicación
$$T(\mathbf{x})=\left(f_0(\mathbf{x}),f_1(\mathbf{x}),\dots,f_{k-1}(\mathbf{x})\right).$$
Como todos los coeficientes son no negativos y cada fila es un promedio dividido por \(m\), la aplicación envía \([0,1]^k\) sobre sí mismo y es monótona coordenada a coordenada. Si comenzamos con
$$\mathbf{p}^{(0)}=(0,\dots,0),\qquad \mathbf{p}^{(t+1)}=T\left(\mathbf{p}^{(t)}\right),$$
obtenemos una sucesión creciente acotada superiormente por \((1,\dots,1)\). Su límite es el punto fijo mínimo del sistema, que coincide con el vector físicamente correcto de probabilidades de extinción.
Los primeros cuatro valores pseudoaleatorios son
$$306,\qquad 3573,\qquad 7404,\qquad 870,$$
y sus residuos módulo \(5\) son
$$1,\qquad 3,\qquad 4,\qquad 0.$$
Así, el tipo \(0\) recibe una regla de residuo \(1\) y una de residuo \(3\), mientras que el tipo \(1\) recibe una de residuo \(4\) y una de residuo \(0\). El sistema queda
$$p_0=\frac{p_0^2+p_1^3}{2},\qquad p_1=\frac{1+p_0p_1}{2}.$$
Partiendo de \((0,0)\), las primeras iteraciones síncronas son
$$\begin{aligned} \mathbf{p}^{(1)}&=\left(0,\frac12\right),\\ \mathbf{p}^{(2)}&=\left(\frac1{16},\frac12\right),\\ \mathbf{p}^{(3)}&=\left(\frac{33}{512},\frac{33}{64}\right). \end{aligned}$$
Si continuamos, obtenemos
$$p_0\approx 0.07243802,$$
en acuerdo con el pequeño caso de comprobación usado por la implementación.
Las implementaciones en C++, Python y Java siguen la misma tubería. Primero generan la sucesión pseudoaleatoria completa de longitud \(km\). Después la dividen en \(k\) bloques consecutivos de tamaño \(m\) y cuentan, para cada tipo, cuántas veces aparece cada una de las cinco clases residuales.
Luego mantienen dos vectores de probabilidades. Durante un barrido sobre todos los tipos, cada valor nuevo se calcula únicamente a partir del vector de la iteración anterior, de modo que la actualización es síncrona y no en el mismo lugar.
Tras cada barrido se registra la mayor diferencia absoluta entre coordenadas viejas y nuevas. La implementación en C++ usa precisión extendida y se detiene por debajo de \(10^{-20}\); las de Python y Java usan precisión doble y se detienen por debajo de \(10^{-15}\). Las tres fijan además un máximo estricto de \(10^6\) iteraciones.
Una vez alcanzada la convergencia, se imprime con ocho decimales la probabilidad de extinción del tipo \(0\).
Generar la sucesión pseudoaleatoria y construir el histograma de cinco conteos para cada tipo cuesta \(O(km)\) tiempo. Si el solucionador de punto fijo necesita \(I\) barridos, la parte iterativa cuesta \(O(Ik)\). En total, el tiempo de ejecución es
$$O(km+Ik).$$
Tal como están escritas, las implementaciones almacenan la sucesión completa de longitud \(km\), los \(5k\) conteos del histograma y dos vectores de probabilidades de longitud \(k\). Por ello, la memoria utilizada es \(O(km)\).
题目中有 \(k\) 种细菌。每一种类型 \(i\) 并不是只有一条繁殖规则，而是拥有 \(m\) 条等概率规则，这些规则来自确定性的伪随机序列
$$r_0=306,\qquad r_{t+1}\equiv r_t^2 \pmod{10007}.$$
第 \(i\) 类细菌对应的规则块是 \(r_{im},r_{im+1},\dots,r_{im+m-1}\)。从一只类型 \(0\) 的细菌开始，后代会继续独立地按照这些规则繁殖。我们要求的是整个族群最终完全灭绝的概率。
这是一个典型的多类型分枝过程。核心简化是先把每个类型的 \(m\) 条规则压缩成五个计数，再把灭绝概率写成 \([0,1]^k\) 上的固定点方程组。
对每个类型 \(i\) 以及每个余数 \(q\in\{0,1,2,3,4\}\)，定义
$$c_{i,q}=\#\left\{j\in\{0,\dots,m-1\}: r_{im+j}\equiv q \pmod{5}\right\}.$$
这些计数满足
$$\sum_{q=0}^4 c_{i,q}=m.$$
一旦知道了这五个数，原始的 \(m\) 条规则就不必逐条保存或逐条分析了。对某个固定类型来说，真正重要的只是五种行为各出现了多少次。
五个余数类分别对应一只类型 \(i\) 细菌的五种可能结果。
余数 \(0\)：不产生任何后代。
余数 \(1\)：产生两只类型 \(i\) 的后代。
余数 \(2\)：产生一只类型 \(2i \bmod k\) 的后代。
余数 \(3\)：产生三只类型 \((i^2+1) \bmod k\) 的后代。
余数 \(4\)：产生一只类型 \(i\) 的后代和一只类型 \((i+1) \bmod k\) 的后代。
如果 \(s_j\) 表示“一只类型 \(j\) 后代最终灭绝”的贡献，那么这五种情况对应的单项式就是
$$1,\qquad s_i^2,\qquad s_{2i \bmod k},\qquad s_{(i^2+1) \bmod k}^3,\qquad s_i s_{(i+1) \bmod k}.$$
这里的平方、立方和乘积都来自独立性：选中某条规则后，它产生的每一条后代谱系都必须灭绝，总概率才是各谱系灭绝概率的乘积。
记 \(\mathbf{s}=(s_0,\dots,s_{k-1})\)。对类型 \(i\) 的 \(m\) 条等概率规则取平均，可得到后代生成多项式
$$f_i(\mathbf{s})=\frac{c_{i,0}+c_{i,1}s_i^2+c_{i,2}s_{2i \bmod k}+c_{i,3}s_{(i^2+1) \bmod k}^3+c_{i,4}s_i s_{(i+1) \bmod k}}{m}.$$
再令 \(p_i\) 表示“从一只类型 \(i\) 的细菌开始，整个过程最终灭绝”的概率。根据分枝过程的标准理论，灭绝向量 \(\mathbf{p}\) 满足
$$p_i=f_i(\mathbf{p})=\frac{c_{i,0}+c_{i,1}p_i^2+c_{i,2}p_{2i \bmod k}+c_{i,3}p_{(i^2+1) \bmod k}^3+c_{i,4}p_i p_{(i+1) \bmod k}}{m}.$$
因此题目要求的答案，就是这个向量的第一项 \(p_0\)。
定义映射
$$T(\mathbf{x})=\left(f_0(\mathbf{x}),f_1(\mathbf{x}),\dots,f_{k-1}(\mathbf{x})\right).$$
由于所有系数都非负，而且每一行都是除以 \(m\) 的平均值，所以 \(T\) 会把 \([0,1]^k\) 映到自身，并且对每个坐标都是单调的。于是从
$$\mathbf{p}^{(0)}=(0,\dots,0),\qquad \mathbf{p}^{(t+1)}=T\left(\mathbf{p}^{(t)}\right)$$
开始迭代，就得到一个逐坐标不减且被 \((1,\dots,1)\) 上界约束的序列。它一定收敛，而且极限正是这个系统的最小固定点，也就是实际有意义的灭绝概率向量。
前四个伪随机值是
$$306,\qquad 3573,\qquad 7404,\qquad 870,$$
模 \(5\) 的余数依次为
$$1,\qquad 3,\qquad 4,\qquad 0.$$
因此类型 \(0\) 恰好得到一条余数为 \(1\) 的规则和一条余数为 \(3\) 的规则；类型 \(1\) 恰好得到一条余数为 \(4\) 的规则和一条余数为 \(0\) 的规则。于是方程组变成
$$p_0=\frac{p_0^2+p_1^3}{2},\qquad p_1=\frac{1+p_0p_1}{2}.$$
从 \((0,0)\) 出发，前几次同步迭代为
$$\begin{aligned} \mathbf{p}^{(1)}&=\left(0,\frac12\right),\\ \mathbf{p}^{(2)}&=\left(\frac1{16},\frac12\right),\\ \mathbf{p}^{(3)}&=\left(\frac{33}{512},\frac{33}{64}\right). \end{aligned}$$
继续迭代可得
$$p_0\approx 0.07243802,$$
这与实现中使用的小规模校验值一致。
C++、Python 和 Java 三个实现采用的流程完全一致。第一步是生成长度为 \(km\) 的完整伪随机序列。第二步是把它切成 \(k\) 个长度为 \(m\) 的连续块，并为每个类型统计五种余数出现的次数。
接着，实现维护两个概率向量。一次扫描所有类型时，每个新值都只依赖上一轮的整个位向量，因此更新是同步的，而不是就地覆盖。扫描结束后，再交换“旧向量”和“新向量”。
每一轮还会记录新旧两次迭代之间的最大绝对坐标差。C++ 实现使用扩展精度，并在该差值低于 \(10^{-20}\) 时停止；Python 与 Java 使用双精度，因此阈值是 \(10^{-15}\)。三者都还设置了 \(10^6\) 次迭代的硬上限。
收敛之后，程序输出类型 \(0\) 的灭绝概率，并保留八位小数。
生成伪随机序列并为每个类型构造五项计数直方图需要 \(O(km)\) 时间。如果固定点迭代一共进行了 \(I\) 轮，那么求解阶段需要 \(O(Ik)\) 时间，所以总时间复杂度为
$$O(km+Ik).$$
按照当前实现方式，程序会显式保存长度为 \(km\) 的完整序列、\(5k\) 个计数以及两个长度为 \(k\) 的概率向量，因此空间复杂度是 \(O(km)\)。
Имеется \(k\) типов бактерий. У типа \(i\) не одна, а \(m\) равновероятных правил размножения, полученных из детерминированной псевдослучайной последовательности
$$r_0=306,\qquad r_{t+1}\equiv r_t^2 \pmod{10007}.$$
Блок \(r_{im},r_{im+1},\dots,r_{im+m-1}\) относится к типу \(i\). Процесс стартует с одной бактерии типа \(0\); затем правила независимо применяются ко всем потомкам. Нужно найти вероятность того, что вся популяция в конечном счете вымрет.
Это многотипный ветвящийся процесс. Главная идея состоит в том, чтобы заменить \(m\) правил каждого типа всего пятью счетчиками остатков. После этого задача сводится к системе уравнений на неподвижную точку в \([0,1]^k\).
Для каждого типа \(i\) и каждого остатка \(q\in\{0,1,2,3,4\}\) положим
$$c_{i,q}=\#\left\{j\in\{0,\dots,m-1\}: r_{im+j}\equiv q \pmod{5}\right\}.$$
Эти величины удовлетворяют равенству
$$\sum_{q=0}^4 c_{i,q}=m.$$
После вычисления пяти счетчиков отдельные значения последовательности больше не нужны. Для типа \(i\) важно только то, сколько раз встречается каждый из пяти вариантов поведения.
Пять классов по модулю \(5\) соответствуют пяти возможным результатам для одной бактерии типа \(i\).
Остаток \(0\): потомков нет.
Остаток \(1\): рождаются две бактерии типа \(i\).
Остаток \(2\): рождается одна бактерия типа \(2i \bmod k\).
Остаток \(3\): рождаются три бактерии типа \((i^2+1) \bmod k\).
Остаток \(4\): рождаются одна бактерия типа \(i\) и одна бактерия типа \((i+1) \bmod k\).
Если \(s_j\) обозначает вклад в вымирание от одного потомка типа \(j\), то этим случаям соответствуют мономы
$$1,\qquad s_i^2,\qquad s_{2i \bmod k},\qquad s_{(i^2+1) \bmod k}^3,\qquad s_i s_{(i+1) \bmod k}.$$
Квадраты, кубы и произведения появляются из независимости родословных: все ветви, рожденные выбранным правилом, должны вымереть одновременно.
Пусть \(\mathbf{s}=(s_0,\dots,s_{k-1})\). Усреднение по \(m\) равновероятным правилам типа \(i\) дает производящий полином
$$f_i(\mathbf{s})=\frac{c_{i,0}+c_{i,1}s_i^2+c_{i,2}s_{2i \bmod k}+c_{i,3}s_{(i^2+1) \bmod k}^3+c_{i,4}s_i s_{(i+1) \bmod k}}{m}.$$
Пусть теперь \(p_i\) означает вероятность вымирания процесса, стартующего с одной бактерии типа \(i\). Тогда стандартная теория ветвящихся процессов дает
$$p_i=f_i(\mathbf{p})=\frac{c_{i,0}+c_{i,1}p_i^2+c_{i,2}p_{2i \bmod k}+c_{i,3}p_{(i^2+1) \bmod k}^3+c_{i,4}p_i p_{(i+1) \bmod k}}{m}.$$
Искомый ответ равен первой координате вектора вымирания \(\mathbf{p}\).
Определим отображение
$$T(\mathbf{x})=\left(f_0(\mathbf{x}),f_1(\mathbf{x}),\dots,f_{k-1}(\mathbf{x})\right).$$
Поскольку все коэффициенты неотрицательны и каждая строка представляет собой среднее, деленное на \(m\), отображение переводит \([0,1]^k\) в себя и монотонно по каждой координате. Поэтому, начиная с
$$\mathbf{p}^{(0)}=(0,\dots,0),\qquad \mathbf{p}^{(t+1)}=T\left(\mathbf{p}^{(t)}\right),$$
мы получаем возрастающую последовательность, ограниченную сверху вектором \((1,\dots,1)\). Она сходится к наименьшей неподвижной точке системы, а это и есть физически правильный вектор вероятностей вымирания.
Первые четыре псевдослучайных значения равны
$$306,\qquad 3573,\qquad 7404,\qquad 870,$$
а их остатки по модулю \(5\) таковы:
$$1,\qquad 3,\qquad 4,\qquad 0.$$
Следовательно, тип \(0\) получает одно правило с остатком \(1\) и одно с остатком \(3\), а тип \(1\) получает одно правило с остатком \(4\) и одно с остатком \(0\). Система принимает вид
$$p_0=\frac{p_0^2+p_1^3}{2},\qquad p_1=\frac{1+p_0p_1}{2}.$$
Если стартовать из \((0,0)\), то первые синхронные итерации равны
$$\begin{aligned} \mathbf{p}^{(1)}&=\left(0,\frac12\right),\\ \mathbf{p}^{(2)}&=\left(\frac1{16},\frac12\right),\\ \mathbf{p}^{(3)}&=\left(\frac{33}{512},\frac{33}{64}\right). \end{aligned}$$
Дальнейшая итерация приводит к значению
$$p_0\approx 0.07243802,$$
что совпадает с малым проверочным случаем, используемым в реализации.
Реализации на C++, Python и Java устроены одинаково. Сначала генерируется вся псевдослучайная последовательность длины \(km\). Затем она разбивается на \(k\) последовательных блоков по \(m\) чисел, и для каждого типа подсчитываются частоты пяти остатков.
После этого поддерживаются два вектора вероятностей. За один проход по всем типам каждая новая координата вычисляется только из предыдущего вектора, так что обновление выполняется синхронно, без перезаписи значений на месте. После прохода векторы меняются ролями.
На каждом проходе также измеряется максимальное абсолютное изменение координаты. Реализация на C++ использует расширенную точность и останавливается при пороге \(10^{-20}\); реализации на Python и Java используют двойную точность и останавливаются при пороге \(10^{-15}\). Во всех трех версиях дополнительно есть жесткий предел в \(10^6\) итераций.
После сходимости печатается вероятность вымирания для типа \(0\) с восемью знаками после запятой.
Построение псевдослучайной последовательности и подсчет пяти счетчиков для каждого типа занимает \(O(km)\) времени. Если итерационный решатель делает \(I\) проходов, то на решение системы уходит \(O(Ik)\) времени. Итого полное время работы равно
$$O(km+Ik).$$
В текущем виде реализации хранят всю последовательность длины \(km\), \(5k\) гистограммных счетчиков и два вектора вероятностей длины \(k\). Поэтому объем памяти равен \(O(km)\).
لدينا \(k\) أنواع من البكتيريا. كل نوع \(i\) لا يملك قاعدة تكاثر واحدة فقط، بل يملك \(m\) قواعد متساوية الاحتمال مشتقة من المتتالية شبه العشوائية الحتمية
$$r_0=306,\qquad r_{t+1}\equiv r_t^2 \pmod{10007}.$$
الكتلة \(r_{im},r_{im+1},\dots,r_{im+m-1}\) تخص النوع \(i\). نبدأ ببكتيريا واحدة من النوع \(0\)، ثم تُطبَّق القواعد بصورة مستقلة على كل السلالات الناتجة، والمطلوب هو احتمال أن ينقرض المجتمع كله في النهاية.
هذه مسألة عملية تفرع متعددة الأنواع. الفكرة الأساسية هي ضغط قواعد كل نوع وعددها \(m\) إلى خمسة عدادات فقط، وبعد ذلك يصبح المطلوب حل نظام نقطة ثابتة على \([0,1]^k\).
لكل نوع \(i\) ولكل باقي \(q\in\{0,1,2,3,4\}\) نعرّف
$$c_{i,q}=\#\left\{j\in\{0,\dots,m-1\}: r_{im+j}\equiv q \pmod{5}\right\}.$$
وهذه القيم تحقق
$$\sum_{q=0}^4 c_{i,q}=m.$$
بعد معرفة هذه الأعداد الخمسة لا نعود بحاجة إلى القيم المفردة للقواعد. ما يهم فعلاً هو عدد مرات ظهور كل سلوك من السلوكات الخمسة الممكنة.
الفئات الخمس لبواقي القسمة على \(5\) تمثل خمس نتائج ممكنة لبكتيريا واحدة من النوع \(i\).
الباقي \(0\): لا يُنتج أي نسل.
الباقي \(1\): يُنتج نسلان من النوع \(i\).
الباقي \(2\): يُنتج نسل واحد من النوع \(2i \bmod k\).
الباقي \(3\): يُنتج ثلاثة أفراد من النوع \((i^2+1) \bmod k\).
الباقي \(4\): يُنتج فرداً من النوع \(i\) وفرداً من النوع \((i+1) \bmod k\).
إذا رمزنا بـ \(s_j\) إلى مساهمة فرد واحد من النوع \(j\) في حدث الانقراض، فإن هذه الحالات الخمس تعطي الحدود
$$1,\qquad s_i^2,\qquad s_{2i \bmod k},\qquad s_{(i^2+1) \bmod k}^3,\qquad s_i s_{(i+1) \bmod k}.$$
الأسس والجداءات هنا ناتجة من الاستقلال: يجب أن تنقرض كل السلالات الناتجة عن القاعدة المختارة لكي يتحقق الانقراض الكلي.
لتكن \(\mathbf{s}=(s_0,\dots,s_{k-1})\). بأخذ متوسط القواعد \(m\) المتساوية الاحتمال للنوع \(i\) نحصل على كثير الحدود المولد
$$f_i(\mathbf{s})=\frac{c_{i,0}+c_{i,1}s_i^2+c_{i,2}s_{2i \bmod k}+c_{i,3}s_{(i^2+1) \bmod k}^3+c_{i,4}s_i s_{(i+1) \bmod k}}{m}.$$
ولنرمز بـ \(p_i\) إلى احتمال انقراض العملية إذا بدأت ببكتيريا واحدة من النوع \(i\). عندئذ تعطي نظرية عمليات التفرع العلاقة
$$p_i=f_i(\mathbf{p})=\frac{c_{i,0}+c_{i,1}p_i^2+c_{i,2}p_{2i \bmod k}+c_{i,3}p_{(i^2+1) \bmod k}^3+c_{i,4}p_i p_{(i+1) \bmod k}}{m}.$$
إذن الجواب المطلوب هو الإحداثي الأول من متجه الانقراض \(\mathbf{p}\).
نعرّف التحويل
$$T(\mathbf{x})=\left(f_0(\mathbf{x}),f_1(\mathbf{x}),\dots,f_{k-1}(\mathbf{x})\right).$$
بما أن جميع المعاملات غير سالبة، وكل صف هو متوسط مقسوم على \(m\)، فإن هذا التحويل يرسل \([0,1]^k\) إلى نفسه وهو أحادي الرتابة في كل إحداثي. لذلك إذا بدأنا من
$$\mathbf{p}^{(0)}=(0,\dots,0),\qquad \mathbf{p}^{(t+1)}=T\left(\mathbf{p}^{(t)}\right),$$
فإننا نحصل على متتالية متزايدة ومحدودة من الأعلى بالمتجه \((1,\dots,1)\). هذه المتتالية تتقارب إلى أصغر نقطة ثابتة للنظام، وهي بالضبط متجه احتمالات الانقراض الصحيح.
أول أربعة حدود من المتتالية شبه العشوائية هي
$$306,\qquad 3573,\qquad 7404,\qquad 870,$$
ومن ثم تكون البواقي modulo \(5\) هي
$$1,\qquad 3,\qquad 4,\qquad 0.$$
إذن النوع \(0\) يحصل على قاعدة واحدة من فئة الباقي \(1\) وقاعدة واحدة من فئة الباقي \(3\)، بينما النوع \(1\) يحصل على قاعدة واحدة من فئة الباقي \(4\) وقاعدة واحدة من فئة الباقي \(0\). لذلك يصبح النظام
$$p_0=\frac{p_0^2+p_1^3}{2},\qquad p_1=\frac{1+p_0p_1}{2}.$$
وبالانطلاق من \((0,0)\)، نحصل في الخطوات المتزامنة الأولى على
$$\begin{aligned} \mathbf{p}^{(1)}&=\left(0,\frac12\right),\\ \mathbf{p}^{(2)}&=\left(\frac1{16},\frac12\right),\\ \mathbf{p}^{(3)}&=\left(\frac{33}{512},\frac{33}{64}\right). \end{aligned}$$
وبمتابعة التكرار نصل إلى
$$p_0\approx 0.07243802,$$
وهو نفس المقدار الصغير الذي تستعمله الشيفرة كحالة تحقق.
تتبع تطبيقات C++ وPython وJava البنية نفسها. أولاً تُولَّد المتتالية شبه العشوائية كاملة بطول \(km\). ثم تُقسَّم إلى \(k\) كتل متتالية طول كل منها \(m\)، وتُحصى لكل نوع أعداد ظهور البواقي الخمسة.
بعد ذلك تحتفظ الشيفرة بمتجهين للاحتمالات. أثناء المرور على جميع الأنواع، يُحسب كل مقدار جديد اعتماداً فقط على متجه الجولة السابقة، ولذلك يكون التحديث متزامناً لا موضعياً. وبعد انتهاء الجولة يتم تبادل المتجهين.
في كل جولة يُسجَّل أيضاً أكبر فرق مطلق بين الإحداثيات القديمة والجديدة. تنفيذ C++ يستخدم دقة موسعة ويتوقف عندما ينخفض هذا الفرق تحت \(10^{-20}\)، أما تنفيذا Python وJava فيستخدمان الدقة المزدوجة ويتوقفان تحت \(10^{-15}\). وفي الحالات الثلاث هناك حد أقصى صارم مقداره \(10^6\) تكرار.
بعد الوصول إلى التقارب، تُطبع قيمة احتمال انقراض النوع \(0\) إلى ثمانية منازل عشرية.
توليد المتتالية شبه العشوائية وبناء عدادات الأنواع الخمسة لكل نوع يحتاج إلى زمن \(O(km)\). وإذا احتاج حل النقطة الثابتة إلى \(I\) جولات، فإن الجزء التكراري يحتاج إلى \(O(Ik)\). لذا يكون الزمن الكلي
$$O(km+Ik).$$
كما هي مكتوبة الآن، تقوم التطبيقات بتخزين المتتالية الكاملة بطول \(km\)، إضافة إلى \(5k\) عداداً ومتجهين للاحتمالات بطول \(k\). لذلك يكون استهلاك الذاكرة \(O(km)\).