Problem Summary

A standard deck has 13 ranks and 4 suits, so each rank appears four times. Cards are drawn without replacement until either two consecutive cards share a rank or the deck runs out. The quantity to compute is the expected total number of draws. Since only rank equality matters, suits never need to be tracked individually.

Mathematical Approach

It is convenient to analyze a generalized deck with \(r\) ranks and \(c\) copies of each rank, then substitute \(r=13\) and \(c=4\) at the end. The key observation is that after at least one card has been drawn, the future depends only on how many copies remain of the last rank and on how the other ranks are distributed by remaining multiplicity.

Step 1: Compress the deck into a rank histogram

After any nonterminal history, call the rank of the most recently drawn card the current rank. Let \(m\in\{0,1,\dots,c-1\}\) be the number of undealt cards of that current rank. For every \(j\in\{0,1,\dots,c\}\), let \(a_j\) be the number of other ranks with exactly \(j\) undealt cards.

Then \((m,a_0,\dots,a_c)\) is a complete Markov state. The identities of the ranks do not matter; only the counts matter, because ranks with the same remaining multiplicity are indistinguishable for every future probability calculation.

For the actual problem \(c=4\), so the histogram has buckets \(a_0,a_1,a_2,a_3,a_4\), and the non-current ranks always satisfy

$$a_0+a_1+a_2+a_3+a_4=12.$$

Step 2: Compute the remaining cards and the stopping probability

If the current state is \((m,\mathbf{a})\), the number of undealt cards is

$$R(m,\mathbf{a})=m+\sum_{j=1}^{c} j a_j.$$

The next draw creates a consecutive pair exactly when it comes from the current rank, which has \(m\) surviving copies. Therefore

$$\Pr(\text{stop on next draw}\mid m,\mathbf{a})=\frac{m}{R(m,\mathbf{a})}.$$

If \(m=0\), the last drawn rank is already exhausted, so the next draw cannot stop the process immediately. The same formula still gives probability \(0\).

Step 3: Describe the non-stopping transitions

Suppose the next card comes from a different rank that currently belongs to bucket \(j\), meaning that this rank has \(j\) undealt copies before the draw. There are \(a_j\) such ranks and \(j a_j\) eligible cards, so

$$\Pr(\text{choose a rank from bucket }j\mid m,\mathbf{a})=\frac{j a_j}{R(m,\mathbf{a})},\qquad 1\le j\le c.$$

After that draw, the newly drawn rank becomes the new current rank and has \(j-1\) cards left. The previous current rank stops being current and re-enters the histogram with \(m\) cards left. Hence the new state is

$$\left(j-1,\ \mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right),$$

where \(\mathbf{e}_k\) is the unit vector of the histogram. This update also handles the case \(m=0\): the exhausted previous current rank is simply recorded in the zero bucket.

Step 4: Write the expectation recurrence

Let \(E(m,\mathbf{a})\) denote the expected number of additional draws from state \((m,\mathbf{a})\). If no cards remain, then

$$E(m,\mathbf{a})=0\qquad\text{when }R(m,\mathbf{a})=0.$$

Otherwise one new card is certainly drawn now, contributing the leading \(1\). If that card matches the current rank, the process stops immediately, so there is no recursive term for that branch. Summing only over the non-stopping branches gives

$$E(m,\mathbf{a})=1+\sum_{j=1}^{c}\frac{j a_j}{R(m,\mathbf{a})}\,E\!\left(j-1,\mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right).$$

This is the full dynamic program. Every transition removes exactly one undealt card, so the recursion always moves toward smaller states.

Step 5: Initialize the real 52-card deck

Before any draw there is no current rank, so the recursion is started after the first card is taken. In a deck with \(r\) ranks and \(c\) copies per rank, the first draw leaves

$$m=c-1,\qquad a_c=r-1,\qquad a_j=0\ \text{for }j\ne c.$$

Therefore the total expected number of draws is

$$1+E\!\left(c-1,\mathbf{a}^{(0)}\right).$$

For the Project Euler deck, \(r=13\) and \(c=4\), so

$$\mathbb{E}[\text{total draws}]=1+E(3,a_0=0,a_1=0,a_2=0,a_3=0,a_4=12).$$

The first transition from this state already explains the simple sanity check

$$\Pr(\text{stop on draw 2})=\frac{3}{51}=\frac{1}{17}.$$

Worked Example: two ranks with two copies each

A very small deck makes the recurrence transparent. Take \(r=2\) and \(c=2\). After the first draw, the state is \(m=1\) and \(a_2=1\), so let

$$X=E(1,a_0=0,a_1=0,a_2=1).$$

There are \(R=3\) undealt cards. With probability \(1/3\) the next card matches the current rank and the process ends after one more draw. With probability \(2/3\) we move to the state \(E(1,a_0=0,a_1=1,a_2=0)\). Hence

$$X=1+\frac{2}{3}E(1,a_0=0,a_1=1,a_2=0).$$

Now let

$$Y=E(1,a_0=0,a_1=1,a_2=0).$$

Then \(R=2\), so

$$Y=1+\frac{1}{2}E(0,a_0=0,a_1=1,a_2=0).$$

Finally, from \(E(0,a_0=0,a_1=1,a_2=0)\) only one card remains, it cannot match the exhausted current rank, and drawing it empties the deck, so

$$E(0,a_0=0,a_1=1,a_2=0)=1.$$

Therefore \(Y=1+\frac12=\frac32\), then \(X=1+\frac23\cdot\frac32=2\), and the full expectation is

$$1+X=3.$$

This matches the small checkpoint used by the implementations and confirms the recurrence.

How the Code Works

The C++, Python, and Java implementations all use the same memoized recursion. They encode the current remainder together with the histogram of remaining rank multiplicities into a compact state key, so each reachable state is evaluated once. For a queried state, the implementation first computes \(R\), returns \(0\) if the deck is empty, and otherwise starts from the mandatory next draw.

It then iterates through the possible buckets \(j=1,\dots,c\). Whenever bucket \(j\) is nonempty, the implementation applies the probability \(\frac{j a_j}{R}\), updates the histogram to reflect the new current rank and the old current rank returning to the bucket structure, and recursively asks for the expectation of the new state. Memoization turns the naive recursion tree into a compact dynamic program over the reachable histogram states.

For the actual input \((r,c)=(13,4)\), the program prints the expectation to eight decimal places. Floating-point arithmetic is sufficient here because the state graph is small and the recurrence is numerically stable for this instance.

Complexity Analysis

For a generalized deck with \(r\) ranks and \(c\) copies per rank, the \(r-1\) non-current ranks can be distributed among the \(c+1\) buckets in at most

$$\binom{r+c-1}{c}$$

ways. The current rank contributes \(c\) possible values of \(m\in\{0,\dots,c-1\}\), so an upper bound on the number of abstract states is

$$|\mathcal{S}|\le c\binom{r+c-1}{c}.$$

Each state considers at most \(c\) transitions, so the memoized running time is \(O(c|\mathcal{S}|)\), and the memory usage is \(O(|\mathcal{S}|)\). For the real deck this gives at most

$$4\binom{16}{4}=7280$$

states, which is tiny.

Footnotes and References

  1. Project Euler problem page: https://projecteuler.net/problem=856
  2. Markov chain: Wikipedia — Markov chain
  3. Expected value: Wikipedia — Expected value
  4. Hypergeometric distribution: Wikipedia — Hypergeometric distribution
  5. Dynamic programming: Wikipedia — Dynamic programming

Problemzusammenfassung

Ein Standardkartenspiel hat 13 Ränge und 4 Farben, also kommt jeder Rang genau viermal vor. Es wird ohne Zurücklegen gezogen, bis entweder zwei aufeinanderfolgende Karten denselben Rang haben oder das Deck leer ist. Gesucht ist der Erwartungswert der gesamten Ziehungszahl. Da nur Ranggleichheit zählt, müssen die Farben nicht einzeln verfolgt werden.

Mathematischer Ansatz

Es ist hilfreich, zunächst ein verallgemeinertes Deck mit \(r\) Rängen und \(c\) Karten pro Rang zu betrachten und erst am Ende \(r=13\) und \(c=4\) einzusetzen. Die Schlüsselerkenntnis lautet: Nach mindestens einer gezogenen Karte hängt die Zukunft nur noch davon ab, wie viele Karten vom zuletzt gezogenen Rang übrig sind und wie die übrigen Ränge nach Restanzahl verteilt sind.

