Each round costs \(t\) units. The payout is \(2^n\) with probability \(2^{-(n+1)}\) for \(n \ge 0\). Starting from bankroll \(B\), we want the probability of being able to keep playing forever, which is the same as never letting the bankroll fall below the amount needed for the next ticket. The bankroll can in principle grow without bound, so the state space is infinite. The implementation therefore replaces the raw infinite recursion by a finite residue model together with a self-similar matrix equation.
The numerical method follows the same probability law in all three implementations. Its key idea is that after shifting the bankroll by \(t-1\), the process becomes translation-invariant in blocks of size \(t-1\).
Let
$$d=t-1,\qquad c=B-d.$$
Here \(c\) is the bankroll measured above the minimum playable threshold. One round changes \(c\) by
$$c' = c + 2^n - t.$$
If \(c' \le 0\), the player can no longer buy another ticket, so that path is ruined.
The exact distribution has infinitely many outcomes \(n=0,1,2,\dots\), but the implementation keeps only
$$0 \le n \lt H,\qquad H=50.$$
The omitted probability mass is
$$\sum_{n=H}^{\infty} 2^{-(n+1)} = 2^{-50},$$
which is far below the requested \(7\)-decimal accuracy.
Every playable shifted bankroll \(c \ge 1\) can be written uniquely as
$$c = k d + r + 1,\qquad k \ge 0,\qquad 0 \le r \lt d.$$
The integer \(k\) is the layer, and \(r\) is the residue inside that layer. This decomposition is useful because
$$k d + r + 1 + 2^n - t = (k-1)d + (r + 2^n).$$
So increasing the starting layer by one simply shifts all future layers by one. The detailed geometry inside each layer stays the same.
Define
$$u_k(r)=\Pr(\text{eventual ruin} \mid c = k d + r + 1).$$
For fixed \(t\), each \(u_k\) is a vector of length \(d\).
Because the process repeats the same pattern from one layer to the next, there is a single \(d \times d\) matrix \(X\) such that
$$u_{k+1} = X u_k,\qquad u_k = X^k u_0.$$
To derive \(X\), start from layer \(1\), where \(c = d + r + 1\). After outcome \(n\),
$$c' = d + r + 1 + 2^n - t = r + 2^n.$$
Write this new state as
$$c' = j_n(r)\,d + \rho_n(r) + 1,$$
with
$$j_n(r)=\left\lfloor\frac{r+2^n-1}{d}\right\rfloor,\qquad \rho_n(r)=(r+2^n-1)\bmod d.$$
If \(j_n(r)=0\), the process lands in the boundary layer. If \(j_n(r)\ge 1\), it lands \(j_n(r)\) layers above the boundary. Therefore the \(r\)-th row of \(X\) satisfies
$$X_{r,*}=\sum_{j_n(r)=0} 2^{-(n+1)} e_{\rho_n(r)}^{\mathsf T}+\sum_{j_n(r)\ge 1} 2^{-(n+1)}\bigl(X^{j_n(r)}\bigr)_{\rho_n(r),*},$$
where \(e_{\rho}^{\mathsf T}\) is the row vector with a \(1\) in column \(\rho\) and zeros elsewhere.
This is a matrix fixed-point equation. The implementations iterate it numerically until successive matrices differ by less than \(10^{-18}\), or until the safety iteration cap is reached.
Now consider the lowest playable layer, where
$$c = r + 1,\qquad 0 \le r \lt d.$$
After outcome \(n\),
$$c' = r + 1 + 2^n - t = r + 2^n - d.$$
If \(c' \le 0\), ruin happens immediately. Otherwise write
$$c' = j'_n(r)\,d + \rho'_n(r) + 1,$$
where
$$j'_n(r)=\left\lfloor\frac{c'-1}{d}\right\rfloor,\qquad \rho'_n(r)=(c'-1)\bmod d.$$
Conditioning on the first round gives
$$u_0(r)=q_r+\sum_{j'_n(r)=0}2^{-(n+1)}u_0\!\left(\rho'_n(r)\right)+\sum_{j'_n(r)\ge 1}2^{-(n+1)}\bigl[X^{j'_n(r)}u_0\bigr]_{\rho'_n(r)},$$
where \(q_r\) is the total probability of immediate ruin from residue \(r\).
All boundary-layer unknowns can therefore be grouped into a finite linear system
$$(I-C)u_0=q,$$
which the implementations solve by Gaussian elimination with pivoting.
For the given starting bankroll \(B\), first convert it to shifted capital
$$c_0=B-d.$$
If \(c_0 \le 0\), the answer is \(0\) because the first ticket cannot even be bought. Otherwise write
$$c_0 = k d + r + 1.$$
The ruin probability is then
$$\bigl[X^k u_0\bigr]_r,$$
so the desired survival probability is
$$\boxed{P_t(B)=1-\bigl[X^k u_0\bigr]_r.}$$
When \(t=2\), we have \(d=1\), so there is only one residue class. The transfer matrix becomes a single scalar \(x\), and the boundary ruin vector becomes a single scalar \(u\).
From layer \(1\), outcome \(n=0\) sends the process to the boundary layer, while \(n\ge 1\) jumps to layer \(2^n-1\). Therefore
$$x=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-1}.$$
For the boundary layer, \(n=0\) is immediate ruin and \(n\ge 1\) sends the process to layer \(2^n-2\), so
$$u=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-2}u.$$
Thus
$$P_2(2)=1-u,\qquad P_2(5)=1-x^3u.$$
The numerical solution gives
$$P_2(2)\approx 0.2522,\qquad P_2(5)\approx 0.6873,$$
which matches the validation checkpoints used by the implementations.
The C++, Python, and Java implementations all follow the same numerical plan. First they precompute the retained payouts \(2^n\) and their probabilities \(2^{-(n+1)}\) for \(0 \le n \lt 50\). Next they enumerate every transition from one reference layer and record only the resulting residue together with the layer jump that must be applied later.
To evaluate many powers of the same transfer matrix efficiently, the implementation builds \(X\), \(X^2\), \(X^4\), and so on, and reconstructs each required \(X^j\) from the binary expansion of \(j\). After the fixed-point iteration for \(X\) stabilizes, the boundary layer is folded into one finite linear system and solved by Gaussian elimination. Finally the initial bankroll is decomposed into layer and residue, the needed matrix power is applied to the boundary vector, and the result is subtracted from \(1\). The Python entry point is only a thin wrapper around the same compiled numerical core.
Let \(d=t-1\). Suppose the reference-layer and boundary transitions involve \(E\) distinct positive layer jumps, and let \(J_{\max}\) be the largest of them. One fixed-point iteration needs the squaring chain up to
$$L=\left\lfloor\log_2 J_{\max}\right\rfloor+1,$$
so its dense-matrix cost is \(O((L+E)d^3)\). Solving the boundary linear system costs \(O(d^3)\). Once the transfer matrix is known, evaluating the answer for a starting bankroll in layer \(k\) needs \(O(d^3 \log k)\) time by binary exponentiation. Memory usage is \(O((L+E)d^2)\). For the actual target \(t=15\), this is entirely manageable because \(d=14\).
In jeder Runde kostet ein Los \(t\) Einheiten. Die Auszahlung ist \(2^n\) mit Wahrscheinlichkeit \(2^{-(n+1)}\) für \(n \ge 0\). Ausgehend von einem Startkapital \(B\) suchen wir die Wahrscheinlichkeit, niemals unter die Grenze zu fallen, bei der das nächste Los nicht mehr gekauft werden kann. Weil das Kapital nach oben unbeschränkt wachsen kann, besitzt das Problem einen unendlichen Zustandsraum. Die Lösung ersetzt diese direkte unendliche Rekursion durch ein endliches Restklassenmodell plus eine selbstähnliche Matrixgleichung.
Alle drei Implementierungen benutzen dieselbe Wahrscheinlichkeitsstruktur. Der entscheidende Trick ist, das Kapital um \(t-1\) zu verschieben; dann wird der Prozess in Blöcken der Größe \(t-1\) translationsinvariant.
Setze
$$d=t-1,\qquad c=B-d.$$
Dabei misst \(c\) das Kapital oberhalb der kleinsten noch spielbaren Grenze. Eine Runde verändert \(c\) zu
$$c' = c + 2^n - t.$$
Falls \(c' \le 0\), kann kein weiteres Los mehr gekauft werden; dieser Pfad ist also ruiniert.
Die exakte Lotterie besitzt unendlich viele Ausgänge \(n=0,1,2,\dots\), aber die Implementierung behält nur
$$0 \le n \lt H,\qquad H=50.$$
Die abgeschnittene Restwahrscheinlichkeit ist
$$\sum_{n=H}^{\infty} 2^{-(n+1)} = 2^{-50},$$
also deutlich kleiner als die geforderte Genauigkeit von \(7\) Nachkommastellen.
Jedes spielbare verschobene Kapital \(c \ge 1\) lässt sich eindeutig schreiben als
$$c = k d + r + 1,\qquad k \ge 0,\qquad 0 \le r \lt d.$$
\(k\) ist die Schicht, \(r\) die Restklasse innerhalb dieser Schicht. Nützlich ist dies wegen
$$k d + r + 1 + 2^n - t = (k-1)d + (r + 2^n).$$
Erhöht man die Startschicht um eins, verschieben sich alle zukünftigen Schichten ebenfalls nur um eins. Die lokale Struktur innerhalb einer Schicht bleibt gleich.
Definiere
$$u_k(r)=\Pr(\text{eventual ruin} \mid c = k d + r + 1).$$
Für festes \(t\) ist jedes \(u_k\) also ein Vektor der Länge \(d\).
Wegen dieser Selbstähnlichkeit gibt es eine einzige \(d \times d\)-Matrix \(X\) mit
$$u_{k+1} = X u_k,\qquad u_k = X^k u_0.$$
Zur Herleitung von \(X\) starte man in Schicht \(1\), also bei \(c = d + r + 1\). Nach Ausgang \(n\) gilt
$$c' = d + r + 1 + 2^n - t = r + 2^n.$$
Schreibe den neuen Zustand als
$$c' = j_n(r)\,d + \rho_n(r) + 1,$$
wobei
$$j_n(r)=\left\lfloor\frac{r+2^n-1}{d}\right\rfloor,\qquad \rho_n(r)=(r+2^n-1)\bmod d.$$
Bei \(j_n(r)=0\) landet der Prozess in der Basisschicht. Bei \(j_n(r)\ge 1\) landet er \(j_n(r)\) Schichten darüber. Deshalb erfüllt die \(r\)-te Zeile von \(X\)
$$X_{r,*}=\sum_{j_n(r)=0} 2^{-(n+1)} e_{\rho_n(r)}^{\mathsf T}+\sum_{j_n(r)\ge 1} 2^{-(n+1)}\bigl(X^{j_n(r)}\bigr)_{\rho_n(r),*},$$
wobei \(e_{\rho}^{\mathsf T}\) der Zeilenvektor mit einer \(1\) in Spalte \(\rho\) und sonst nur Nullen ist.
Das ist eine Matrix-Fixpunktgleichung. Die Implementierungen iterieren sie numerisch, bis sich zwei aufeinanderfolgende Matrizen eintragsweise um weniger als \(10^{-18}\) unterscheiden oder die Sicherheitsgrenze der Iterationen erreicht ist.
Nun betrachten wir die niedrigste noch spielbare Schicht mit
$$c = r + 1,\qquad 0 \le r \lt d.$$
Nach Ausgang \(n\) ist
$$c' = r + 1 + 2^n - t = r + 2^n - d.$$
Falls \(c' \le 0\), tritt der Ruin sofort ein. Andernfalls schreiben wir
$$c' = j'_n(r)\,d + \rho'_n(r) + 1,$$
mit
$$j'_n(r)=\left\lfloor\frac{c'-1}{d}\right\rfloor,\qquad \rho'_n(r)=(c'-1)\bmod d.$$
Durch Konditionieren über die erste Runde erhält man
$$u_0(r)=q_r+\sum_{j'_n(r)=0}2^{-(n+1)}u_0\!\left(\rho'_n(r)\right)+\sum_{j'_n(r)\ge 1}2^{-(n+1)}\bigl[X^{j'_n(r)}u_0\bigr]_{\rho'_n(r)},$$
wobei \(q_r\) die Gesamtwahrscheinlichkeit des unmittelbaren Ruins aus Restklasse \(r\) ist.
Alle unbekannten Werte der Basisschicht bilden damit ein endliches lineares System
$$(I-C)u_0=q,$$
das per Gauß-Elimination mit Pivotisierung gelöst wird.
Für das gegebene Startkapital \(B\) berechnet man zuerst das verschobene Kapital
$$c_0=B-d.$$
Ist \(c_0 \le 0\), dann ist die Antwort \(0\), weil nicht einmal das erste Los gekauft werden kann. Andernfalls schreibe
$$c_0 = k d + r + 1.$$
Die Ruinwahrscheinlichkeit ist dann
$$\bigl[X^k u_0\bigr]_r,$$
also die gesuchte Überlebenswahrscheinlichkeit
$$\boxed{P_t(B)=1-\bigl[X^k u_0\bigr]_r.}$$
Für \(t=2\) gilt \(d=1\), also gibt es nur eine einzige Restklasse. Die Transfermatrix reduziert sich auf einen Skalar \(x\), und der Randvektor des Ruins auf einen Skalar \(u\).
Aus Schicht \(1\) führt \(n=0\) in die Basisschicht, während \(n\ge 1\) in die Schicht \(2^n-1\) springt. Daher
$$x=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-1}.$$
Für die Basisschicht bedeutet \(n=0\) sofortigen Ruin, und \(n\ge 1\) springt in die Schicht \(2^n-2\). Also
$$u=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-2}u.$$
Damit ist
$$P_2(2)=1-u,\qquad P_2(5)=1-x^3u.$$
Die numerische Lösung liefert
$$P_2(2)\approx 0.2522,\qquad P_2(5)\approx 0.6873,$$
genau die Kontrollwerte, die in den Implementierungen verwendet werden.
Die C++-, Python- und Java-Implementierungen folgen demselben numerischen Plan. Zuerst werden die berücksichtigten Auszahlungen \(2^n\) und die zugehörigen Wahrscheinlichkeiten \(2^{-(n+1)}\) für \(0 \le n \lt 50\) vorab berechnet. Danach werden alle Übergänge aus einer Referenzschicht aufgelistet; gespeichert werden nur die Zielrestklasse und der nötige Schichtsprung.
Um viele Potenzen derselben Transfermatrix effizient auszuwerten, baut die Implementierung \(X\), \(X^2\), \(X^4\) usw. auf und setzt jede benötigte Potenz \(X^j\) aus der Binärdarstellung von \(j\) zusammen. Sobald die Fixpunktiteration für \(X\) stabil ist, werden alle Übergänge aus der Basisschicht in ein endliches lineares System eingearbeitet und per Gauß-Elimination gelöst. Am Ende wird das Startkapital in Schicht und Restklasse zerlegt, die passende Matrixpotenz angewandt und das Ergebnis von \(1\) abgezogen. Der Python-Einstiegspunkt ist dabei nur eine schmale Hülle um denselben kompilierten numerischen Kern.
Sei \(d=t-1\). Angenommen, die Übergänge aus Referenz- und Randebene enthalten \(E\) verschiedene positive Schichtsprünge, und \(J_{\max}\) sei der größte davon. Eine Fixpunktiteration benötigt die Quadrierungskette bis
$$L=\left\lfloor\log_2 J_{\max}\right\rfloor+1,$$
und damit für dichte Matrizen \(O((L+E)d^3)\) Zeit. Das Lösen des Randsystems kostet \(O(d^3)\). Ist die Transfermatrix bekannt, so kostet die Auswertung für ein Startkapital in Schicht \(k\) mittels schneller Potenzierung \(O(d^3 \log k)\). Der Speicherverbrauch liegt bei \(O((L+E)d^2)\). Für den eigentlichen Zielwert \(t=15\) ist das problemlos, weil dann \(d=14\) ist.
Her turda bir bilet satın almak için \(t\) birim ödenir. Ödeme, \(n \ge 0\) için \(2^{-(n+1)}\) olasılıkla \(2^n\) şeklindedir. Başlangıç sermayesi \(B\) iken amaç, oyunu sonsuza kadar sürdürebilme olasılığını bulmaktır; bu da bir sonraki bileti alamayacak seviyenin altına hiç düşmemekle aynıdır. Sermaye teorik olarak sınırsız büyüyebildiği için durum uzayı sonsuzdur. Çözüm, ham sonsuz özyinelemeyi sonlu bir artık-sınıf modeli ve onun üzerine kurulu öz-benzer bir matris denklemiyle değiştirir.
Üç dildeki uygulamalar aynı olasılık modelini kullanır. Temel fikir, sermayeyi \(t-1\) kadar kaydırınca sürecin \(t-1\) uzunluklu bloklarda öteleme bakımından aynı davranmasıdır.
Şöyle tanımlayalım:
$$d=t-1,\qquad c=B-d.$$
Burada \(c\), sermayenin oynanabilir en alt eşik üzerindeki payını gösterir. Bir turun ardından
$$c' = c + 2^n - t$$
olur. Eğer \(c' \le 0\) ise oyuncu artık yeni bilet alamaz; bu yol iflasla biter.
Gerçek dağılım sonsuz sayıda sonuç içerir, ancak uygulama yalnızca
$$0 \le n \lt H,\qquad H=50$$
aralığını tutar. Atılan kuyruk olasılığı
$$\sum_{n=H}^{\infty} 2^{-(n+1)} = 2^{-50}$$
olduğundan, istenen \(7\) basamaklı doğruluk için tamamen ihmal edilebilir düzeydedir.
Her oynanabilir kaydırılmış sermaye \(c \ge 1\) tek şekilde
$$c = k d + r + 1,\qquad k \ge 0,\qquad 0 \le r \lt d$$
biçiminde yazılır. Burada \(k\) katmanı, \(r\) ise o katman içindeki artık sınıfını gösterir. Bunun yararı şu eşitlikten gelir:
$$k d + r + 1 + 2^n - t = (k-1)d + (r + 2^n).$$
Yani başlangıç katmanını bir artırmak, gelecekteki tüm katmanları da yalnızca bir kaydırır; katman içindeki ince yapı değişmez.
Şimdi
$$u_k(r)=\Pr(\text{eventual ruin} \mid c = k d + r + 1)$$
tanımını yapalım. Sabit \(t\) için her \(u_k\), uzunluğu \(d\) olan bir vektördür.
Bu öz-benzerlik nedeniyle tek bir \(d \times d\) boyutlu \(X\) matrisi vardır ve
$$u_{k+1} = X u_k,\qquad u_k = X^k u_0$$
sağlanır.
\(X\)'i türetmek için \(c = d + r + 1\) olan \(1\). katmandan başlayalım. \(n\) sonucu sonrası
$$c' = d + r + 1 + 2^n - t = r + 2^n$$
elde edilir. Bunu
$$c' = j_n(r)\,d + \rho_n(r) + 1$$
şeklinde yazalım; burada
$$j_n(r)=\left\lfloor\frac{r+2^n-1}{d}\right\rfloor,\qquad \rho_n(r)=(r+2^n-1)\bmod d.$$
Eğer \(j_n(r)=0\) ise süreç sınır katmanına iner. Eğer \(j_n(r)\ge 1\) ise sınırın \(j_n(r)\) kat üstüne çıkar. Bu yüzden \(X\)'in \(r\). satırı
$$X_{r,*}=\sum_{j_n(r)=0} 2^{-(n+1)} e_{\rho_n(r)}^{\mathsf T}+\sum_{j_n(r)\ge 1} 2^{-(n+1)}\bigl(X^{j_n(r)}\bigr)_{\rho_n(r),*}$$
denklemini sağlar. Burada \(e_{\rho}^{\mathsf T}\), yalnızca \(\rho\) sütununda \(1\) olan satır vektörüdür.
Bu bir matris sabit-nokta denklemidir. Uygulamalar bunu ardışık iki matris arasındaki en büyük fark \(10^{-18}\)'in altına inene kadar ya da güvenlik amaçlı iterasyon sınırına ulaşılana kadar nümerik olarak yineler.
Şimdi en alt oynanabilir katmana bakalım:
$$c = r + 1,\qquad 0 \le r \lt d.$$
\(n\) sonucu sonrası
$$c' = r + 1 + 2^n - t = r + 2^n - d$$
olur. Eğer \(c' \le 0\) ise iflas hemen gerçekleşir. Aksi halde
$$c' = j'_n(r)\,d + \rho'_n(r) + 1$$
ve
$$j'_n(r)=\left\lfloor\frac{c'-1}{d}\right\rfloor,\qquad \rho'_n(r)=(c'-1)\bmod d$$
yazılır.
İlk tur üzerinden koşullandırınca
$$u_0(r)=q_r+\sum_{j'_n(r)=0}2^{-(n+1)}u_0\!\left(\rho'_n(r)\right)+\sum_{j'_n(r)\ge 1}2^{-(n+1)}\bigl[X^{j'_n(r)}u_0\bigr]_{\rho'_n(r)}$$
elde edilir. Burada \(q_r\), \(r\) artık sınıfından anında iflas olasılığıdır.
Böylece tüm bilinmeyen sınır katmanı değerleri sonlu bir lineer sisteme toplanır:
$$(I-C)u_0=q.$$
Bu sistem Gauss eliminasyonu ve pivotlama ile çözülür.
Verilen başlangıç sermayesi \(B\) için önce kaydırılmış sermayeyi hesaplarız:
$$c_0=B-d.$$
Eğer \(c_0 \le 0\) ise ilk bilet bile alınamayacağı için cevap \(0\)'dır. Aksi halde
$$c_0 = k d + r + 1$$
yazarız. İflas olasılığı
$$\bigl[X^k u_0\bigr]_r$$
olduğundan, aranan hayatta kalma olasılığı
$$\boxed{P_t(B)=1-\bigl[X^k u_0\bigr]_r.}$$
\(t=2\) olduğunda \(d=1\) olur; dolayısıyla yalnızca tek bir artık sınıfı vardır. Aktarım matrisi tek skaler \(x\)'e, sınır iflas vektörü de tek skaler \(u\)'ya iner.
\(1\). katmandan başlayınca \(n=0\) sonucu sınır katmanına iner, \(n\ge 1\) ise \(2^n-1\). katmana sıçrar. Bu nedenle
$$x=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-1}.$$
Sınır katmanında \(n=0\) anında iflas, \(n\ge 1\) ise \(2^n-2\). katmana geçiş demektir. Bu yüzden
$$u=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-2}u.$$
Dolayısıyla
$$P_2(2)=1-u,\qquad P_2(5)=1-x^3u.$$
Nümerik çözüm
$$P_2(2)\approx 0.2522,\qquad P_2(5)\approx 0.6873$$
değerlerini verir; bunlar uygulamalardaki doğrulama kontrol noktalarıyla aynıdır.
C++, Python ve Java uygulamaları aynı nümerik planı izler. Önce \(0 \le n \lt 50\) için saklanan ödemeler \(2^n\) ve bunların olasılıkları \(2^{-(n+1)}\) önceden hesaplanır. Ardından tek bir referans katmandan çıkılabilecek tüm geçişler listelenir; yalnızca hedef artık sınıfı ve sonradan uygulanacak katman sıçraması saklanır.
Aynı aktarım matrisinin çok sayıda kuvvetini verimli hesaplamak için uygulama \(X\), \(X^2\), \(X^4\) gibi kuvvetleri oluşturur ve gerekli her \(X^j\)'yi \(j\)'nin ikili açılımından yeniden kurar. \(X\) için sabit-nokta iterasyonu kararlı hale geldikten sonra sınır katmanından daha yüksek katmanlara giden tüm katkılar tek bir sonlu lineer sisteme katlanır ve Gauss eliminasyonu ile çözülür. Son aşamada başlangıç sermayesi katman ve artık sınıfına ayrılır, uygun matris kuvveti uygulanır ve sonuç \(1\)'den çıkarılır. Python tarafı ise aynı derlenmiş nümerik çekirdeği çalıştıran ince bir sarmalayıcıdır.
\(d=t-1\) olsun. Referans katmanı ve sınır katmanı geçişlerinde \(E\) farklı pozitif katman sıçraması görülsün ve bunların en büyüğü \(J_{\max}\) olsun. Bir sabit-nokta iterasyonu,
$$L=\left\lfloor\log_2 J_{\max}\right\rfloor+1$$
uzunluğuna kadar kareleme zinciri gerektirir; bu yüzden yoğun matris maliyeti \(O((L+E)d^3)\) olur. Sınırdaki lineer sistemi çözmek \(O(d^3)\) zaman ister. Aktarım matrisi elde edildikten sonra, başlangıç sermayesi \(k\). katmandaysa sonucu ikili üs alma ile hesaplamak \(O(d^3 \log k)\) sürer. Bellek kullanımı \(O((L+E)d^2)\)'dir. Asıl hedef olan \(t=15\) için \(d=14\) olduğundan yöntem rahatlıkla uygulanabilir.
En cada ronda se paga \(t\) por un billete. El premio es \(2^n\) con probabilidad \(2^{-(n+1)}\) para \(n \ge 0\). Partiendo de un capital inicial \(B\), queremos la probabilidad de poder seguir jugando para siempre, es decir, de no caer nunca por debajo del importe necesario para comprar el siguiente billete. Como el capital puede crecer sin límite, el espacio de estados es infinito. La solución sustituye esa recursión infinita por un modelo finito de residuos junto con una ecuación matricial autosimilar.
Las tres implementaciones siguen la misma ley probabilística. La idea central es desplazar el capital en \(t-1\); después de ese cambio, el proceso se repite por bloques de tamaño \(t-1\).
Definimos
$$d=t-1,\qquad c=B-d.$$
Aquí \(c\) mide cuánto capital hay por encima del umbral mínimo para poder seguir jugando. Después de una ronda,
$$c' = c + 2^n - t.$$
Si \(c' \le 0\), el jugador ya no puede comprar otro billete y esa trayectoria queda arruinada.
La lotería exacta tiene infinitos resultados \(n=0,1,2,\dots\), pero la implementación conserva solo
$$0 \le n \lt H,\qquad H=50.$$
La masa de probabilidad descartada es
$$\sum_{n=H}^{\infty} 2^{-(n+1)} = 2^{-50},$$
muy por debajo de la precisión decimal pedida.
Cada capital desplazado jugable \(c \ge 1\) se escribe de forma única como
$$c = k d + r + 1,\qquad k \ge 0,\qquad 0 \le r \lt d.$$
El entero \(k\) es la capa y \(r\) es el residuo dentro de esa capa. Esto funciona porque
$$k d + r + 1 + 2^n - t = (k-1)d + (r + 2^n).$$
Subir una capa en el capital inicial solo desplaza una capa todas las posiciones futuras; la geometría interna permanece igual.
Definimos
$$u_k(r)=\Pr(\text{eventual ruin} \mid c = k d + r + 1).$$
Para un \(t\) fijo, cada \(u_k\) es un vector de longitud \(d\).
Como el patrón se repite de capa en capa, existe una única matriz \(d \times d\) llamada \(X\) tal que
$$u_{k+1} = X u_k,\qquad u_k = X^k u_0.$$
Para obtener \(X\), empezamos en la capa \(1\), donde \(c = d + r + 1\). Tras el resultado \(n\),
$$c' = d + r + 1 + 2^n - t = r + 2^n.$$
Escribimos el nuevo estado como
$$c' = j_n(r)\,d + \rho_n(r) + 1,$$
con
$$j_n(r)=\left\lfloor\frac{r+2^n-1}{d}\right\rfloor,\qquad \rho_n(r)=(r+2^n-1)\bmod d.$$
Si \(j_n(r)=0\), caemos en la capa frontera. Si \(j_n(r)\ge 1\), aterrizamos \(j_n(r)\) capas por encima. Por tanto, la fila \(r\) de \(X\) cumple
$$X_{r,*}=\sum_{j_n(r)=0} 2^{-(n+1)} e_{\rho_n(r)}^{\mathsf T}+\sum_{j_n(r)\ge 1} 2^{-(n+1)}\bigl(X^{j_n(r)}\bigr)_{\rho_n(r),*},$$
donde \(e_{\rho}^{\mathsf T}\) es el vector fila con un \(1\) en la columna \(\rho\) y ceros en las demás.
Esta es una ecuación matricial de punto fijo. Las implementaciones la iteran numéricamente hasta que dos matrices consecutivas difieren menos de \(10^{-18}\) entrada a entrada, o hasta alcanzar el tope de seguridad de iteraciones.
Ahora tomamos la capa jugable más baja, donde
$$c = r + 1,\qquad 0 \le r \lt d.$$
Después del resultado \(n\),
$$c' = r + 1 + 2^n - t = r + 2^n - d.$$
Si \(c' \le 0\), la ruina es inmediata. En caso contrario, escribimos
$$c' = j'_n(r)\,d + \rho'_n(r) + 1,$$
donde
$$j'_n(r)=\left\lfloor\frac{c'-1}{d}\right\rfloor,\qquad \rho'_n(r)=(c'-1)\bmod d.$$
Condicionando sobre la primera ronda se obtiene
$$u_0(r)=q_r+\sum_{j'_n(r)=0}2^{-(n+1)}u_0\!\left(\rho'_n(r)\right)+\sum_{j'_n(r)\ge 1}2^{-(n+1)}\bigl[X^{j'_n(r)}u_0\bigr]_{\rho'_n(r)},$$
donde \(q_r\) es la probabilidad total de ruina inmediata desde el residuo \(r\).
Agrupando todos los incógnitos de la capa frontera aparece un sistema lineal finito:
$$(I-C)u_0=q.$$
Las implementaciones lo resuelven con eliminación gaussiana y pivoteo.
Para el capital inicial \(B\), primero se convierte a capital desplazado:
$$c_0=B-d.$$
Si \(c_0 \le 0\), la respuesta es \(0\) porque ni siquiera puede comprarse el primer billete. En otro caso, se escribe
$$c_0 = k d + r + 1.$$
La probabilidad de ruina es
$$\bigl[X^k u_0\bigr]_r,$$
y por tanto la probabilidad buscada de supervivencia es
$$\boxed{P_t(B)=1-\bigl[X^k u_0\bigr]_r.}$$
Cuando \(t=2\), se tiene \(d=1\), así que solo existe una clase residual. La matriz de transferencia se reduce a un escalar \(x\), y el vector de ruina de frontera a un escalar \(u\).
Desde la capa \(1\), el resultado \(n=0\) lleva a la capa frontera, mientras que \(n\ge 1\) salta a la capa \(2^n-1\). Por ello
$$x=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-1}.$$
En la frontera, \(n=0\) significa ruina inmediata y \(n\ge 1\) lleva a la capa \(2^n-2\), de modo que
$$u=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-2}u.$$
Así,
$$P_2(2)=1-u,\qquad P_2(5)=1-x^3u.$$
La solución numérica da
$$P_2(2)\approx 0.2522,\qquad P_2(5)\approx 0.6873,$$
exactamente los valores de control usados por las implementaciones.
Las implementaciones en C++, Python y Java siguen el mismo esquema numérico. Primero precalculan los premios retenidos \(2^n\) y sus probabilidades \(2^{-(n+1)}\) para \(0 \le n \lt 50\). Después enumeran todas las transiciones salientes de una capa de referencia y almacenan solo el residuo de llegada y el salto de capa que habrá que aplicar más adelante.
Para evaluar con eficiencia muchas potencias de la misma matriz de transferencia, la implementación construye \(X\), \(X^2\), \(X^4\), etc., y recompone cada potencia \(X^j\) a partir de la expansión binaria de \(j\). Una vez estabilizada la iteración de punto fijo para \(X\), todas las contribuciones de la capa frontera a capas superiores se pliegan en un único sistema lineal finito y se resuelven por eliminación gaussiana. Por último, el capital inicial se descompone en capa y residuo, se aplica la potencia de matriz necesaria al vector de frontera y se resta el resultado de \(1\). La entrada en Python es solo un envoltorio ligero alrededor del mismo núcleo numérico compilado.
Sea \(d=t-1\). Supongamos que las transiciones desde la capa de referencia y desde la frontera contienen \(E\) saltos positivos de capa distintos, y sea \(J_{\max}\) el mayor de ellos. Una iteración de punto fijo necesita la cadena de cuadrados hasta
$$L=\left\lfloor\log_2 J_{\max}\right\rfloor+1,$$
por lo que el coste con matrices densas es \(O((L+E)d^3)\). Resolver el sistema lineal de la frontera cuesta \(O(d^3)\). Una vez conocida la matriz de transferencia, evaluar la respuesta para una banca inicial en la capa \(k\) requiere \(O(d^3 \log k)\) usando exponenciación binaria. La memoria es \(O((L+E)d^2)\). Para el objetivo real \(t=15\), esto es muy cómodo porque \(d=14\).
每一轮都要先支付 \(t\) 的票价，然后以概率 \(2^{-(n+1)}\) 获得奖金 \(2^n\)，其中 \(n \ge 0\)。给定初始资金 \(B\)，题目要求的是“永远还能继续买下一张票”的概率，也就是资金永远不会跌到无法继续游戏的水平之下。由于资金理论上可以无限增大，直接把资金值当作状态会得到一个无限状态马尔可夫过程。实现采用的办法不是硬算这个无限递推，而是把状态拆成有限个余数类，再配合一个具有自相似结构的矩阵方程来求解。
三种语言的实现使用的是同一套数值模型。核心观察是：把资金整体平移 \(t-1\) 之后，过程会按长度为 \(t-1\) 的区块重复出现，因此可以把无限状态压缩成“层号 + 余数”的形式。
定义
$$d=t-1,\qquad c=B-d.$$
这里的 \(c\) 表示资金相对于“还能继续买票的最低门槛”高出了多少。经过一轮之后，新的平移资金满足
$$c' = c + 2^n - t.$$
如果 \(c' \le 0\)，就说明下一张票已经买不起了，因此这条路径立刻算作破产。
真实的分布包含无限多个结果 \(n=0,1,2,\dots\)，但实现只保留
$$0 \le n \lt H,\qquad H=50.$$
被截掉的总概率是
$$\sum_{n=H}^{\infty} 2^{-(n+1)} = 2^{-50},$$
大约只有 \(8.88\times 10^{-16}\)，远小于题目最终要求的 \(7\) 位小数精度，因此这种截断不会影响最终答案。
任何仍然可玩的平移资金 \(c \ge 1\) 都可以唯一写成
$$c = k d + r + 1,\qquad k \ge 0,\qquad 0 \le r \lt d.$$
其中 \(k\) 表示第几层，\(r\) 表示该层内部的余数位置。这样写的好处在于
$$k d + r + 1 + 2^n - t = (k-1)d + (r + 2^n).$$
也就是说，初始层号整体增加 \(1\)，未来轨迹的层号也只是整体增加 \(1\)，层内结构本身并不会改变。这正是自相似性的来源。
定义
$$u_k(r)=\Pr(\text{eventual ruin} \mid c = k d + r + 1).$$
对于固定的票价 \(t\)，每个 \(u_k\) 都是一个长度为 \(d\) 的向量，记录该层每个余数位置的破产概率。
由于每层都重复同一种几何结构，所以存在一个统一的 \(d \times d\) 矩阵 \(X\)，使得
$$u_{k+1} = X u_k,\qquad u_k = X^k u_0.$$
要导出 \(X\)，先看第 \(1\) 层，也就是
$$c = d + r + 1.$$
若本轮结果是 \(n\)，则
$$c' = d + r + 1 + 2^n - t = r + 2^n.$$
把这个新状态再写成层号和余数的形式：
$$c' = j_n(r)\,d + \rho_n(r) + 1,$$
其中
$$j_n(r)=\left\lfloor\frac{r+2^n-1}{d}\right\rfloor,\qquad \rho_n(r)=(r+2^n-1)\bmod d.$$
如果 \(j_n(r)=0\)，那么新状态落在最底层；如果 \(j_n(r)\ge 1\)，那么新状态落在离底层 \(j_n(r)\) 层高的位置。因此，矩阵 \(X\) 的第 \(r\) 行满足
$$X_{r,*}=\sum_{j_n(r)=0} 2^{-(n+1)} e_{\rho_n(r)}^{\mathsf T}+\sum_{j_n(r)\ge 1} 2^{-(n+1)}\bigl(X^{j_n(r)}\bigr)_{\rho_n(r),*},$$
这里 \(e_{\rho}^{\mathsf T}\) 表示只有第 \(\rho\) 列为 \(1\)、其余列都为 \(0\) 的行向量。
这就是一个矩阵固定点方程。实现中通过数值迭代来逼近它：不断用当前矩阵去生成右边的新矩阵，直到相邻两次迭代的最大元素差小于 \(10^{-18}\)，或者达到安全的迭代上限为止。
最低可玩层对应
$$c = r + 1,\qquad 0 \le r \lt d.$$
在这一层中，若结果为 \(n\)，则
$$c' = r + 1 + 2^n - t = r + 2^n - d.$$
如果 \(c' \le 0\)，那就是立刻破产。否则仍可写成
$$c' = j'_n(r)\,d + \rho'_n(r) + 1,$$
其中
$$j'_n(r)=\left\lfloor\frac{c'-1}{d}\right\rfloor,\qquad \rho'_n(r)=(c'-1)\bmod d.$$
按第一轮分类后，就有
$$u_0(r)=q_r+\sum_{j'_n(r)=0}2^{-(n+1)}u_0\!\left(\rho'_n(r)\right)+\sum_{j'_n(r)\ge 1}2^{-(n+1)}\bigl[X^{j'_n(r)}u_0\bigr]_{\rho'_n(r)},$$
这里 \(q_r\) 表示从余数 \(r\) 出发“一步直接破产”的总概率。
把所有最低层未知量合并起来，就得到一个有限维线性系统
$$(I-C)u_0=q.$$
实现通过带主元选取的高斯消元来解这个系统。
给定初始资金 \(B\) 后，先转成平移资金
$$c_0=B-d.$$
若 \(c_0 \le 0\)，说明连第一张票都买不起，答案直接是 \(0\)。否则把它写成
$$c_0 = k d + r + 1.$$
对应的破产概率就是
$$\bigl[X^k u_0\bigr]_r,$$
所以要求的“不破产”概率为
$$\boxed{P_t(B)=1-\bigl[X^k u_0\bigr]_r.}$$
当 \(t=2\) 时，\(d=1\)，因此只有一个余数类。此时传递矩阵退化成一个标量 \(x\)，最低层破产向量也退化成一个标量 \(u\)。这正好把矩阵模型压缩成一维方程。
从第 \(1\) 层出发时，结果 \(n=0\) 会落到最低层，而 \(n\ge 1\) 会跳到第 \(2^n-1\) 层，因此有
$$x=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-1}.$$
对于最低层，\(n=0\) 表示立刻破产，\(n\ge 1\) 会跳到第 \(2^n-2\) 层，因此
$$u=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-2}u.$$
于是
$$P_2(2)=1-u,\qquad P_2(5)=1-x^3u.$$
数值求解得到
$$P_2(2)\approx 0.2522,\qquad P_2(5)\approx 0.6873,$$
与实现中使用的校验值完全一致。
C++、Python 和 Java 三个实现都遵循同一套数值流程。首先预计算保留下来的奖金 \(2^n\) 以及相应概率 \(2^{-(n+1)}\)，其中 \(0 \le n \lt 50\)。然后枚举一个参考层中所有可能的一步转移，只记录目标余数和需要跨越的层数，而不是保存整个无限状态空间。
为了高效求出许多不同的矩阵幂，实现会先构造 \(X\)、\(X^2\)、\(X^4\) 等二次幂，再根据指数的二进制展开拼出所需的每个 \(X^j\)。固定点矩阵收敛以后，最低层通往更高层的所有贡献都会被折叠进一个有限线性系统，再用高斯消元解出。最后把初始资金分解为层号和余数，应用相应的矩阵幂到最低层向量上，再用 \(1\) 减去破产概率。Python 版本本身只是一个轻量包装，它调用的是同样的编译后数值核心。
记 \(d=t-1\)。假设参考层和最低层转移中一共出现了 \(E\) 个不同的正层跳跃，最大者为 \(J_{\max}\)。一次固定点迭代需要构造到
$$L=\left\lfloor\log_2 J_{\max}\right\rfloor+1$$
为止的平方链，因此稠密矩阵部分的代价是 \(O((L+E)d^3)\)。最低层线性系统的求解需要 \(O(d^3)\)。当传递矩阵已经求出之后，对处在第 \(k\) 层的初始资金求答案，可通过二进制快速幂在 \(O(d^3 \log k)\) 时间内完成。内存使用量为 \(O((L+E)d^2)\)。对于目标参数 \(t=15\)，这里的 \(d=14\)，所以实际计算规模很小。
В каждом раунде билет стоит \(t\). Выплата равна \(2^n\) с вероятностью \(2^{-(n+1)}\) при \(n \ge 0\). При начальном капитале \(B\) требуется найти вероятность того, что игрок сможет продолжать игру бесконечно долго, то есть никогда не опустится ниже суммы, необходимой для покупки следующего билета. Поскольку капитал может расти без ограничения, прямое описание процесса приводит к бесконечному пространству состояний. Поэтому решение заменяет бесконечную рекурсию конечной моделью по остаткам и самоподобным матричным уравнением.
Во всех трех реализациях используется одна и та же вероятностная модель. Главная идея состоит в том, чтобы сдвинуть капитал на \(t-1\); после этого процесс повторяется блоками длины \(t-1\), и бесконечное пространство состояний удобно разбивается на слои.
Положим
$$d=t-1,\qquad c=B-d.$$
Величина \(c\) показывает, насколько капитал находится выше минимального playable-порога. После одного раунда имеем
$$c' = c + 2^n - t.$$
Если \(c' \le 0\), следующий билет купить уже нельзя, и такая траектория считается разорением.
Точное распределение содержит бесконечно много исходов \(n=0,1,2,\dots\), но реализация сохраняет только
$$0 \le n \lt H,\qquad H=50.$$
Отброшенная суммарная вероятность равна
$$\sum_{n=H}^{\infty} 2^{-(n+1)} = 2^{-50},$$
что намного меньше требуемой точности в \(7\) десятичных знаков.
Любой playable-сдвинутый капитал \(c \ge 1\) единственным образом представим в виде
$$c = k d + r + 1,\qquad k \ge 0,\qquad 0 \le r \lt d.$$
Число \(k\) есть номер слоя, а \(r\) задает положение внутри слоя. Это удобно, потому что
$$k d + r + 1 + 2^n - t = (k-1)d + (r + 2^n).$$
Если увеличить начальный слой на единицу, все будущие слои тоже просто сдвигаются на единицу; внутренняя структура переходов не меняется.
Определим
$$u_k(r)=\Pr(\text{eventual ruin} \mid c = k d + r + 1).$$
При фиксированном \(t\) каждый \(u_k\) является вектором длины \(d\).
Благодаря самоподобию существует одна матрица \(d \times d\), обозначим ее через \(X\), такая что
$$u_{k+1} = X u_k,\qquad u_k = X^k u_0.$$
Чтобы получить \(X\), рассмотрим слой \(1\), то есть состояние
$$c = d + r + 1.$$
После исхода \(n\) получаем
$$c' = d + r + 1 + 2^n - t = r + 2^n.$$
Запишем его в стандартной форме
$$c' = j_n(r)\,d + \rho_n(r) + 1,$$
где
$$j_n(r)=\left\lfloor\frac{r+2^n-1}{d}\right\rfloor,\qquad \rho_n(r)=(r+2^n-1)\bmod d.$$
Если \(j_n(r)=0\), процесс попадает в граничный слой. Если \(j_n(r)\ge 1\), он попадает на \(j_n(r)\) слоев выше границы. Поэтому строка \(r\) матрицы \(X\) удовлетворяет равенству
$$X_{r,*}=\sum_{j_n(r)=0} 2^{-(n+1)} e_{\rho_n(r)}^{\mathsf T}+\sum_{j_n(r)\ge 1} 2^{-(n+1)}\bigl(X^{j_n(r)}\bigr)_{\rho_n(r),*},$$
где \(e_{\rho}^{\mathsf T}\) есть строка с единицей в столбце \(\rho\) и нулями в остальных местах.
Это матричное уравнение на неподвижную точку. Реализации решают его численно итерациями, пока максимальное изменение по элементам не станет меньше \(10^{-18}\), либо пока не будет достигнут защитный предел числа итераций.
Теперь возьмем самый нижний playable-слой:
$$c = r + 1,\qquad 0 \le r \lt d.$$
После исхода \(n\)
$$c' = r + 1 + 2^n - t = r + 2^n - d.$$
Если \(c' \le 0\), разорение наступает сразу. В противном случае записываем
$$c' = j'_n(r)\,d + \rho'_n(r) + 1,$$
где
$$j'_n(r)=\left\lfloor\frac{c'-1}{d}\right\rfloor,\qquad \rho'_n(r)=(c'-1)\bmod d.$$
Условие по первому раунду дает
$$u_0(r)=q_r+\sum_{j'_n(r)=0}2^{-(n+1)}u_0\!\left(\rho'_n(r)\right)+\sum_{j'_n(r)\ge 1}2^{-(n+1)}\bigl[X^{j'_n(r)}u_0\bigr]_{\rho'_n(r)},$$
где \(q_r\) обозначает полную вероятность немедленного разорения из остатка \(r\).
Все неизвестные граничного слоя собираются в конечную линейную систему
$$(I-C)u_0=q,$$
которая решается методом Гаусса с выбором главного элемента.
Для заданного начального капитала \(B\) сначала вычисляем сдвинутый капитал
$$c_0=B-d.$$
Если \(c_0 \le 0\), ответ равен \(0\), потому что первый билет купить нельзя. Иначе записываем
$$c_0 = k d + r + 1.$$
Тогда вероятность разорения равна
$$\bigl[X^k u_0\bigr]_r,$$
а искомая вероятность никогда не разориться есть
$$\boxed{P_t(B)=1-\bigl[X^k u_0\bigr]_r.}$$
Если \(t=2\), то \(d=1\), и остается только один класс остатка. Тогда матрица переноса превращается в один скаляр \(x\), а граничный вектор разорения в один скаляр \(u\).
Из слоя \(1\) исход \(n=0\) переводит процесс в граничный слой, а \(n\ge 1\) уносит его в слой \(2^n-1\). Поэтому
$$x=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-1}.$$
Для граничного слоя \(n=0\) означает немедленное разорение, а \(n\ge 1\) ведет в слой \(2^n-2\), значит
$$u=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-2}u.$$
Следовательно,
$$P_2(2)=1-u,\qquad P_2(5)=1-x^3u.$$
Численное решение дает
$$P_2(2)\approx 0.2522,\qquad P_2(5)\approx 0.6873,$$
что совпадает с контрольными значениями, используемыми в реализациях.
Реализации на C++, Python и Java следуют одному и тому же численному плану. Сначала они предварительно вычисляют сохраненные выплаты \(2^n\) и вероятности \(2^{-(n+1)}\) для \(0 \le n \lt 50\). Затем перебираются все переходы из одного опорного слоя; сохраняются только целевой остаток и величина прыжка по слоям.
Чтобы эффективно получать много степеней одной и той же матрицы переноса, реализация строит \(X\), \(X^2\), \(X^4\) и так далее, а каждую нужную степень \(X^j\) собирает по двоичному разложению числа \(j\). После стабилизации итераций для \(X\) все переходы из граничного слоя на более высокие слои сворачиваются в одну конечную линейную систему, которая решается методом Гаусса. В конце начальный капитал раскладывается на слой и остаток, к граничному вектору применяется нужная степень матрицы, и результат вычитается из \(1\). Python-обертка только запускает тот же скомпилированный численный механизм.
Пусть \(d=t-1\). Предположим, что переходы из опорного и граничного слоев содержат \(E\) различных положительных прыжков по слоям, а \(J_{\max}\) есть максимальный из них. Одна итерация неподвижной точки требует цепочки квадратов до
$$L=\left\lfloor\log_2 J_{\max}\right\rfloor+1,$$
поэтому стоимость плотных матричных операций равна \(O((L+E)d^3)\). Решение граничной линейной системы занимает \(O(d^3)\). После нахождения матрицы переноса ответ для стартового капитала в слое \(k\) вычисляется за \(O(d^3 \log k)\) с помощью быстрого возведения в степень. Память имеет порядок \(O((L+E)d^2)\). Для целевого случая \(t=15\) имеем \(d=14\), так что размер вычислений невелик.
في كل جولة يدفع اللاعب ثمن تذكرة مقداره \(t\). ويكون العائد \(2^n\) باحتمال \(2^{-(n+1)}\) لكل \(n \ge 0\). إذا كان رأس المال الابتدائي هو \(B\)، فنحن نريد احتمال أن يبقى اللاعب قادرًا على الاستمرار إلى الأبد، أي ألا يهبط رأس المال أبدًا تحت الحد الأدنى اللازم لشراء التذكرة التالية. وبما أن رأس المال قد يكبر بلا حد، فإن فضاء الحالات لا نهائي. لذلك تستبدل الخوارزمية هذا الوصف اللانهائي المباشر بنموذج منتهٍ يعتمد على البواقي، مع معادلة مصفوفية ذات بنية ذاتية التشابه.
تعتمد تطبيقات C++ وPython وJava على النموذج الاحتمالي نفسه. والفكرة الأساسية هي إزاحة رأس المال بمقدار \(t-1\)، وعندها يصبح السلوك متكررًا على كتل طولها \(t-1\)، مما يسمح بضغط الفضاء اللانهائي إلى «طبقة + باقي».
نعرّف
$$d=t-1,\qquad c=B-d.$$
وهنا تمثل \(c\) مقدار رأس المال فوق أصغر حد يسمح بالاستمرار في اللعب. بعد جولة واحدة يصبح
$$c' = c + 2^n - t.$$
إذا كان \(c' \le 0\)، فلن يستطيع اللاعب شراء تذكرة أخرى، وتُعد هذه المسار حالة إفلاس.
التوزيع الحقيقي يحتوي عددًا لا نهائيًا من القيم \(n=0,1,2,\dots\)، لكن التنفيذ يحتفظ فقط بالقيم
$$0 \le n \lt H,\qquad H=50.$$
أما كتلة الاحتمال المحذوفة فهي
$$\sum_{n=H}^{\infty} 2^{-(n+1)} = 2^{-50},$$
وهي أصغر بكثير من مستوى الدقة المطلوب عند \(7\) منازل عشرية.
كل رأس مال مُزاح وقابل للعب بحيث \(c \ge 1\) يمكن كتابته بشكل وحيد على الصورة
$$c = k d + r + 1,\qquad k \ge 0,\qquad 0 \le r \lt d.$$
العدد \(k\) هو رقم الطبقة، و\(r\) هو الباقي داخل تلك الطبقة. وهذه الكتابة مفيدة لأن
$$k d + r + 1 + 2^n - t = (k-1)d + (r + 2^n).$$
أي إن زيادة الطبقة الابتدائية بمقدار واحد لا تغيّر البنية المحلية، بل تزيح جميع الطبقات المستقبلية بمقدار واحد فقط.
نعرّف
$$u_k(r)=\Pr(\text{eventual ruin} \mid c = k d + r + 1).$$
وعند ثبات \(t\)، يكون كل \(u_k\) متجهًا طوله \(d\).
بسبب هذا التشابه الذاتي، توجد مصفوفة واحدة بحجم \(d \times d\) نرمز لها بـ \(X\) بحيث
$$u_{k+1} = X u_k,\qquad u_k = X^k u_0.$$
لاشتقاق \(X\)، نبدأ من الطبقة \(1\)، أي من الحالة
$$c = d + r + 1.$$
بعد نتيجة \(n\) نحصل على
$$c' = d + r + 1 + 2^n - t = r + 2^n.$$
ثم نكتب الحالة الجديدة بالشكل
$$c' = j_n(r)\,d + \rho_n(r) + 1,$$
حيث
$$j_n(r)=\left\lfloor\frac{r+2^n-1}{d}\right\rfloor,\qquad \rho_n(r)=(r+2^n-1)\bmod d.$$
إذا كان \(j_n(r)=0\)، فإننا نهبط إلى الطبقة الحدّية. وإذا كان \(j_n(r)\ge 1\)، فإننا نصل إلى طبقة أعلى بمقدار \(j_n(r)\) فوق الحد. لذلك تحقق السطر \(r\) من المصفوفة \(X\) العلاقة
$$X_{r,*}=\sum_{j_n(r)=0} 2^{-(n+1)} e_{\rho_n(r)}^{\mathsf T}+\sum_{j_n(r)\ge 1} 2^{-(n+1)}\bigl(X^{j_n(r)}\bigr)_{\rho_n(r),*},$$
حيث \(e_{\rho}^{\mathsf T}\) هو متجه صفّي يحتوي \(1\) في العمود \(\rho\) وصفرًا في غيره.
هذه معادلة نقطة ثابتة لمصفوفة. وتقوم التطبيقات بحلها عدديًا بالتكرار حتى يصبح أكبر فرق عنصرًا بعنصر بين مصفوفتين متتاليتين أقل من \(10^{-18}\)، أو حتى الوصول إلى حد أمان لعدد التكرارات.
الآن ننظر إلى أدنى طبقة ما زالت تسمح باللعب، حيث
$$c = r + 1,\qquad 0 \le r \lt d.$$
بعد النتيجة \(n\) يصبح
$$c' = r + 1 + 2^n - t = r + 2^n - d.$$
إذا كان \(c' \le 0\)، فالإفلاس فوري. وإلا نكتب
$$c' = j'_n(r)\,d + \rho'_n(r) + 1,$$
مع
$$j'_n(r)=\left\lfloor\frac{c'-1}{d}\right\rfloor,\qquad \rho'_n(r)=(c'-1)\bmod d.$$
وبأخذ الاحتمال الشرطي على الجولة الأولى نحصل على
$$u_0(r)=q_r+\sum_{j'_n(r)=0}2^{-(n+1)}u_0\!\left(\rho'_n(r)\right)+\sum_{j'_n(r)\ge 1}2^{-(n+1)}\bigl[X^{j'_n(r)}u_0\bigr]_{\rho'_n(r)},$$
حيث \(q_r\) هو احتمال الإفلاس الفوري انطلاقًا من الباقي \(r\).
وعند جمع مجاهيل الطبقة الحدّية كلها نحصل على نظام خطي منتهٍ:
$$(I-C)u_0=q.$$
وتحله التطبيقات بحذف غاوسي مع اختيار محور مناسب.
بالنسبة إلى رأس المال الابتدائي \(B\)، نحسب أولًا رأس المال المزاح
$$c_0=B-d.$$
إذا كان \(c_0 \le 0\)، فالإجابة هي \(0\) لأن أول تذكرة نفسها غير ممكنة. وإلا نكتب
$$c_0 = k d + r + 1.$$
وعندئذ يكون احتمال الإفلاس
$$\bigl[X^k u_0\bigr]_r,$$
ومن ثم فإن احتمال عدم الإفلاس أبدًا هو
$$\boxed{P_t(B)=1-\bigl[X^k u_0\bigr]_r.}$$
إذا كان \(t=2\)، فإن \(d=1\)، وبالتالي لا يبقى إلا صنف باقي واحد. عندها تختزل مصفوفة الانتقال إلى عدد وحيد \(x\)، ويختزل متجه الإفلاس الحدّي إلى عدد وحيد \(u\).
من الطبقة \(1\)، تؤدي النتيجة \(n=0\) إلى الطبقة الحدّية، بينما تؤدي كل نتيجة \(n\ge 1\) إلى الطبقة \(2^n-1\). لذلك
$$x=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-1}.$$
أما في الطبقة الحدّية، فإن \(n=0\) يعني إفلاسًا فوريًا، و\(n\ge 1\) ينقل العملية إلى الطبقة \(2^n-2\)، ومن ثم
$$u=\frac12+\sum_{n=1}^{H-1}2^{-(n+1)}x^{2^n-2}u.$$
وبالتالي
$$P_2(2)=1-u,\qquad P_2(5)=1-x^3u.$$
ويعطي الحل العددي
$$P_2(2)\approx 0.2522,\qquad P_2(5)\approx 0.6873,$$
وهي نفس قيم التحقق المستخدمة داخل التطبيقات.
تتبع تطبيقات C++ وPython وJava الخطة العددية نفسها. أولًا تُحسب مسبقًا الجوائز المحتفَظ بها \(2^n\) واحتمالاتها \(2^{-(n+1)}\) من أجل \(0 \le n \lt 50\). بعد ذلك تُحصى جميع الانتقالات الخارجة من طبقة مرجعية واحدة، ويُخزَّن فقط الباقي الهدف ومقدار القفزة بين الطبقات.
ولكي تُحسب قوى كثيرة للمصفوفة نفسها بكفاءة، تبني الخوارزمية السلسلة \(X\)، \(X^2\)، \(X^4\)، ... ثم تعيد تركيب كل قوة مطلوبة \(X^j\) من التمثيل الثنائي للأس \(j\). وبعد أن تستقر عملية التكرار الخاصة بالنقطة الثابتة، تُطوى كل مساهمات الطبقة الحدّية إلى الطبقات الأعلى داخل نظام خطي منتهٍ، ثم يُحل بحذف غاوسي. وفي المرحلة الأخيرة يُفكك رأس المال الابتدائي إلى طبقة وباقٍ، وتُطبق قوة المصفوفة المناسبة على المتجه الحدّي، ثم يُطرح الناتج من \(1\). أما مدخل Python فهو مجرد غلاف خفيف يشغّل النواة العددية المترجمة نفسها.
ليكن \(d=t-1\). افترض أن انتقالات الطبقة المرجعية والطبقة الحدّية تولد \(E\) قفزات موجبة مختلفة بين الطبقات، وأن أكبرها هو \(J_{\max}\). عندئذ تحتاج كل خطوة من تكرار النقطة الثابتة إلى سلسلة تربيع حتى
$$L=\left\lfloor\log_2 J_{\max}\right\rfloor+1,$$
ولذلك تكون كلفة المصفوفات الكثيفة \(O((L+E)d^3)\). أما حل النظام الخطي الحدّي فيكلف \(O(d^3)\). وبعد معرفة مصفوفة الانتقال، فإن تقييم الجواب لرأس مال ابتدائي يقع في الطبقة \(k\) يحتاج \(O(d^3 \log k)\) باستخدام الرفع الثنائي للأس. ويبلغ استهلاك الذاكرة \(O((L+E)d^2)\). وفي الحالة الهدف \(t=15\) يكون \(d=14\)، ولهذا يبقى الحساب صغيرًا ومريحًا.