Schritt 1: Das Deck als Rang-Histogramm komprimieren

Nach jeder nichtterminalen Historie nennen wir den Rang der zuletzt gezogenen Karte den aktuellen Rang. Sei \(m\in\{0,1,\dots,c-1\}\) die Zahl der noch nicht gezogenen Karten dieses aktuellen Rangs. Für jedes \(j\in\{0,1,\dots,c\}\) sei \(a_j\) die Anzahl der anderen Ränge, von denen genau \(j\) Karten übrig sind.

Dann ist \((m,a_0,\dots,a_c)\) ein vollständiger Markov-Zustand. Welche konkreten Ränge das sind, spielt keine Rolle; entscheidend sind nur die Häufigkeiten, weil Ränge mit derselben Restanzahl für alle zukünftigen Wahrscheinlichkeiten austauschbar sind.

Im eigentlichen Problem ist \(c=4\). Das Histogramm besitzt also die Körbe \(a_0,a_1,a_2,a_3,a_4\), und für die nicht aktuellen Ränge gilt immer

$$a_0+a_1+a_2+a_3+a_4=12.$$

Schritt 2: Restkartenzahl und Stoppwahrscheinlichkeit bestimmen

Im Zustand \((m,\mathbf{a})\) ist die Zahl der noch nicht gezogenen Karten

$$R(m,\mathbf{a})=m+\sum_{j=1}^{c} j a_j.$$

Der nächste Zug erzeugt genau dann sofort ein benachbartes Paar, wenn er vom aktuellen Rang kommt, von dem noch \(m\) Karten existieren. Daher

$$\Pr(\text{Stopp im nächsten Zug}\mid m,\mathbf{a})=\frac{m}{R(m,\mathbf{a})}.$$

Falls \(m=0\), ist der zuletzt gezogene Rang bereits erschöpft, also kann der nächste Zug nicht sofort stoppen. Dieselbe Formel liefert dann korrekt die Wahrscheinlichkeit \(0\).

Schritt 3: Die nicht stoppenden Übergänge beschreiben

Angenommen, die nächste Karte stammt von einem anderen Rang aus Korb \(j\), also von einem Rang mit \(j\) Restkarten vor dem Zug. Es gibt \(a_j\) solche Ränge und somit \(j a_j\) passende Karten. Deshalb gilt

$$\Pr(\text{Wahl eines Rangs aus Korb }j\mid m,\mathbf{a})=\frac{j a_j}{R(m,\mathbf{a})},\qquad 1\le j\le c.$$

Nach diesem Zug wird der neu getroffene Rang zum aktuellen Rang und hat noch \(j-1\) Karten. Der bisher aktuelle Rang wird wieder gewöhnlich und kehrt mit \(m\) Restkarten ins Histogramm zurück. Der neue Zustand lautet also

$$\left(j-1,\ \mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right),$$

wobei \(\mathbf{e}_k\) den Einheitsvektor des Histogramms bezeichnet. Auch der Fall \(m=0\) ist damit abgedeckt: Der alte aktuelle Rang landet dann einfach im Null-Korb.

Schritt 4: Die Erwartungsrekursion aufstellen

Sei \(E(m,\mathbf{a})\) die erwartete Zahl zusätzlicher Züge aus dem Zustand \((m,\mathbf{a})\). Wenn keine Karten mehr übrig sind, dann gilt

$$E(m,\mathbf{a})=0\qquad\text{falls }R(m,\mathbf{a})=0.$$

Andernfalls wird jetzt sicher genau eine weitere Karte gezogen; das erklärt die führende \(1\). Falls diese Karte den aktuellen Rang trifft, endet der Prozess sofort, daher gibt es für diesen Erfolgszweig keinen Rekursionsterm. Es bleiben nur die nicht stoppenden Äste:

$$E(m,\mathbf{a})=1+\sum_{j=1}^{c}\frac{j a_j}{R(m,\mathbf{a})}\,E\!\left(j-1,\mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right).$$

Das ist bereits das vollständige dynamische Programm. Jeder Übergang entfernt genau eine noch nicht gezogene Karte, also bewegt sich die Rekursion immer zu kleineren Zuständen.

Schritt 5: Das reale 52-Karten-Deck initialisieren

Vor dem ersten Zug gibt es noch keinen aktuellen Rang. Deshalb startet die Rekursion erst nach der ersten gezogenen Karte. Für ein Deck mit \(r\) Rängen und \(c\) Karten pro Rang bleibt dann

$$m=c-1,\qquad a_c=r-1,\qquad a_j=0\ \text{für }j\ne c.$$

Somit ist die gesuchte gesamte Erwartung

$$1+E\!\left(c-1,\mathbf{a}^{(0)}\right).$$

Für das Project-Euler-Deck mit \(r=13\) und \(c=4\) folgt

$$\mathbb{E}[\text{Gesamtzahl der Ziehungen}]=1+E(3,a_0=0,a_1=0,a_2=0,a_3=0,a_4=12).$$

Schon der erste Übergang aus diesem Zustand erklärt die einfache Plausibilitätsprüfung

$$\Pr(\text{Stopp bei Zug 2})=\frac{3}{51}=\frac{1}{17}.$$

Durchgerechnetes Beispiel: zwei Ränge mit je zwei Karten

Ein sehr kleines Deck macht die Rekursion transparent. Nehmen wir \(r=2\) und \(c=2\). Nach der ersten Karte ist der Zustand \(m=1\) und \(a_2=1\). Setzen wir

$$X=E(1,a_0=0,a_1=0,a_2=1).$$

Es bleiben \(R=3\) Karten. Mit Wahrscheinlichkeit \(1/3\) trifft die nächste Karte den aktuellen Rang und der Prozess endet nach genau einem weiteren Zug. Mit Wahrscheinlichkeit \(2/3\) gelangen wir in den Zustand \(E(1,a_0=0,a_1=1,a_2=0)\). Also

$$X=1+\frac{2}{3}E(1,a_0=0,a_1=1,a_2=0).$$

Definieren wir nun

$$Y=E(1,a_0=0,a_1=1,a_2=0).$$

Dann ist \(R=2\) und damit

$$Y=1+\frac{1}{2}E(0,a_0=0,a_1=1,a_2=0).$$

Im Zustand \(E(0,a_0=0,a_1=1,a_2=0)\) bleibt nur noch eine Karte. Sie kann den erschöpften aktuellen Rang nicht treffen, und nach diesem Zug ist das Deck leer. Deshalb gilt

$$E(0,a_0=0,a_1=1,a_2=0)=1.$$

Also ist \(Y=1+\frac12=\frac32\), daraus \(X=1+\frac23\cdot\frac32=2\), und die gesamte Erwartung lautet

$$1+X=3.$$

Genau dieser kleine Kontrollwert wird auch von den Implementierungen bestätigt.

Wie der Code funktioniert

Die C++-, Python- und Java-Implementierungen verwenden dieselbe memoisiert rekursive Idee. Der aktuelle Restwert und das Histogramm der verbleibenden Rangmultiplizitäten werden zu einem kompakten Zustandsschlüssel kodiert, sodass jeder erreichbare Zustand nur einmal ausgewertet wird. Für einen Zustand berechnet die Implementierung zunächst \(R\), gibt \(0\) zurück, falls das Deck leer ist, und zählt ansonsten den zwingend folgenden nächsten Zug.

Danach werden die möglichen Körbe \(j=1,\dots,c\) durchlaufen. Immer wenn Korb \(j\) nicht leer ist, verwendet die Implementierung die Wahrscheinlichkeit \(\frac{j a_j}{R}\), aktualisiert das Histogramm für den neuen aktuellen Rang und den zurückkehrenden alten aktuellen Rang und ruft die Erwartung des Folgezustands rekursiv ab. Durch Memoisierung wird aus dem naiven Rekursionsbaum ein kleines dynamisches Programm über die erreichbaren Histogrammzustände.

Für die eigentliche Eingabe \((r,c)=(13,4)\) wird der Erwartungswert mit acht Nachkommastellen ausgegeben. Gleitkommaarithmetik reicht hier aus, weil der Zustandsgraph klein und die Rekursion numerisch gutartig ist.

Komplexitätsanalyse

Für ein verallgemeinertes Deck mit \(r\) Rängen und \(c\) Karten pro Rang lassen sich die \(r-1\) nicht aktuellen Ränge auf höchstens

$$\binom{r+c-1}{c}$$

verschiedene Arten auf die \(c+1\) Körbe verteilen. Der aktuelle Rang liefert zusätzlich \(c\) mögliche Werte \(m\in\{0,\dots,c-1\}\). Also gilt als obere Schranke für die Zahl abstrakter Zustände

$$|\mathcal{S}|\le c\binom{r+c-1}{c}.$$

Pro Zustand werden höchstens \(c\) Übergänge betrachtet. Damit beträgt die memoisiert Laufzeit \(O(c|\mathcal{S}|)\), und der Speicherbedarf ist \(O(|\mathcal{S}|)\). Für das reale Deck ergibt sich höchstens

$$4\binom{16}{4}=7280$$

Zustände, also sehr wenig.

Fußnoten und Referenzen

  1. Project-Euler-Problemseite: https://projecteuler.net/problem=856
  2. Markov-Kette: Wikipedia — Markov chain
  3. Erwartungswert: Wikipedia — Expected value
  4. Hypergeometrische Verteilung: Wikipedia — Hypergeometric distribution
  5. Dynamische Programmierung: Wikipedia — Dynamic programming

Problem Özeti

Standart bir destede 13 rütbe vardır ve her rütbeden 4 kart bulunur. Kartlar geri koymadan çekilir; ardışık iki kartın rütbesi aynı olduğunda süreç hemen durur, aksi halde deste bitene kadar devam eder. Hesaplanması istenen büyüklük toplam çekiş sayısının beklenen değeridir. Sadece rütbe eşitliği önemli olduğu için renk veya tür bilgisi ayrı ayrı izlenmez.

Matematiksel Yaklaşım

Önce \(r\) farklı rütbeli ve her rütbeden \(c\) kopyalı genelleştirilmiş bir deste düşünmek, sonra en sonda \(r=13\) ve \(c=4\) koymak uygundur. Temel gözlem şudur: En az bir kart çekildikten sonra geleceği belirleyen tek şey, son çekilen rütbeden kaç kart kaldığı ve diğer rütbelerin kalan kart sayılarına göre nasıl dağıldığıdır.

Adım 1: Desteyi rütbe histogramına indirgeme

Henüz bitmemiş herhangi bir geçmişten sonra, en son çekilen kartın rütbesine mevcut rütbe diyelim. \(m\in\{0,1,\dots,c-1\}\), bu mevcut rütbeden destede kalan kart sayısı olsun. Her \(j\in\{0,1,\dots,c\}\) için \(a_j\), diğer rütbeler arasında tam olarak \(j\) kartı kalan rütbe sayısını göstersin.

Böylece \((m,a_0,\dots,a_c)\) tam bir Markov durumudur. Hangi somut rütbenin hangi kovada olduğunun önemi yoktur; sadece sayılar önemlidir, çünkü aynı kalan çokluğa sahip rütbeler gelecekteki tüm olasılıklar açısından eşdeğerdir.

Gerçek problemde \(c=4\) olduğu için kovalar \(a_0,a_1,a_2,a_3,a_4\) olur ve mevcut olmayan 12 rütbe için her zaman

$$a_0+a_1+a_2+a_3+a_4=12$$

eşitliği geçerlidir.

Adım 2: Kalan kart sayısı ve durma olasılığı

\((m,\mathbf{a})\) durumunda destede kalan toplam kart sayısı

$$R(m,\mathbf{a})=m+\sum_{j=1}^{c} j a_j$$

şeklindedir. Bir sonraki çekiş ancak mevcut rütbeden gelirse ardışık eşleşme üretir; bu rütbeden de \(m\) adet kart kalmıştır. Dolayısıyla

$$\Pr(\text{bir sonraki çekişte durma}\mid m,\mathbf{a})=\frac{m}{R(m,\mathbf{a})}.$$

Eğer \(m=0\) ise son çekilen rütbe zaten tükenmiştir; bu yüzden bir sonraki kart süreci anında bitiremez. Aynı formül bu durumda doğru olarak \(0\) verir.

Adım 3: Durmayan geçişlerin yapısı

Şimdi sonraki kartın, çekişten önce \(j\) kartı kalmış başka bir rütbeden geldiğini düşünelim. Böyle \(a_j\) tane rütbe ve toplam \(j a_j\) uygun kart vardır. Bu nedenle

$$\Pr(\text{\(j\) kovasından bir rütbe seçilmesi}\mid m,\mathbf{a})=\frac{j a_j}{R(m,\mathbf{a})},\qquad 1\le j\le c.$$

Bu çekişten sonra yeni çekilen rütbe mevcut rütbe olur ve ondan \(j-1\) kart kalır. Eski mevcut rütbe artık mevcut değildir ve \(m\) kartıyla yeniden histogram içine döner. Yeni durum bu yüzden

$$\left(j-1,\ \mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right)$$

olur; burada \(\mathbf{e}_k\) histogramın birim vektörüdür. \(m=0\) durumu da aynı formülle kapsanır; tükenmiş eski mevcut rütbe sadece sıfır kovasına yazılır.

Adım 4: Beklenti bağıntısını kurma

\(E(m,\mathbf{a})\), \((m,\mathbf{a})\) durumundan itibaren beklenen ek çekiş sayısı olsun. Hiç kart kalmadıysa

$$E(m,\mathbf{a})=0\qquad\text{eğer }R(m,\mathbf{a})=0.$$

Aksi halde şimdi kesin olarak bir kart daha çekilecektir; baştaki \(1\) bunun içindir. O kart mevcut rütbe ile eşleşirse süreç hemen biter, bu yüzden başarı dalı için ayrıca özyineli bir terim yoktur. Sadece durmayan dalları toplarsak

$$E(m,\mathbf{a})=1+\sum_{j=1}^{c}\frac{j a_j}{R(m,\mathbf{a})}\,E\!\left(j-1,\mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right)$$

elde edilir. Bu, dinamik programın tamamıdır. Her geçiş desteden bir kart eksilttiği için özyineleme her zaman daha küçük durumlara gider.

Adım 5: Gerçek 52 kartlık destenin başlangıcı

İlk çekişten önce henüz mevcut bir rütbe yoktur; bu yüzden özyineleme ilk kart çekildikten sonra başlatılır. \(r\) rütbeli ve her rütbeden \(c\) kopyalı destede ilk çekiş sonrası

$$m=c-1,\qquad a_c=r-1,\qquad a_j=0\ \text{\(j\ne c\) için}$$

olur. Demek ki toplam beklenen çekiş sayısı

$$1+E\!\left(c-1,\mathbf{a}^{(0)}\right)$$

şeklindedir. Project Euler destesi için \(r=13\) ve \(c=4\) olduğundan

$$\mathbb{E}[\text{toplam çekiş}]=1+E(3,a_0=0,a_1=0,a_2=0,a_3=0,a_4=12)$$

yazılır. Bu başlangıç durumunun ilk adımı, şu basit kontrolü hemen açıklar:

$$\Pr(\text{2. çekişte durma})=\frac{3}{51}=\frac{1}{17}.$$

Çözümlü Örnek: iki rütbe ve her birinden iki kopya

Küçük bir deste bağıntıyı çok daha görünür kılar. \(r=2\) ve \(c=2\) alalım. İlk karttan sonra durum \(m=1\) ve \(a_2=1\) olur. Şimdi

$$X=E(1,a_0=0,a_1=0,a_2=1)$$

olsun. Destede \(R=3\) kart kalmıştır. Olasılık \(1/3\) ile sonraki kart mevcut rütbeden gelir ve süreç bir ek çekişten sonra biter. Olasılık \(2/3\) ile \(E(1,a_0=0,a_1=1,a_2=0)\) durumuna geçilir. Dolayısıyla

$$X=1+\frac{2}{3}E(1,a_0=0,a_1=1,a_2=0).$$

Şimdi

$$Y=E(1,a_0=0,a_1=1,a_2=0)$$

olsun. Bu durumda \(R=2\) olduğundan

$$Y=1+\frac{1}{2}E(0,a_0=0,a_1=1,a_2=0).$$

\(E(0,a_0=0,a_1=1,a_2=0)\) durumunda yalnızca bir kart kalmıştır. Bu kart tükenmiş mevcut rütbe ile eşleşemez ve çekildikten sonra deste boşalır. Bu yüzden

$$E(0,a_0=0,a_1=1,a_2=0)=1.$$

Böylece \(Y=1+\frac12=\frac32\), ardından \(X=1+\frac23\cdot\frac32=2\) bulunur. Toplam beklenti ise

$$1+X=3$$

olur. Bu değer, uygulamaların kullandığı küçük kontrol örneğiyle tam uyumludur.

Kod Nasıl Çalışır

C++, Python ve Java uygulamaları aynı memoization tabanlı özyinelemeyi kullanır. Mevcut rütbenin kalan kart sayısı ile diğer rütbelerin histogramı, her ulaşılabilir durumun yalnızca bir kez hesaplanacağı kompakt bir durum anahtarına dönüştürülür. Bir durum sorgulandığında uygulama önce \(R\)'yi hesaplar, deste boşsa \(0\) döndürür, değilse zorunlu olan bir sonraki çekişi sayarak başlar.

Daha sonra \(j=1,\dots,c\) kovaları dolaşılır. Kova \(j\) boş değilse uygulama \(\frac{j a_j}{R}\) olasılığını kullanır, yeni mevcut rütbeyi ve eski mevcut rütbenin histogram içine geri dönmesini yansıtan durumu kurar ve yeni durumun beklentisini özyineli olarak ister. Memoization sayesinde naif özyineleme ağacı, ulaşılabilir histogram durumları üzerinde küçük bir dinamik programa dönüşür.

Gerçek girdi \((r,c)=(13,4)\) için program sonucu sekiz ondalık basamakla yazar. Durum grafiği küçük olduğu için kayan nokta hesabı bu örnekte yeterlidir.

Karmaşıklık Analizi

\(r\) rütbeli ve her rütbeden \(c\) kopyalı genel bir destede, mevcut olmayan \(r-1\) rütbe en fazla

$$\binom{r+c-1}{c}$$

farklı şekilde \(c+1\) kovaya dağıtılabilir. Mevcut rütbe de \(m\in\{0,\dots,c-1\}\) olmak üzere \(c\) farklı değer alır. Bu yüzden soyut durum sayısı için üst sınır

$$|\mathcal{S}|\le c\binom{r+c-1}{c}$$

olur. Her durumda en fazla \(c\) geçiş incelendiğinden memoization'lı çalışma süresi \(O(c|\mathcal{S}|)\), bellek kullanımı ise \(O(|\mathcal{S}|)\) olur. Gerçek deste için bu üst sınır

$$4\binom{16}{4}=7280$$

durumdur; yani oldukça küçüktür.

Dipnotlar ve Kaynaklar

  1. Project Euler problem sayfası: https://projecteuler.net/problem=856
  2. Markov zinciri: Wikipedia — Markov chain
  3. Beklenen değer: Wikipedia — Expected value
  4. Hipergeometrik dağılım: Wikipedia — Hypergeometric distribution
  5. Dinamik programlama: Wikipedia — Dynamic programming

Resumen del Problema

Una baraja estándar tiene 13 rangos y 4 palos, así que cada rango aparece cuatro veces. Se extraen cartas sin reemplazo hasta que dos cartas consecutivas tengan el mismo rango o hasta que se agote la baraja. Lo que se busca es el valor esperado del número total de extracciones. Como solo importa la igualdad de rango, no hace falta seguir los palos por separado.

Enfoque Matemático

Conviene estudiar primero una baraja generalizada con \(r\) rangos y \(c\) copias por rango, y solo al final sustituir \(r=13\) y \(c=4\). La observación central es que, una vez extraída al menos una carta, el futuro depende únicamente de cuántas copias quedan del último rango visto y de cómo se distribuyen los demás rangos según su multiplicidad restante.

Paso 1: Comprimir la baraja en un histograma de rangos

Después de cualquier historia no terminal, llamemos rango actual al rango de la carta extraída más recientemente. Sea \(m\in\{0,1,\dots,c-1\}\) el número de cartas no repartidas de ese rango actual. Para cada \(j\in\{0,1,\dots,c\}\), sea \(a_j\) el número de otros rangos que tienen exactamente \(j\) cartas no repartidas.

Entonces \((m,a_0,\dots,a_c)\) es un estado de Markov completo. La identidad concreta de los rangos ya no importa; solo importan los conteos, porque los rangos con la misma multiplicidad restante son indistinguibles para todas las probabilidades futuras.

En el problema real \(c=4\), así que los cubos del histograma son \(a_0,a_1,a_2,a_3,a_4\), y los rangos no actuales siempre satisfacen

$$a_0+a_1+a_2+a_3+a_4=12.$$

Paso 2: Calcular las cartas restantes y la probabilidad de detenerse

Si el estado actual es \((m,\mathbf{a})\), el número de cartas no repartidas es

$$R(m,\mathbf{a})=m+\sum_{j=1}^{c} j a_j.$$

La siguiente extracción produce una pareja consecutiva exactamente cuando proviene del rango actual, del que sobreviven \(m\) copias. Por tanto

$$\Pr(\text{detenerse en la siguiente extracción}\mid m,\mathbf{a})=\frac{m}{R(m,\mathbf{a})}.$$

Si \(m=0\), el último rango visto ya está agotado y la siguiente carta no puede detener el proceso de inmediato. La misma fórmula sigue dando correctamente la probabilidad \(0\).

Paso 3: Describir las transiciones que no detienen el proceso

Supongamos que la siguiente carta viene de un rango distinto que antes de la extracción pertenecía al cubo \(j\), es decir, tenía \(j\) copias sin repartir. Existen \(a_j\) rangos de ese tipo y, en total, \(j a_j\) cartas elegibles. Por ello

$$\Pr(\text{elegir un rango del cubo }j\mid m,\mathbf{a})=\frac{j a_j}{R(m,\mathbf{a})},\qquad 1\le j\le c.$$

Tras esa extracción, el rango recién visto pasa a ser el nuevo rango actual y le quedan \(j-1\) cartas. El antiguo rango actual deja de ser actual y vuelve al histograma con \(m\) cartas restantes. Así, el nuevo estado es

$$\left(j-1,\ \mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right),$$

donde \(\mathbf{e}_k\) es el vector unidad del histograma. Esto también cubre el caso \(m=0\): el antiguo rango actual, ya agotado, simplemente entra en el cubo cero.

Paso 4: Escribir la recurrencia del valor esperado

Sea \(E(m,\mathbf{a})\) el número esperado de extracciones adicionales a partir del estado \((m,\mathbf{a})\). Si ya no quedan cartas, entonces

$$E(m,\mathbf{a})=0\qquad\text{cuando }R(m,\mathbf{a})=0.$$

En otro caso, ahora se extrae con certeza una carta más, y eso explica el término inicial \(1\). Si esa carta coincide con el rango actual, el proceso termina inmediatamente, así que esa rama no aporta término recursivo. Sumando solo las ramas que no terminan, obtenemos

$$E(m,\mathbf{a})=1+\sum_{j=1}^{c}\frac{j a_j}{R(m,\mathbf{a})}\,E\!\left(j-1,\mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right).$$

Ese es todo el programa dinámico. Cada transición elimina exactamente una carta no repartida, de modo que la recursión siempre avanza hacia estados más pequeños.

Paso 5: Inicializar la baraja real de 52 cartas

Antes de la primera extracción todavía no existe un rango actual, así que la recursión se inicia después de sacar la primera carta. En una baraja con \(r\) rangos y \(c\) copias por rango, esa primera extracción deja

$$m=c-1,\qquad a_c=r-1,\qquad a_j=0\ \text{para }j\ne c.$$

Por tanto, el valor esperado total es

$$1+E\!\left(c-1,\mathbf{a}^{(0)}\right).$$

Para la baraja del problema, \(r=13\) y \(c=4\), así que

$$\mathbb{E}[\text{extracciones totales}]=1+E(3,a_0=0,a_1=0,a_2=0,a_3=0,a_4=12).$$

La primera transición desde este estado ya explica la comprobación elemental

$$\Pr(\text{detenerse en la extracción 2})=\frac{3}{51}=\frac{1}{17}.$$

Ejemplo trabajado: dos rangos con dos copias cada uno

Una baraja muy pequeña vuelve transparente la recurrencia. Tomemos \(r=2\) y \(c=2\). Después de la primera carta, el estado es \(m=1\) y \(a_2=1\). Definamos

$$X=E(1,a_0=0,a_1=0,a_2=1).$$

Quedan \(R=3\) cartas. Con probabilidad \(1/3\), la siguiente carta coincide con el rango actual y el proceso termina tras una extracción adicional. Con probabilidad \(2/3\), pasamos al estado \(E(1,a_0=0,a_1=1,a_2=0)\). Luego

$$X=1+\frac{2}{3}E(1,a_0=0,a_1=1,a_2=0).$$

Sea ahora

$$Y=E(1,a_0=0,a_1=1,a_2=0).$$

Entonces \(R=2\), y por tanto

$$Y=1+\frac{1}{2}E(0,a_0=0,a_1=1,a_2=0).$$

Desde \(E(0,a_0=0,a_1=1,a_2=0)\) queda una sola carta. No puede coincidir con el rango actual agotado y, al extraerla, la baraja se vacía. Así que

$$E(0,a_0=0,a_1=1,a_2=0)=1.$$

De aquí sale \(Y=1+\frac12=\frac32\), luego \(X=1+\frac23\cdot\frac32=2\), y la esperanza total es

$$1+X=3.$$

Este valor coincide con la pequeña comprobación que usan las implementaciones.

Cómo Funciona el Código

Las implementaciones en C++, Python y Java usan la misma recursión con memoización. El resto del rango actual y el histograma de multiplicidades restantes se codifican en una clave compacta, de modo que cada estado alcanzable se evalúa una sola vez. Para un estado dado, la implementación calcula primero \(R\), devuelve \(0\) si la baraja ya está vacía y, en caso contrario, empieza contando la siguiente extracción obligatoria.

Después recorre los posibles cubos \(j=1,\dots,c\). Siempre que el cubo \(j\) no esté vacío, aplica la probabilidad \(\frac{j a_j}{R}\), actualiza el histograma para reflejar el nuevo rango actual y el retorno del antiguo rango actual a la estructura de cubos, y pide recursivamente la esperanza del nuevo estado. La memoización convierte el árbol recursivo ingenuo en un programa dinámico pequeño sobre los histogramas realmente alcanzables.

Para la entrada real \((r,c)=(13,4)\), el programa imprime la esperanza con ocho decimales. La aritmética en coma flotante es suficiente aquí porque el grafo de estados es pequeño y la recurrencia está bien condicionada.

Análisis de Complejidad

En una baraja general con \(r\) rangos y \(c\) copias por rango, los \(r-1\) rangos no actuales pueden repartirse entre los \(c+1\) cubos de a lo sumo

$$\binom{r+c-1}{c}$$

formas. El rango actual aporta además \(c\) posibles valores \(m\in\{0,\dots,c-1\}\), así que un límite superior para el número de estados abstractos es

$$|\mathcal{S}|\le c\binom{r+c-1}{c}.$$

Cada estado examina como máximo \(c\) transiciones, por lo que el tiempo total con memoización es \(O(c|\mathcal{S}|)\), y la memoria usada es \(O(|\mathcal{S}|)\). En la baraja real esto da, como máximo,

$$4\binom{16}{4}=7280$$

estados, una cantidad muy pequeña.

Notas y Referencias

  1. Página del problema: https://projecteuler.net/problem=856
  2. Cadena de Markov: Wikipedia — Markov chain
  3. Valor esperado: Wikipedia — Expected value
  4. Distribución hipergeométrica: Wikipedia — Hypergeometric distribution
  5. Programación dinámica: Wikipedia — Dynamic programming

问题概述

一副标准扑克牌有 13 个点数,每个点数有 4 张牌。现在从中不放回抽牌,直到出现“相邻两张牌点数相同”为止;如果始终没有出现这种情况,就一直抽到牌堆耗尽。要求计算总抽牌数的期望。由于题目只关心点数是否相同,并不关心花色,所以花色信息可以完全忽略。

数学方法

为了把结构讲清楚,先考虑一个推广模型:共有 \(r\) 个点数,每个点数各有 \(c\) 张牌。最后再把 \(r=13\)、\(c=4\) 代回原题。核心观察是:只要已经抽出至少一张牌,那么未来的概率分布只取决于“上一张牌所属点数还剩多少张”以及“其他点数按照剩余张数如何分组”,而不取决于更早之前抽牌的详细顺序。

步骤 1:把牌堆压缩成点数直方图状态

在任意一个尚未结束的历史之后,把最近一次抽到的点数称为当前点数。记 \(m\in\{0,1,\dots,c-1\}\) 为该当前点数在未抽牌中还剩多少张。对于每个 \(j\in\{0,1,\dots,c\}\),记 \(a_j\) 为“除当前点数以外,其余点数中恰好还剩 \(j\) 张牌的点数个数”。

这样,\((m,a_0,\dots,a_c)\) 就构成了一个完整的 Markov 状态。原因是:未来只需要知道各类点数还剩多少张,同一剩余张数的点数在未来概率上完全等价,所以不必记录它们的具体身份。

对原题而言,\(c=4\),所以状态中的直方图桶就是 \(a_0,a_1,a_2,a_3,a_4\)。除当前点数以外,一共有 12 个点数,因此始终满足

$$a_0+a_1+a_2+a_3+a_4=12.$$

步骤 2:计算剩余牌数与下一步直接停止的概率

在状态 \((m,\mathbf{a})\) 中,剩余未抽的牌总数为

$$R(m,\mathbf{a})=m+\sum_{j=1}^{c} j a_j.$$

下一次抽牌会立即形成相邻同点数,恰好当且仅当抽到的是当前点数中的某一张,而当前点数还剩 \(m\) 张。所以有

$$\Pr(\text{下一抽即停止}\mid m,\mathbf{a})=\frac{m}{R(m,\mathbf{a})}.$$

如果 \(m=0\),说明当前点数已经被抽空了,那么下一张牌不可能立刻和上一张牌构成同点数相邻对。上式在这种情况下自然给出概率 \(0\)。

步骤 3:写出不会立即停止时的状态转移

假设下一张牌不是当前点数,而是来自某个“抽牌前还有 \(j\) 张剩余牌”的其他点数。这样的点数共有 \(a_j\) 个,对应的可抽牌张数一共是 \(j a_j\),因此

$$\Pr(\text{从桶 }j\text{ 选中一个点数}\mid m,\mathbf{a})=\frac{j a_j}{R(m,\mathbf{a})},\qquad 1\le j\le c.$$

抽完这张牌后,新抽到的那个点数就变成新的当前点数,并且它的剩余张数从 \(j\) 变成 \(j-1\)。原来的当前点数不再是当前点数,它会以“还剩 \(m\) 张”的身份重新并入直方图。因此新状态为

$$\left(j-1,\ \mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right),$$

其中 \(\mathbf{e}_k\) 表示直方图中的第 \(k\) 个单位向量。这个写法也自动覆盖了 \(m=0\) 的情形:旧的当前点数已经耗尽,于是它会进入零张剩余的桶。

步骤 4:建立期望的递推方程

令 \(E(m,\mathbf{a})\) 表示从状态 \((m,\mathbf{a})\) 出发,未来还要再抽多少张牌的期望值。如果已经没有牌了,那么

$$E(m,\mathbf{a})=0\qquad\text{当 }R(m,\mathbf{a})=0.$$

否则,现在必然还要再抽一张牌,所以方程前面一定有一个 \(1\)。如果这张牌恰好来自当前点数,那么过程就在这一抽后立刻终止,因此“成功分支”不会再带来递推项。于是只需对那些没有立刻终止的分支求和:

$$E(m,\mathbf{a})=1+\sum_{j=1}^{c}\frac{j a_j}{R(m,\mathbf{a})}\,E\!\left(j-1,\mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right).$$

这就是完整的动态规划方程。每一步转移都会让剩余牌数减少 1,因此递归一定会朝着更小的状态前进,不会出现循环依赖。

步骤 5:代入原题的初始状态

在第一张牌抽出之前,还不存在“当前点数”,所以递推自然从“已经抽完第一张牌之后”的状态开始。对于一般的 \(r\) 点数、每点数 \(c\) 张的牌堆,第一张牌抽出后有

$$m=c-1,\qquad a_c=r-1,\qquad a_j=0\ \text{当 }j\ne c.$$

因此总抽牌数的期望就是

$$1+E\!\left(c-1,\mathbf{a}^{(0)}\right).$$

把原题参数 \(r=13\)、\(c=4\) 代入,得到

$$\mathbb{E}[\text{总抽牌数}]=1+E(3,a_0=0,a_1=0,a_2=0,a_3=0,a_4=12).$$

从这个初始状态出发,第二张牌立刻让过程终止的概率就是

$$\Pr(\text{第 2 张牌时停止})=\frac{3}{51}=\frac{1}{17},$$

这也是一个很直接的正确性检查。

例子:两个点数、每个点数两张牌

用一个很小的模型可以把递推算得一清二楚。取 \(r=2\)、\(c=2\)。抽出第一张牌后,当前点数还剩 1 张,另一个点数还剩 2 张,因此状态为 \(m=1\)、\(a_2=1\)。记

$$X=E(1,a_0=0,a_1=0,a_2=1).$$

此时剩余牌数 \(R=3\)。下一张牌以概率 \(1/3\) 来自当前点数,于是再抽 1 张就结束;以概率 \(2/3\) 来自另一个点数,转入状态 \(E(1,a_0=0,a_1=1,a_2=0)\)。所以

$$X=1+\frac{2}{3}E(1,a_0=0,a_1=1,a_2=0).$$

再记

$$Y=E(1,a_0=0,a_1=1,a_2=0).$$

此时 \(R=2\),于是

$$Y=1+\frac{1}{2}E(0,a_0=0,a_1=1,a_2=0).$$

而在状态 \(E(0,a_0=0,a_1=1,a_2=0)\) 中,只剩最后一张牌。它不可能与已经耗尽的当前点数匹配,抽完后牌堆清空,因此

$$E(0,a_0=0,a_1=1,a_2=0)=1.$$

从而得到 \(Y=1+\frac12=\frac32\),再代回得

$$X=1+\frac23\cdot\frac32=2.$$

最后总期望为

$$1+X=3.$$

这个小例子与实现中的检验值完全一致,也清楚说明了递推方程为什么正确。

代码如何实现

C++、Python 和 Java 实现采用的是同一套“记忆化递归”思路。它们把当前点数剩余张数与其他点数的剩余张数直方图编码成一个紧凑的状态键,这样每个可达状态只会被计算一次。对于某个状态,程序先计算 \(R\);若牌堆已经为空,就返回 \(0\);否则先把“接下来必抽的一张牌”记入答案。

随后程序遍历所有可能的桶 \(j=1,\dots,c\)。只要第 \(j\) 个桶非空,就按概率 \(\frac{j a_j}{R}\) 进入相应后继状态,并把新当前点数与旧当前点数回到直方图中的变化准确更新。记忆化把原本会指数膨胀的递归树压缩成一个很小的状态图上的动态规划。

对于原题的 \((r,c)=(13,4)\),程序最后输出八位小数。由于状态数很少,而且递推结构稳定,这里使用浮点数就足够了。

复杂度分析

对于一般的 \(r\) 点数、每点数 \(c\) 张的牌堆,除当前点数外还有 \(r-1\) 个点数,它们被分配到 \(c+1\) 个桶中的方式最多有

$$\binom{r+c-1}{c}$$

种。当前点数的剩余量 \(m\) 还可以取 \(0,1,\dots,c-1\) 共 \(c\) 个值,所以抽象状态总数满足上界

$$|\mathcal{S}|\le c\binom{r+c-1}{c}.$$

每个状态最多考察 \(c\) 条转移,因此记忆化后的总时间复杂度为 \(O(c|\mathcal{S}|)\),空间复杂度为 \(O(|\mathcal{S}|)\)。对于真实牌堆,这个上界是

$$4\binom{16}{4}=7280,$$

规模非常小。

注释与参考资料

  1. Project Euler 题目页面:https://projecteuler.net/problem=856
  2. Markov 链:Wikipedia — Markov chain
  3. 期望值:Wikipedia — Expected value
  4. 超几何分布:Wikipedia — Hypergeometric distribution
  5. 动态规划:Wikipedia — Dynamic programming

Краткое описание задачи

В стандартной колоде 13 рангов, и каждый ранг встречается 4 раза. Карты тянутся без возвращения до тех пор, пока либо две соседние вытянутые карты не окажутся одного ранга, либо колода не закончится. Требуется найти математическое ожидание общего числа вытянутых карт. Поскольку важна только равенство рангов, масти можно вообще не отслеживать.

Математический подход

Удобно сначала разобрать обобщенную колоду с \(r\) рангами и \(c\) картами каждого ранга, а затем подставить \(r=13\) и \(c=4\). Главная идея состоит в том, что после того как уже вытянута хотя бы одна карта, будущее определяется только количеством оставшихся карт последнего ранга и тем, как остальные ранги распределены по числу оставшихся копий.

Шаг 1: Сжимаем колоду до гистограммы рангов

После любой незавершенной истории назовем ранг последней вытянутой карты текущим рангом. Пусть \(m\in\{0,1,\dots,c-1\}\) — число еще не вытянутых карт этого текущего ранга. Для каждого \(j\in\{0,1,\dots,c\}\) пусть \(a_j\) обозначает количество всех остальных рангов, у которых осталось ровно \(j\) карт.

Тогда \((m,a_0,\dots,a_c)\) является полным марковским состоянием. Конкретные названия рангов уже не нужны; важны только количества, потому что ранги с одинаковой остаточной кратностью полностью эквивалентны для дальнейших вероятностей.

В исходной задаче \(c=4\), значит, используются корзины \(a_0,a_1,a_2,a_3,a_4\), а для 12 нетекущих рангов всегда выполнено

$$a_0+a_1+a_2+a_3+a_4=12.$$

Шаг 2: Число оставшихся карт и вероятность немедленной остановки

В состоянии \((m,\mathbf{a})\) число еще не вытянутых карт равно

$$R(m,\mathbf{a})=m+\sum_{j=1}^{c} j a_j.$$

Следующая вытяжка немедленно создает соседнюю пару ровно тогда, когда приходит карта текущего ранга, а таких карт осталось \(m\). Поэтому

$$\Pr(\text{остановка на следующем ходе}\mid m,\mathbf{a})=\frac{m}{R(m,\mathbf{a})}.$$

Если \(m=0\), то последний увиденный ранг уже исчерпан, и следующая карта не может завершить процесс сразу. Та же формула корректно дает вероятность \(0\).

Шаг 3: Переходы, которые не завершают процесс

Предположим, что следующая карта приходит не из текущего ранга, а из некоторого другого ранга, который до вытяжки лежал в корзине \(j\), то есть имел \(j\) оставшихся карт. Таких рангов \(a_j\), а подходящих карт всего \(j a_j\). Значит,

$$\Pr(\text{выбрать ранг из корзины }j\mid m,\mathbf{a})=\frac{j a_j}{R(m,\mathbf{a})},\qquad 1\le j\le c.$$

После этой вытяжки новый ранг становится текущим и у него остается \(j-1\) карт. Старый текущий ранг перестает быть текущим и возвращается в гистограмму с \(m\) оставшимися картами. Следовательно, новое состояние имеет вид

$$\left(j-1,\ \mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right),$$

где \(\mathbf{e}_k\) — единичный вектор гистограммы. Это автоматически покрывает и случай \(m=0\): исчерпанный прежний текущий ранг просто попадает в нулевую корзину.

Шаг 4: Рекуррентное уравнение для ожидания

Обозначим через \(E(m,\mathbf{a})\) ожидаемое число дополнительных вытяжек из состояния \((m,\mathbf{a})\). Если карт уже не осталось, то

$$E(m,\mathbf{a})=0\qquad\text{при }R(m,\mathbf{a})=0.$$

Иначе сейчас обязательно произойдет еще одна вытяжка, что и дает начальный член \(1\). Если эта карта совпала с текущим рангом, процесс сразу заканчивается, поэтому для успешной ветви дополнительного рекурсивного слагаемого нет. Остается суммировать только по неостанавливающим ветвям:

$$E(m,\mathbf{a})=1+\sum_{j=1}^{c}\frac{j a_j}{R(m,\mathbf{a})}\,E\!\left(j-1,\mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right).$$

Это и есть все динамическое программирование. Каждый переход уменьшает число оставшихся карт на единицу, поэтому рекурсия всегда идет к меньшим состояниям.

Шаг 5: Начальное состояние реальной колоды

До первой вытяжки текущего ранга еще нет, поэтому рекурсию естественно запускать после того, как первая карта уже вынута. В колоде с \(r\) рангами и \(c\) копиями каждого ранга после первой карты остается

$$m=c-1,\qquad a_c=r-1,\qquad a_j=0\ \text{для }j\ne c.$$

Следовательно, искомое полное ожидание равно

$$1+E\!\left(c-1,\mathbf{a}^{(0)}\right).$$

Для колоды Project Euler, где \(r=13\) и \(c=4\), получаем

$$\mathbb{E}[\text{общее число вытяжек}]=1+E(3,a_0=0,a_1=0,a_2=0,a_3=0,a_4=12).$$

Уже первый шаг из этого состояния дает простую проверку

$$\Pr(\text{остановка на 2-й карте})=\frac{3}{51}=\frac{1}{17}.$$

Разобранный пример: два ранга по две карты

На маленькой колоде рекурсия видна особенно ясно. Возьмем \(r=2\) и \(c=2\). После первой карты состояние равно \(m=1\) и \(a_2=1\). Обозначим

$$X=E(1,a_0=0,a_1=0,a_2=1).$$

Осталось \(R=3\) карты. С вероятностью \(1/3\) следующая карта совпадает с текущим рангом, и процесс заканчивается после еще одной вытяжки. С вероятностью \(2/3\) мы переходим в состояние \(E(1,a_0=0,a_1=1,a_2=0)\). Поэтому

$$X=1+\frac{2}{3}E(1,a_0=0,a_1=1,a_2=0).$$

Теперь положим

$$Y=E(1,a_0=0,a_1=1,a_2=0).$$

Тогда \(R=2\), и значит

$$Y=1+\frac{1}{2}E(0,a_0=0,a_1=1,a_2=0).$$

В состоянии \(E(0,a_0=0,a_1=1,a_2=0)\) остается только одна карта. Она не может совпасть с уже исчерпанным текущим рангом, а после ее вытяжки колода пустеет, поэтому

$$E(0,a_0=0,a_1=1,a_2=0)=1.$$

Отсюда \(Y=1+\frac12=\frac32\), затем \(X=1+\frac23\cdot\frac32=2\), а полное ожидание равно

$$1+X=3.$$

Именно это маленькое контрольное значение подтверждают реализации.

Как Работает Код

Реализации на C++, Python и Java используют одну и ту же идею: рекурсию с мемоизацией. Остаток текущего ранга и гистограмма кратностей остальных рангов кодируются в компактный ключ состояния, поэтому каждый достижимый узел вычисляется только один раз. Для данного состояния программа сначала считает \(R\), возвращает \(0\), если колода уже пуста, и иначе добавляет обязательную следующую вытяжку.

Затем перебираются корзины \(j=1,\dots,c\). Если корзина \(j\) непуста, программа использует вероятность \(\frac{j a_j}{R}\), обновляет описание состояния так, чтобы новый текущий ранг и старый текущий ранг были учтены правильно, и рекурсивно получает ожидание дочернего состояния. Мемоизация превращает наивное дерево рекурсии в компактную динамику по реально достижимым гистограммам.

Для входа \((r,c)=(13,4)\) программа печатает ответ с точностью до восьми знаков после запятой. Арифметики с плавающей точкой здесь достаточно, потому что граф состояний очень мал.

Анализ сложности

В общей колоде с \(r\) рангами и \(c\) копиями каждого ранга \(r-1\) нетекущих рангов можно разложить по \(c+1\) корзинам не более чем

$$\binom{r+c-1}{c}$$

способами. Текущий ранг добавляет еще \(c\) возможных значений \(m\in\{0,\dots,c-1\}\), так что число абстрактных состояний ограничено сверху величиной

$$|\mathcal{S}|\le c\binom{r+c-1}{c}.$$

Каждое состояние рассматривает не более \(c\) переходов, поэтому общее время работы с мемоизацией равно \(O(c|\mathcal{S}|)\), а память — \(O(|\mathcal{S}|)\). Для реальной колоды это дает максимум

$$4\binom{16}{4}=7280$$

состояний, то есть совсем немного.

Примечания и Ссылки

  1. Страница задачи Project Euler: https://projecteuler.net/problem=856
  2. Цепь Маркова: Wikipedia — Markov chain
  3. Математическое ожидание: Wikipedia — Expected value
  4. Гипергеометрическое распределение: Wikipedia — Hypergeometric distribution
  5. Динамическое программирование: Wikipedia — Dynamic programming

ملخص المسألة

تتكون الرزمة القياسية من 13 رتبة، ولكل رتبة 4 بطاقات. نسحب البطاقات دون إرجاع حتى تظهر بطاقتان متتاليتان من الرتبة نفسها، أو حتى تنفد الرزمة تمامًا. المطلوب هو القيمة المتوقعة لعدد البطاقات المسحوبة كليًا. وبما أن الشرط يعتمد على الرتبة فقط، فلا حاجة إلى تتبع الأصناف منفردة.

المنهج الرياضي

من الأنسب أولًا دراسة رزمة معممة تحتوي على \(r\) رتب و\(c\) نسخ من كل رتبة، ثم نعوض في النهاية بـ \(r=13\) و\(c=4\). الفكرة الأساسية هي أن المستقبل، بعد سحب بطاقة واحدة على الأقل، لا يعتمد على التسلسل الكامل للسحوبات السابقة، بل يعتمد فقط على عدد البطاقات المتبقية من رتبة آخر بطاقة سُحبت، وعلى كيفية توزع الرتب الأخرى بحسب عدد النسخ المتبقية منها.

الخطوة 1: ضغط حالة الرزمة إلى هيستوغرام للرتب

بعد أي تاريخ لم ينته بعد، نسمي رتبة آخر بطاقة سُحبت الرتبة الحالية. ولتكن \(m\in\{0,1,\dots,c-1\}\) هي عدد البطاقات غير المسحوبة المتبقية من هذه الرتبة الحالية. ولكل \(j\in\{0,1,\dots,c\}\)، ليكن \(a_j\) عدد الرتب الأخرى التي بقي منها بالضبط \(j\) بطاقات.

إذًا تصبح \((m,a_0,\dots,a_c)\) حالة ماركوف كاملة. هوية الرتب نفسها لا تهم؛ ما يهم فقط هو الأعداد، لأن أي رتبتين لهما العدد نفسه من النسخ المتبقية تكونان متماثلتين من حيث جميع الاحتمالات المستقبلية.

في المسألة الأصلية لدينا \(c=4\)، ولذلك تكون سلال الهيستوغرام هي \(a_0,a_1,a_2,a_3,a_4\)، وتبقى الرتب غير الحالية دومًا محققة للعلاقة

$$a_0+a_1+a_2+a_3+a_4=12.$$

الخطوة 2: حساب عدد البطاقات المتبقية واحتمال التوقف في السحبة التالية

إذا كانت الحالة الحالية هي \((m,\mathbf{a})\)، فإن عدد البطاقات غير المسحوبة يساوي

$$R(m,\mathbf{a})=m+\sum_{j=1}^{c} j a_j.$$

تتكون الزوجية المتجاورة في السحبة التالية فقط إذا كانت البطاقة التالية من الرتبة الحالية نفسها، وهذه الرتبة ما زال منها \(m\) بطاقات. لذلك

$$\Pr(\text{stop next}\mid m,\mathbf{a})=\frac{m}{R(m,\mathbf{a})}.$$

إذا كان \(m=0\)، فهذا يعني أن الرتبة الحالية قد استنفدت بالفعل، ومن ثم لا يمكن للسحبة التالية أن توقف العملية مباشرة. والصيغة نفسها تعطي عندئذ القيمة الصحيحة وهي \(0\).

الخطوة 3: وصف الانتقالات التي لا توقف العملية

افترض أن البطاقة التالية لا تأتي من الرتبة الحالية، بل من رتبة أخرى كانت قبل السحب موجودة في السلة \(j\)، أي كان بقي منها \(j\) بطاقات. يوجد \(a_j\) رتب من هذا النوع، وبالتالي يوجد \(j a_j\) بطاقات مؤهلة. ومن ثم

$$\Pr(\text{choose bucket }j\mid m,\mathbf{a})=\frac{j a_j}{R(m,\mathbf{a})},\qquad 1\le j\le c.$$

بعد هذه السحبة تصبح الرتبة الجديدة هي الرتبة الحالية، ويتبقى منها \(j-1\) بطاقات. أما الرتبة الحالية القديمة فتتوقف عن كونها حالية وتعود إلى الهيستوغرام ومعها \(m\) بطاقات متبقية. لذلك تكون الحالة الجديدة

$$\left(j-1,\ \mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right),$$

حيث \(\mathbf{e}_k\) هو متجه الوحدة في الهيستوغرام. وهذه الصيغة تشمل أيضًا الحالة \(m=0\): فالرتبة الحالية القديمة المستنفدة تُسجَّل ببساطة في سلة الصفر.

الخطوة 4: كتابة علاقة التوقع

لنرمز بـ \(E(m,\mathbf{a})\) إلى العدد المتوقع من السحوبات الإضافية انطلاقًا من الحالة \((m,\mathbf{a})\). إذا لم تبق أي بطاقة، فلدينا

$$E(m,\mathbf{a})=0\qquad\text{if }R(m,\mathbf{a})=0.$$

وإلا فإن بطاقة واحدة أخرى ستُسحب يقينًا الآن، وهذا ما يفسر الحد الابتدائي \(1\). وإذا طابقت هذه البطاقة الرتبة الحالية فإن العملية تتوقف فورًا، لذلك لا يوجد حد تكراري لذلك الفرع الناجح. وبجمع الفروع التي لا تتوقف مباشرة نحصل على

$$E(m,\mathbf{a})=1+\sum_{j=1}^{c}\frac{j a_j}{R(m,\mathbf{a})}\,E\!\left(j-1,\mathbf{a}-\mathbf{e}_j+\mathbf{e}_m\right).$$

هذه هي صيغة البرمجة الديناميكية كاملة. كل انتقال يزيل بطاقة واحدة من الرزمة المتبقية، ولذا فإن الاستدعاء التكراري يتجه دائمًا نحو حالات أصغر.

الخطوة 5: تهيئة الرزمة الحقيقية ذات 52 بطاقة

قبل السحبة الأولى لا توجد رتبة حالية بعد، لذلك تبدأ العلاقة بعد سحب البطاقة الأولى. في رزمة فيها \(r\) رتب و\(c\) نسخ من كل رتبة، تترك السحبة الأولى الحالة

$$m=c-1,\qquad a_c=r-1,\qquad a_j=0\ \text{for }j\ne c.$$

وبالتالي فإن التوقع الكلي يساوي

$$1+E\!\left(c-1,\mathbf{a}^{(0)}\right).$$

وبالتعويض في مسألة Project Euler بـ \(r=13\) و\(c=4\) نحصل على

$$\mathbb{E}[\text{total draws}]=1+E(3,a_0=0,a_1=0,a_2=0,a_3=0,a_4=12).$$

ومن هذه الحالة الابتدائية نحصل فورًا على فحص بسيط للصحة:

$$\Pr(\text{stop on draw 2})=\frac{3}{51}=\frac{1}{17}.$$

مثال محلول: رتبتان ولكل رتبة بطاقتان

الرزمة الصغيرة تجعل العلاقة واضحة جدًا. خذ \(r=2\) و\(c=2\). بعد سحب البطاقة الأولى تصبح الحالة \(m=1\) و\(a_2=1\). لنكتب

$$X=E(1,a_0=0,a_1=0,a_2=1).$$

يبقى \(R=3\) بطاقات. باحتمال \(1/3\) تأتي البطاقة التالية من الرتبة الحالية نفسها، فتتوقف العملية بعد سحبة إضافية واحدة. وباحتمال \(2/3\) ننتقل إلى الحالة \(E(1,a_0=0,a_1=1,a_2=0)\). إذن

$$X=1+\frac{2}{3}E(1,a_0=0,a_1=1,a_2=0).$$

ولنضع الآن

$$Y=E(1,a_0=0,a_1=1,a_2=0).$$

حينها يكون \(R=2\)، ولذلك

$$Y=1+\frac{1}{2}E(0,a_0=0,a_1=1,a_2=0).$$

أما الحالة \(E(0,a_0=0,a_1=1,a_2=0)\) ففيها تبقى بطاقة واحدة فقط. لا يمكن لهذه البطاقة أن تطابق الرتبة الحالية المستنفدة، وبعد سحبها تصبح الرزمة فارغة، ولذلك

$$E(0,a_0=0,a_1=1,a_2=0)=1.$$

ومن ثم \(Y=1+\frac12=\frac32\)، وبعدها

$$X=1+\frac23\cdot\frac32=2.$$

إذن التوقع الكلي هو

$$1+X=3.$$

وهذا يطابق مثال التحقق الصغير الذي تستخدمه التطبيقات.

كيف تعمل الشيفرة

تستخدم تطبيقات C++ وPython وJava الفكرة نفسها: علاقة عودية مع حفظ النتائج المحسوبة مسبقًا. يتم ترميز عدد البطاقات المتبقية من الرتبة الحالية مع هيستوغرام الرتب الأخرى في مفتاح حالة مضغوط، وبذلك تُحسب كل حالة قابلة للوصول مرة واحدة فقط. عند تقييم حالة ما، تحسب الشيفرة أولًا \(R\)، وتعيد \(0\) إذا كانت الرزمة قد فرغت، وإلا تبدأ من السحبة التالية الإلزامية.

بعد ذلك تمر الشيفرة على السلال الممكنة \(j=1,\dots,c\). فإذا كانت السلة \(j\) غير فارغة، استُخدم الاحتمال \(\frac{j a_j}{R}\)، ثم حُدثت الحالة بحيث تمثل الرتبة الحالية الجديدة وعودة الرتبة الحالية القديمة إلى الهيستوغرام، ثم استُدعي التوقع للحالة التالية. هكذا يتحول شجر الاستدعاء التكراري الخام إلى برنامج ديناميكي صغير على فضاء الحالات القابلة للوصول فعلًا.

في الإدخال الحقيقي \((r,c)=(13,4)\) تطبع البرامج النتيجة حتى ثماني منازل عشرية. والحساب العشري كاف هنا لأن مخطط الحالات صغير جدًا.

تحليل التعقيد

في الرزمة العامة ذات \(r\) رتب و\(c\) نسخ لكل رتبة، يمكن توزيع الرتب \(r-1\) غير الحالية على \(c+1\) سلال بحد أقصى قدره

$$\binom{r+c-1}{c}$$

من الأشكال. وتضيف الرتبة الحالية \(c\) قيم ممكنة لـ \(m\in\{0,\dots,c-1\}\)، ولذلك نحصل على الحد الأعلى التالي لعدد الحالات المجردة:

$$|\mathcal{S}|\le c\binom{r+c-1}{c}.$$

وكل حالة تنظر في ما لا يزيد على \(c\) انتقالات، لذا فإن زمن التنفيذ مع الحفظ يساوي \(O(c|\mathcal{S}|)\)، أما الذاكرة فتساوي \(O(|\mathcal{S}|)\). وفي الرزمة الحقيقية نحصل على حد أعلى مقداره

$$4\binom{16}{4}=7280$$

حالة فقط، وهو عدد صغير جدًا.

هوامش ومراجع

  1. صفحة مسألة Project Euler: https://projecteuler.net/problem=856
  2. سلسلة ماركوف: Wikipedia — Markov chain
  3. القيمة المتوقعة: Wikipedia — Expected value
  4. التوزيع فوق الهندسي: Wikipedia — Hypergeometric distribution
  5. البرمجة الديناميكية: Wikipedia — Dynamic programming