Problem Summary

Consider the lattice square \(\{0,\dots,N\}^2\), which contains \((N+1)^2\) grid points. Each point is selected independently with probability

$$p=\frac{1}{N+1},\qquad q=1-p.$$

Let \(H\) be the convex hull of the selected points. The goal is to compute \(\mathbb{E}[\operatorname{Area}(H)]\). A direct enumeration of all \(2^{(N+1)^2}\) subsets is impossible, so the solution rewrites the expected area as a sum over possible hull edges and then evaluates that sum direction by direction.

Mathematical Approach

The central idea is that area can be decomposed into oriented edge contributions. Once that is done, linearity of expectation allows us to analyze one candidate hull edge at a time and sum the results.

Step 1: Rewrite area as a sum over oriented hull edges

If the hull vertices in counterclockwise order are \(v_0,\dots,v_{m-1}\), then the shoelace formula gives

$$2\,\operatorname{Area}(H)=\sum_{t=0}^{m-1} v_t\times v_{t+1},\qquad v_t\times v_{t+1}=x_t y_{t+1}-y_t x_{t+1}.$$

Taking expectation and using linearity, we may sum over all ordered lattice-point pairs \((u,v)\):

$$2\,\mathbb{E}[\operatorname{Area}(H)]=\sum_{u,v}(u\times v)\,\mathbb{P}(u\to v\text{ appears as a counterclockwise hull edge}).$$

So the problem becomes a probability question: for a fixed oriented segment \(u\to v\), when does it contribute as an actual edge of the convex hull?

Step 2: Group candidate edges by primitive direction and supporting line

Any lattice segment has direction proportional to a primitive vector \(d=(dx,dy)\) with \(\gcd(|dx|,|dy|)=1\). Fix such a direction. All lattice points on a line parallel to \(d\) satisfy

$$s(x,y)=dy\,x-dx\,y=\text{constant}.$$

Therefore the implementation groups candidate edges by primitive direction and by the value of \(s\), which identifies one supporting line. On a fixed line, ordering the lattice points by repeated steps of \(d\) gives

$$P_0,P_1,\dots,P_{k-1}.$$

Every possible hull edge with that slope is some pair \(P_i\to P_j\) with \(0\le i\lt j\lt k\).

Step 3: Compute the probability that one pair is the hull edge

Fix one supporting line and one pair \(P_i\to P_j\). For this segment to be the hull edge whose interior lies on the left, four independent conditions must hold:

$$\text{(a) }P_i\text{ and }P_j\text{ are selected},$$

$$\text{(b) every point strictly to the right of the line is unselected},$$

$$\text{(c) every collinear point before }P_i\text{ or after }P_j\text{ is unselected},$$

$$\text{(d) at least one point strictly to the left of the line is selected}.$$

If \(L\) is the number of lattice points strictly on the left, \(R\) the number strictly on the right, and

$$O=i+(k-1-j)$$

the number of collinear points outside the segment, then independence gives

$$\mathbb{P}(P_i\to P_j\text{ is that hull edge})=p^2 q^{R} q^{O}(1-q^{L}).$$

Points lying between \(P_i\) and \(P_j\) on the same line are allowed to be selected; they remain on the boundary segment and do not change the area contribution. The factor \(1-q^L\) also removes the degenerate case in which all selected points are collinear and the area is zero.

Step 4: Sum all pairs on one supporting line

For a fixed line, \(L\) and \(R\) depend only on the line, while \(O=i+(k-1-j)\) depends on the chosen pair. Hence the total expected oriented contribution from one line is

$$p^2 q^{R}(1-q^{L})\sum_{0\le i\lt j\lt k}(P_i\times P_j)\,q^{i}q^{k-1-j}.$$

This is exactly why the implementation loops over all ordered pairs on a line: the cross product \(P_i\times P_j\) supplies the geometric term, and the powers of \(q\) encode the probability that these two points are the extreme selected points on that line.

Step 5: Sum over all primitive directions

Now sum the previous expression over every supporting line of every primitive direction. Each genuine convex-hull edge is counted exactly once: its slope determines the primitive direction, and the counterclockwise orientation is the unique orientation for which the hull lies on the left and the exterior lies on the right.

Therefore

$$\boxed{\mathbb{E}[\operatorname{Area}(H)]=\frac12\sum_{\text{primitive }d}\ \sum_{\ell\parallel d} p^2 q^{R(\ell)}(1-q^{L(\ell)})\sum_{0\le i\lt j\lt k_\ell}(P_i\times P_j)\,q^{i}q^{k_\ell-1-j}.}$$

This is the formula evaluated by the program.

Worked Example: \(N=1\)

For \(N=1\), the grid is the unit square with four corner points and

$$p=\frac12.$$

The hull has area \(1\) only when all four corners are selected, which happens with probability \(1/16\). It has area \(1/2\) when exactly three of the four corners are selected, and there are four such configurations, each with probability \(1/16\). All remaining configurations have area \(0\).

Hence

$$\mathbb{E}[\operatorname{Area}(H)]=1\cdot\frac{1}{16}+\frac12\cdot 4\cdot\frac{1}{16}=\frac{3}{16}=0.18750,$$

which matches the checkpoint used by the implementation.

How the Code Works

The C++, Python, and Java implementations all follow the same pipeline. They first compute \(p\), \(q\), and a table of powers \(q^0,q^1,\dots,q^{(N+1)^2}\), because every probability factor is assembled from these values. Next they enumerate every primitive direction \((dx,dy)\) inside the square.

For one direction, the implementation scans all lattice points and evaluates \(s(x,y)=dy\,x-dx\,y\). This simultaneously gives the population of each parallel line and identifies the boundary point from which that line should be traversed so that every line is visited exactly once. Prefix sums over the line populations then give, for each line, how many points lie strictly to its left and strictly to its right.

After that, each line is walked in order, its points are stored as \(P_0,\dots,P_{k-1}\), and every pair \(P_i,P_j\) contributes

$$(P_i\times P_j)\,q^i q^{k-1-j}$$

to the line sum. Multiplying by the line factor

$$p^2 q^R(1-q^L)$$

turns that geometric sum into an expected oriented-area contribution. Adding all directions and dividing by \(2\) yields the final expected area. Implementations with native parallel support split the direction set across workers, but the mathematical result is identical.

Complexity Analysis

The number of primitive directions with \(|dx|,|dy|\le N\) is \(O(N^2)\). For a fixed direction, building line statistics over the \((N+1)^2\) lattice points costs \(O(N^2)\), while the pair accumulation on one line of length \(k\) costs \(O(k^2)\). Summed over all lines of that direction, \(\sum k=(N+1)^2\) and the worst case is \(\sum k^2=O(N^3)\), so one direction costs \(O(N^3)\) time in the worst case.

Consequently the total work is \(O(N^5)\). The memory usage is \(O(N^2)\): the power table, the per-direction line counts, and the temporary storage for a single line are all quadratic or smaller. Parallel execution reduces wall-clock time but does not change the asymptotic bounds.

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=514
  2. The oriented-area identity for polygons, often called the shoelace formula.
  3. Convex hulls, supporting lines, and the fact that a hull edge is determined by one empty half-plane.
  4. Primitive lattice directions, i.e. coprime integer step vectors, which enumerate all lattice slopes without duplication.

Problemzusammenfassung

Betrachtet wird das Gitterquadrat \(\{0,\dots,N\}^2\) mit \((N+1)^2\) Gitterpunkten. Jeder Punkt wird unabhängig mit der Wahrscheinlichkeit

$$p=\frac{1}{N+1},\qquad q=1-p$$

ausgewählt. Sei \(H\) die konvexe Hülle der gewählten Punkte. Gesucht ist \(\mathbb{E}[\operatorname{Area}(H)]\). Die direkte Betrachtung aller \(2^{(N+1)^2}\) Punktmengen ist unmöglich, daher wird der Erwartungswert der Fläche als Summe über mögliche Hüllkanten umgeschrieben und dann richtungsweise ausgewertet.

Mathematischer Ansatz

Der zentrale Gedanke ist, die Fläche in Beiträge orientierter Kanten zu zerlegen. Danach erlaubt die Linearität des Erwartungswerts, jede mögliche Hüllkante einzeln zu analysieren und alle Beiträge zu addieren.

Schritt 1: Fläche als Summe orientierter Hüllkanten schreiben

Seien \(v_0,\dots,v_{m-1}\) die Hüllpunkte in gegen den Uhrzeigersinn laufender Reihenfolge. Dann liefert die Gaußsche Trapezformel

$$2\,\operatorname{Area}(H)=\sum_{t=0}^{m-1} v_t\times v_{t+1},\qquad v_t\times v_{t+1}=x_t y_{t+1}-y_t x_{t+1}.$$

Nimmt man den Erwartungswert und benutzt Linearität, darf man über alle geordneten Punktpaare \((u,v)\) summieren:

$$2\,\mathbb{E}[\operatorname{Area}(H)]=\sum_{u,v}(u\times v)\,\mathbb{P}(u\to v\text{ erscheint als Hüllkante gegen den Uhrzeigersinn}).$$

Damit wird aus dem Flächenproblem ein Wahrscheinlichkeitsproblem: Für wann ist ein festes orientiertes Segment \(u\to v\) tatsächlich eine Kante der konvexen Hülle?

Schritt 2: Nach primitiver Richtung und Stützgerade gruppieren

Jedes Gittersegment besitzt eine Richtung, die proportional zu einem primitiven Vektor \(d=(dx,dy)\) mit \(\gcd(|dx|,|dy|)=1\) ist. Fixieren wir eine solche Richtung. Alle Gitterpunkte auf einer zu \(d\) parallelen Geraden erfüllen

$$s(x,y)=dy\,x-dx\,y=\text{konstant}.$$

Die Implementierung gruppiert Kandidaten daher zuerst nach der primitiven Richtung und dann nach dem Wert von \(s\), der genau eine Stützgerade festlegt. Auf einer festen Geraden erhält man durch wiederholte Schritte mit \(d\) die geordnete Liste

$$P_0,P_1,\dots,P_{k-1}.$$

Jede mögliche Hüllkante mit dieser Steigung ist dann ein Paar \(P_i\to P_j\) mit \(0\le i\lt j\lt k\).

Schritt 3: Wahrscheinlichkeit eines konkreten Kantenpaares

Fixiere eine Stützgerade und ein Paar \(P_i\to P_j\). Damit dieses Segment die Hüllkante mit Innerem auf der linken Seite ist, müssen vier unabhängige Bedingungen gelten:

$$\text{(a) }P_i\text{ und }P_j\text{ sind ausgewählt},$$

$$\text{(b) jeder Punkt strikt rechts der Geraden ist nicht ausgewählt},$$

$$\text{(c) jeder kollineare Punkt vor }P_i\text{ oder nach }P_j\text{ ist nicht ausgewählt},$$

$$\text{(d) mindestens ein Punkt strikt links der Geraden ist ausgewählt}.$$

Seien \(L\) die Zahl der Punkte strikt links, \(R\) die Zahl der Punkte strikt rechts und

$$O=i+(k-1-j)$$

die Zahl der kollinearen Punkte außerhalb des Segments. Dann folgt aus der Unabhängigkeit

$$\mathbb{P}(P_i\to P_j\text{ ist genau diese Hüllkante})=p^2 q^{R} q^{O}(1-q^{L}).$$

Punkte zwischen \(P_i\) und \(P_j\) auf derselben Geraden dürfen durchaus gewählt sein; sie liegen dann nur auf dem Randsegment und verändern den Flächenbeitrag nicht. Der Faktor \(1-q^L\) schließt außerdem den degenerierten Fall aus, dass alle gewählten Punkte kollinear sind und die Fläche somit \(0\) ist.

Schritt 4: Alle Paare einer Stützgeraden aufsummieren

Für eine feste Gerade hängen \(L\) und \(R\) nur von dieser Geraden ab, während \(O=i+(k-1-j)\) vom gewählten Paar abhängt. Der gesamte orientierte Erwartungsbeitrag einer Geraden ist deshalb

$$p^2 q^{R}(1-q^{L})\sum_{0\le i\lt j\lt k}(P_i\times P_j)\,q^{i}q^{k-1-j}.$$

Genau deshalb durchläuft die Implementierung alle geordneten Paare auf einer Geraden: Das Kreuzprodukt \(P_i\times P_j\) liefert den geometrischen Anteil, und die Potenzen von \(q\) codieren die Wahrscheinlichkeit, dass diese beiden Punkte die extremen ausgewählten Punkte auf dieser Geraden sind.

Schritt 5: Über alle primitiven Richtungen summieren

Nun wird der vorige Ausdruck über alle Stützgeraden aller primitiven Richtungen summiert. Jede echte Kante der konvexen Hülle wird genau einmal erfasst: Ihre Steigung bestimmt die primitive Richtung, und die Orientierung gegen den Uhrzeigersinn ist genau diejenige, bei der die Hülle links und das Äußere rechts liegt.

Daher gilt

$$\boxed{\mathbb{E}[\operatorname{Area}(H)]=\frac12\sum_{\text{primitive }d}\ \sum_{\ell\parallel d} p^2 q^{R(\ell)}(1-q^{L(\ell)})\sum_{0\le i\lt j\lt k_\ell}(P_i\times P_j)\,q^{i}q^{k_\ell-1-j}.}$$

Genau diese Formel wertet das Programm aus.

Durchgerechnetes Beispiel: \(N=1\)

Für \(N=1\) besteht das Gitter aus den vier Ecken des Einheitsquadrats und

$$p=\frac12.$$

Die Hülle hat genau dann Fläche \(1\), wenn alle vier Ecken gewählt werden; das geschieht mit Wahrscheinlichkeit \(1/16\). Fläche \(1/2\) entsteht, wenn genau drei der vier Ecken gewählt werden. Davon gibt es vier Konfigurationen, jede mit Wahrscheinlichkeit \(1/16\). Alle übrigen Fälle haben Fläche \(0\).

Also

$$\mathbb{E}[\operatorname{Area}(H)]=1\cdot\frac{1}{16}+\frac12\cdot 4\cdot\frac{1}{16}=\frac{3}{16}=0.18750,$$

genau der Kontrollwert der Implementierung.

Wie der Code Arbeitet

Die C++-, Python- und Java-Implementierungen folgen derselben Pipeline. Zuerst berechnen sie \(p\), \(q\) und eine Tabelle der Potenzen \(q^0,q^1,\dots,q^{(N+1)^2}\), denn jede Wahrscheinlichkeitskomponente wird aus diesen Werten zusammengesetzt. Danach werden alle primitiven Richtungen \((dx,dy)\) innerhalb des Quadrats aufgelistet.

Für eine feste Richtung durchläuft die Implementierung alle Gitterpunkte und berechnet \(s(x,y)=dy\,x-dx\,y\). Dadurch erhält sie zugleich die Anzahl der Punkte auf jeder parallelen Geraden und den Randpunkt, von dem aus diese Gerade einmalig traversiert werden soll. Präfixsummen über die Geradenbesetzungen liefern dann für jede Gerade die Zahl der Punkte strikt links bzw. strikt rechts.

Anschließend wird jede Gerade in Reihenfolge durchlaufen, ihre Punkte werden als \(P_0,\dots,P_{k-1}\) gespeichert, und jedes Paar \(P_i,P_j\) trägt

$$(P_i\times P_j)\,q^i q^{k-1-j}$$

zur Geradensumme bei. Multipliziert man noch den Geradenfaktor

$$p^2 q^R(1-q^L),$$

erhält man daraus einen erwarteten orientierten Flächenbeitrag. Die Summe über alle Richtungen und die abschließende Division durch \(2\) ergeben den gesuchten Erwartungswert. Laufzeitoptimierte Implementierungen teilen die Richtungsmenge auf mehrere Arbeiter auf, das mathematische Ergebnis bleibt aber identisch.

Komplexitätsanalyse

Die Zahl primitiver Richtungen mit \(|dx|,|dy|\le N\) ist \(O(N^2)\). Für eine feste Richtung kostet der Aufbau der Geradenstatistik über die \((N+1)^2\) Gitterpunkte \(O(N^2)\), während die Paarakkumulation auf einer Geraden der Länge \(k\) Aufwand \(O(k^2)\) hat. Über alle Geraden dieser Richtung gilt \(\sum k=(N+1)^2\), und im ungünstigsten Fall ist \(\sum k^2=O(N^3)\). Damit kostet eine Richtung im Worst Case \(O(N^3)\) Zeit.

Folglich beträgt der Gesamtaufwand \(O(N^5)\). Der Speicherbedarf ist \(O(N^2)\): Potenztabelle, richtungsabhängige Geradenzähler und temporäre Speicherung einer einzelnen Geraden sind sämtlich quadratisch oder kleiner. Parallelisierung verbessert die Wandzeit, ändert aber nicht die asymptotischen Schranken.

Fußnoten und Referenzen

  1. Problemseite: https://projecteuler.net/problem=514
  2. Die orientierte Flächenformel für Polygone, meist als Shoelace-Formel oder Gaußsche Trapezformel bekannt.
  3. Konvexe Hüllen, Stützgeraden und die Charakterisierung einer Hüllkante durch eine leere Halbebene.
  4. Primitive Gitterrichtungen, also teilerfremde ganzzahlige Schrittvektoren, die alle Gittersteigungen ohne Duplikate erfassen.

Problem Özeti

\(\{0,\dots,N\}^2\) kafes karesi üzerinde toplam \((N+1)^2\) nokta vardır. Her nokta bağımsız olarak

$$p=\frac{1}{N+1},\qquad q=1-p$$

olasılığıyla seçilir. \(H\), seçilen noktaların konveks kabuğu olsun. Amaç \(\mathbb{E}[\operatorname{Area}(H)]\) değerini hesaplamaktır. Tüm \(2^{(N+1)^2}\) alt kümeleri doğrudan saymak mümkün olmadığı için çözüm, beklenen alanı olası kabuk kenarlarının toplamına dönüştürür ve bunu yönlere ayırarak değerlendirir.

Matematiksel Yaklaşım

Temel fikir, alanı yönlü kenar katkılarına ayırmaktır. Bu yapıldıktan sonra beklentinin doğrusallığı sayesinde her aday kabuk kenarı ayrı ayrı incelenebilir.

Adım 1: Alanı yönlü kabuk kenarlarının toplamı olarak yaz

Kabuk köşeleri saat yönünün tersine \(v_0,\dots,v_{m-1}\) sırasıyla verilmiş olsun. Shoelace formülü şunu verir:

$$2\,\operatorname{Area}(H)=\sum_{t=0}^{m-1} v_t\times v_{t+1},\qquad v_t\times v_{t+1}=x_t y_{t+1}-y_t x_{t+1}.$$

Beklenti alınır ve doğrusallık kullanılırsa tüm sıralı nokta çiftleri \((u,v)\) üzerinden toplayabiliriz:

$$2\,\mathbb{E}[\operatorname{Area}(H)]=\sum_{u,v}(u\times v)\,\mathbb{P}(u\to v\text{ saat yönünün tersine bir kabuk kenarıdır}).$$

Böylece problem bir olasılık problemine dönüşür: Sabit bir yönlü \(u\to v\) parçası hangi koşullarda gerçekten konveks kabuğun kenarı olur?

Adım 2: Aday kenarları primitif yön ve taşıyıcı doğruya göre grupla

Her kafes parçasının yönü, \(\gcd(|dx|,|dy|)=1\) koşulunu sağlayan bir primitif vektör \(d=(dx,dy)\) ile orantılıdır. Böyle bir yönü sabitleyelim. \(d\)'ye paralel aynı doğru üzerindeki tüm kafes noktaları

$$s(x,y)=dy\,x-dx\,y=\text{sabit}$$

eşitliğini sağlar. Bu yüzden implementasyon aday kenarları önce primitif yöne, sonra da taşıyıcı doğruyu belirleyen \(s\) değerine göre gruplar. Sabit bir doğru üzerindeki noktalar, \(d\) doğrultusunda ilerlenerek

$$P_0,P_1,\dots,P_{k-1}$$

şeklinde sıralanır. Bu eğime sahip her olası kabuk kenarı da \(0\le i\lt j\lt k\) olacak biçimde bir \(P_i\to P_j\) çiftidir.

Adım 3: Bir çiftin gerçekten kabuk kenarı olma olasılığını hesapla

Bir taşıyıcı doğru ve bir \(P_i\to P_j\) çifti sabit olsun. Bu parçanın iç bölgesi solda kalan kabuk kenarı olması için dört bağımsız koşul gerekir:

$$\text{(a) }P_i\text{ ve }P_j\text{ seçilmiş olmalı},$$

$$\text{(b) doğrunun sağındaki bütün noktalar seçilmemiş olmalı},$$

$$\text{(c) aynı doğru üzerinde }P_i\text{'den önce ve }P_j\text{'den sonra kalan noktalar seçilmemiş olmalı},$$

$$\text{(d) doğrunun solunda en az bir seçilmiş nokta bulunmalı}.$$

Doğrunun solundaki nokta sayısı \(L\), sağındaki nokta sayısı \(R\) ve segment dışında aynı doğru üzerindeki kollinear nokta sayısı

$$O=i+(k-1-j)$$

olsun. Bağımsızlıktan

$$\mathbb{P}(P_i\to P_j\text{ bu kabuk kenarıdır})=p^2 q^{R} q^{O}(1-q^{L})$$

elde edilir. \(P_i\) ile \(P_j\) arasındaki aynı doğru üzerindeki noktaların seçilmesine izin verilir; bunlar sınır parçası üzerinde kalır ve alan katkısını değiştirmez. \(1-q^L\) çarpanı ayrıca tüm seçili noktaların tek bir doğru üzerinde kaldığı ve alanın sıfır olduğu dejenerik durumu dışlar.

Adım 4: Aynı taşıyıcı doğru üzerindeki tüm çiftleri topla

Sabit bir doğru için \(L\) ve \(R\) yalnızca doğruya bağlıdır; \(O=i+(k-1-j)\) ise seçilen çifte bağlıdır. Dolayısıyla tek bir doğrunun toplam yönlü beklenti katkısı

$$p^2 q^{R}(1-q^{L})\sum_{0\le i\lt j\lt k}(P_i\times P_j)\,q^{i}q^{k-1-j}$$

olur. İmplementasyondaki çift döngüsü tam olarak bunu hesaplar: \(P_i\times P_j\) geometrik terimi verir, \(q\) kuvvetleri ise bu iki noktanın kendi doğruları üzerindeki uç seçili noktalar olma olasılığını kodlar.

Adım 5: Tüm primitif yönler üzerinde topla

Şimdi bir önceki ifadeyi tüm primitif yönlerin tüm taşıyıcı doğruları üzerinde toplayalım. Konveks kabuğun her gerçek kenarı tam bir kez sayılır: eğim primitif yönü belirler, saat yönünün tersine yönlendirme ise kabuğun solda, dış bölgenin sağda kaldığı tek yönlendirmedir.

Böylece

$$\boxed{\mathbb{E}[\operatorname{Area}(H)]=\frac12\sum_{\text{primitive }d}\ \sum_{\ell\parallel d} p^2 q^{R(\ell)}(1-q^{L(\ell)})\sum_{0\le i\lt j\lt k_\ell}(P_i\times P_j)\,q^{i}q^{k_\ell-1-j}.}$$

Programın değerlendirdiği formül tam olarak budur.

Çözümlü Örnek: \(N=1\)

\(N=1\) için kafes, birim karenin dört köşesinden oluşur ve

$$p=\frac12.$$

Kabuk alanı yalnızca dört köşenin hepsi seçildiğinde \(1\) olur; bunun olasılığı \(1/16\)'dır. Tam olarak üç köşe seçildiğinde alan \(1/2\) olur. Böyle dört durum vardır ve her birinin olasılığı da \(1/16\)'dır. Geriye kalan tüm durumlarda alan \(0\)'dır.

Bu nedenle

$$\mathbb{E}[\operatorname{Area}(H)]=1\cdot\frac{1}{16}+\frac12\cdot 4\cdot\frac{1}{16}=\frac{3}{16}=0.18750,$$

ki bu değer implementasyondaki kontrol noktasıyla aynıdır.

Kod Nasıl Çalışır

C++, Python ve Java implementasyonları aynı işlem hattını izler. Önce \(p\), \(q\) ve \(q^0,q^1,\dots,q^{(N+1)^2}\) kuvvetleri hazırlanır; çünkü tüm olasılık terimleri bunlardan kurulur. Sonra kare içindeki tüm primitif yönler \((dx,dy)\) listelenir.

Tek bir yön için implementasyon bütün kafes noktalarını tarar ve \(s(x,y)=dy\,x-dx\,y\) değerini hesaplar. Bu işlem aynı anda her paralel doğrunun kaç nokta içerdiğini verir ve her doğrunun yalnızca bir kez gezilmesi için başlangıç sınır noktasını belirler. Doğru nüfuslarının prefix toplamları kullanılarak her doğru için soldaki ve sağdaki nokta sayıları bulunur.

Daha sonra her doğru sıra ile yürünür, noktalar \(P_0,\dots,P_{k-1}\) olarak tutulur ve her \(P_i,P_j\) çifti

$$(P_i\times P_j)\,q^i q^{k-1-j}$$

katkısını üretir. Bunun

$$p^2 q^R(1-q^L)$$

ile çarpılması, geometrik toplamı yönlü beklenen alan katkısına dönüştürür. Tüm yönler toplanıp sonuç \(2\)'ye bölününce son beklenen alan elde edilir. Yerel paralelleştirme destekleyen dillerde yönler işçilere bölünür, fakat hesaplanan matematiksel değer değişmez.

Karmaşıklık Analizi

\(|dx|,|dy|\le N\) koşulunu sağlayan primitif yön sayısı \(O(N^2)\)'dir. Sabit bir yön için \((N+1)^2\) kafes noktasından doğru istatistiklerini çıkarmak \(O(N^2)\) zaman alır; uzunluğu \(k\) olan tek bir doğru üzerindeki çift toplama ise \(O(k^2)\) maliyetlidir. O yönün tüm doğruları için \(\sum k=(N+1)^2\) olur ve en kötü durumda \(\sum k^2=O(N^3)\)'tür. Dolayısıyla tek bir yönün maliyeti en kötü durumda \(O(N^3)\)'tür.

Buna göre toplam iş miktarı \(O(N^5)\) olur. Bellek kullanımı \(O(N^2)\)'dir: kuvvet tablosu, yöne bağlı doğru sayımları ve tek bir doğru için geçici depolama birlikte en fazla kuadratik büyüklüktedir. Paralelleştirme çalışma süresini kısaltır ama asimptotik sınırları değiştirmez.

Dipnotlar ve Kaynaklar

  1. Problem sayfası: https://projecteuler.net/problem=514
  2. Çokgenlerde yönlü alan özdeşliği; yaygın adıyla shoelace formülü.
  3. Konveks kabuklar, taşıyıcı doğrular ve bir kabuk kenarının boş bir yarı-düzlem ile belirlenmesi.
  4. Tekrarsız bütün kafes eğimlerini veren, aralarında asal tam sayı adım vektörleri yani primitif yönler.

Resumen del Problema

Consideramos el cuadrado de la retícula \(\{0,\dots,N\}^2\), que contiene \((N+1)^2\) puntos de la cuadrícula. Cada punto se selecciona de manera independiente con probabilidad

$$p=\frac{1}{N+1},\qquad q=1-p.$$

Sea \(H\) la envolvente convexa de los puntos seleccionados. Hay que calcular \(\mathbb{E}[\operatorname{Area}(H)]\). Enumerar directamente los \(2^{(N+1)^2}\) subconjuntos es inviable, así que la solución reescribe el área esperada como una suma sobre posibles aristas del casco y luego evalúa esa suma dirección por dirección.

Enfoque Matemático

La idea central es descomponer el área en contribuciones de aristas orientadas. Hecho eso, la linealidad de la esperanza permite estudiar una arista candidata cada vez y sumar todos los aportes.

Paso 1: Reescribir el área como suma de aristas orientadas del casco

Si los vértices del casco en orden antihorario son \(v_0,\dots,v_{m-1}\), la fórmula del cordón dice que

$$2\,\operatorname{Area}(H)=\sum_{t=0}^{m-1} v_t\times v_{t+1},\qquad v_t\times v_{t+1}=x_t y_{t+1}-y_t x_{t+1}.$$

Tomando esperanza y usando linealidad, podemos sumar sobre todos los pares ordenados \((u,v)\) de puntos de la retícula:

$$2\,\mathbb{E}[\operatorname{Area}(H)]=\sum_{u,v}(u\times v)\,\mathbb{P}(u\to v\text{ aparece como arista antihoraria del casco}).$$

Así, el problema se convierte en una cuestión probabilística: para un segmento orientado fijo \(u\to v\), ¿cuándo contribuye realmente como arista de la envolvente convexa?

Paso 2: Agrupar aristas candidatas por dirección primitiva y recta soporte

Cualquier segmento de la retícula tiene una dirección proporcional a un vector primitivo \(d=(dx,dy)\) con \(\gcd(|dx|,|dy|)=1\). Fijemos una de esas direcciones. Todos los puntos de una recta paralela a \(d\) satisfacen

$$s(x,y)=dy\,x-dx\,y=\text{constante}.$$

Por eso la implementación agrupa los candidatos primero por la dirección primitiva y después por el valor de \(s\), que identifica una recta soporte concreta. En una recta fija, al ordenar los puntos avanzando por pasos de \(d\), obtenemos

$$P_0,P_1,\dots,P_{k-1}.$$

Cualquier arista posible con esa pendiente es entonces algún par \(P_i\to P_j\) con \(0\le i\lt j\lt k\).

Paso 3: Calcular la probabilidad de que un par sea la arista del casco

Fijemos una recta soporte y un par \(P_i\to P_j\). Para que ese segmento sea la arista del casco cuyo interior queda a la izquierda, deben cumplirse cuatro condiciones independientes:

$$\text{(a) }P_i\text{ y }P_j\text{ están seleccionados},$$

$$\text{(b) todos los puntos estrictamente a la derecha de la recta no están seleccionados},$$

$$\text{(c) todos los puntos colineales antes de }P_i\text{ o después de }P_j\text{ no están seleccionados},$$

$$\text{(d) al menos un punto estrictamente a la izquierda de la recta está seleccionado}.$$

Si \(L\) es el número de puntos estrictamente a la izquierda, \(R\) el número de puntos estrictamente a la derecha y

$$O=i+(k-1-j)$$

el número de puntos colineales fuera del segmento, la independencia da

$$\mathbb{P}(P_i\to P_j\text{ es esa arista del casco})=p^2 q^{R} q^{O}(1-q^{L}).$$

Los puntos situados entre \(P_i\) y \(P_j\) en la misma recta pueden estar seleccionados; siguen perteneciendo al mismo segmento de borde y no cambian el aporte de área. El factor \(1-q^L\) también elimina el caso degenerado en el que todos los puntos seleccionados son colineales y el área vale \(0\).

Paso 4: Sumar todos los pares de una recta soporte

En una recta fija, \(L\) y \(R\) dependen solo de la recta, mientras que \(O=i+(k-1-j)\) depende del par elegido. Por tanto, la contribución orientada esperada de una sola recta es

$$p^2 q^{R}(1-q^{L})\sum_{0\le i\lt j\lt k}(P_i\times P_j)\,q^{i}q^{k-1-j}.$$

Eso explica el doble bucle sobre parejas en la implementación: el producto cruzado \(P_i\times P_j\) aporta la parte geométrica, y las potencias de \(q\) codifican la probabilidad de que esos dos puntos sean los extremos seleccionados en su recta.

Paso 5: Sumar sobre todas las direcciones primitivas

Ahora se suma la expresión anterior sobre todas las rectas soporte de todas las direcciones primitivas. Cada arista real de la envolvente convexa aparece exactamente una vez: su pendiente fija la dirección primitiva, y la orientación antihoraria es la única en la que el casco queda a la izquierda y el exterior a la derecha.

En consecuencia,

$$\boxed{\mathbb{E}[\operatorname{Area}(H)]=\frac12\sum_{\text{primitive }d}\ \sum_{\ell\parallel d} p^2 q^{R(\ell)}(1-q^{L(\ell)})\sum_{0\le i\lt j\lt k_\ell}(P_i\times P_j)\,q^{i}q^{k_\ell-1-j}.}$$

Esa es exactamente la fórmula que evalúa el programa.

Ejemplo trabajado: \(N=1\)

Cuando \(N=1\), la retícula es el cuadrado unidad con sus cuatro esquinas y

$$p=\frac12.$$

La envolvente tiene área \(1\) solo cuando se seleccionan las cuatro esquinas, lo que ocurre con probabilidad \(1/16\). Tiene área \(1/2\) cuando se seleccionan exactamente tres de las cuatro esquinas; hay cuatro configuraciones así, cada una con probabilidad \(1/16\). En todos los demás casos el área es \(0\).

Por tanto,

$$\mathbb{E}[\operatorname{Area}(H)]=1\cdot\frac{1}{16}+\frac12\cdot 4\cdot\frac{1}{16}=\frac{3}{16}=0.18750,$$

que coincide con el punto de control usado por la implementación.

Cómo Funciona el Código

Las implementaciones en C++, Python y Java siguen la misma secuencia. Primero calculan \(p\), \(q\) y una tabla de potencias \(q^0,q^1,\dots,q^{(N+1)^2}\), porque todos los factores probabilísticos se construyen a partir de esos valores. Después enumeran todas las direcciones primitivas \((dx,dy)\) dentro del cuadrado.

Para una dirección dada, la implementación recorre todos los puntos de la retícula y evalúa \(s(x,y)=dy\,x-dx\,y\). Con ello obtiene al mismo tiempo cuántos puntos hay en cada recta paralela e identifica el punto de frontera desde el que debe recorrerse cada recta para visitarla exactamente una vez. Las sumas prefijo de las poblaciones de recta permiten saber, para cada recta, cuántos puntos quedan estrictamente a la izquierda y cuántos estrictamente a la derecha.

Luego se recorre cada recta en orden, se guardan sus puntos como \(P_0,\dots,P_{k-1}\), y cada pareja \(P_i,P_j\) aporta

$$(P_i\times P_j)\,q^i q^{k-1-j}$$

a la suma de esa recta. Al multiplicar por el factor de recta

$$p^2 q^R(1-q^L)$$

esa suma geométrica se convierte en una contribución esperada de área orientada. Al sumar todas las direcciones y dividir entre \(2\), se obtiene el área esperada final. Las implementaciones con soporte nativo de paralelización reparten las direcciones entre varios trabajadores, pero el resultado matemático es el mismo.

Análisis de Complejidad

El número de direcciones primitivas con \(|dx|,|dy|\le N\) es \(O(N^2)\). Para una dirección fija, construir las estadísticas de rectas sobre los \((N+1)^2\) puntos cuesta \(O(N^2)\), mientras que acumular parejas en una recta de longitud \(k\) cuesta \(O(k^2)\). Sumando sobre todas las rectas de esa dirección se tiene \(\sum k=(N+1)^2\), y en el peor caso \(\sum k^2=O(N^3)\). Por ello, una dirección cuesta \(O(N^3)\) tiempo en el peor caso.

En consecuencia, el trabajo total es \(O(N^5)\). La memoria usada es \(O(N^2)\): la tabla de potencias, los contadores de recta por dirección y el almacenamiento temporal de una recta son cuadráticos o menores. La paralelización reduce el tiempo de reloj, pero no cambia las cotas asintóticas.

Notas y Referencias

  1. Página del problema: https://projecteuler.net/problem=514
  2. La identidad de área orientada para polígonos, conocida habitualmente como fórmula del cordón.
  3. Envolventes convexas, rectas soporte y el hecho de que una arista del casco queda determinada por una semirrecta vacía.
  4. Direcciones primitivas de la retícula, es decir, vectores enteros coprimos que enumeran todas las pendientes sin duplicación.

问题概述

考虑格点正方形 \(\{0,\dots,N\}^2\),其中一共有 \((N+1)^2\) 个格点。每个点彼此独立地以概率

$$p=\frac{1}{N+1},\qquad q=1-p$$

被选中。记 \(H\) 为所有被选中点的凸包。题目要求计算 \(\mathbb{E}[\operatorname{Area}(H)]\)。如果直接枚举全部 \(2^{(N+1)^2}\) 个点集,规模完全无法接受,因此解法必须把“期望面积”改写成“所有可能凸包边的贡献之和”,再按方向逐类统计。

数学方法

核心思想是把面积拆成有向边的贡献。这样一来,利用期望的线性性质,就可以把问题化为“某一条给定线段以多大概率成为凸包边”这一类局部概率计算。

步骤 1:把面积改写成有向凸包边之和

设凸包顶点按逆时针顺序为 \(v_0,\dots,v_{m-1}\)。鞋带公式给出

$$2\,\operatorname{Area}(H)=\sum_{t=0}^{m-1} v_t\times v_{t+1},\qquad v_t\times v_{t+1}=x_t y_{t+1}-y_t x_{t+1}.$$

对上式取期望,并利用线性性质,就可以把求和改写成对所有有序格点对 \((u,v)\) 的求和:

$$2\,\mathbb{E}[\operatorname{Area}(H)]=\sum_{u,v}(u\times v)\,\mathbb{P}(u\to v\text{ 作为逆时针凸包边出现}).$$

因此,原问题转化为:对于一个固定的有向线段 \(u\to v\),它在什么条件下会真正成为凸包的一条边?

步骤 2:按原始方向和支撑直线分组

任何格点线段的方向都与某个原始向量 \(d=(dx,dy)\) 平行,其中 \(\gcd(|dx|,|dy|)=1\)。固定这样一个方向。所有与 \(d\) 平行、并且位于同一条直线上的格点都满足

$$s(x,y)=dy\,x-dx\,y=\text{常数}.$$

所以实现中先按原始方向分组,再按 \(s\) 的取值分组;这个 \(s\) 值恰好标识了一条支撑直线。在一条固定的直线上,沿着方向 \(d\) 依次前进,可以把线上格点写成

$$P_0,P_1,\dots,P_{k-1}.$$

于是,具有该斜率的任何候选凸包边都可以表示成某个 \(0\le i\lt j\lt k\) 的点对 \(P_i\to P_j\)。

步骤 3:计算一对端点成为凸包边的概率

现在固定一条支撑直线和一个候选点对 \(P_i\to P_j\)。如果这条线段要成为“内部在左侧”的凸包边,就必须同时满足四个彼此独立的条件:

$$\text{(a) }P_i\text{ 和 }P_j\text{ 都被选中},$$

$$\text{(b) 直线右侧的所有点都没有被选中},$$

$$\text{(c) 与该线段共线但位于 }P_i\text{ 之前或 }P_j\text{ 之后的点都没有被选中},$$

$$\text{(d) 直线左侧至少有一个点被选中}.$$

设直线左侧格点数为 \(L\),右侧格点数为 \(R\),再设线段之外、但仍与之共线的格点数为

$$O=i+(k-1-j).$$

由于所有格点独立选取,上述事件的概率就是

$$\mathbb{P}(P_i\to P_j\text{ 是这条凸包边})=p^2 q^{R} q^{O}(1-q^{L}).$$

需要注意的是,位于 \(P_i\) 和 \(P_j\) 之间的共线格点可以自由选中,因为它们仍然落在同一条边界线段上,不会改变面积贡献。因子 \(1-q^L\) 同时也排除了“所有被选中点都共线,凸包面积为零”的退化情形。

步骤 4:对同一条支撑直线上的所有端点对求和

对于一条固定直线,\(L\) 和 \(R\) 只取决于这条直线本身,而 \(O=i+(k-1-j)\) 取决于选了哪一对端点。因此,这条直线的总有向期望贡献为

$$p^2 q^{R}(1-q^{L})\sum_{0\le i\lt j\lt k}(P_i\times P_j)\,q^{i}q^{k-1-j}.$$

这正是实现里那层“枚举同一直线上所有点对”的原因:\(P_i\times P_j\) 给出几何上的有向面积项,而 \(q\) 的幂次则精确表达了这两个点恰好是该直线上最外侧被选中点的概率。

步骤 5:对所有原始方向求和

把上面的表达式再对所有原始方向、所有对应的支撑直线求和。每一条真正的凸包边都会被恰好计算一次:它的斜率唯一确定一个原始方向,而逆时针方向则唯一对应“凸包在左侧、外部在右侧”的取向。

于是得到

$$\boxed{\mathbb{E}[\operatorname{Area}(H)]=\frac12\sum_{\text{primitive }d}\ \sum_{\ell\parallel d} p^2 q^{R(\ell)}(1-q^{L(\ell)})\sum_{0\le i\lt j\lt k_\ell}(P_i\times P_j)\,q^{i}q^{k_\ell-1-j}.}$$

程序计算的正是这个公式。

例子:\(N=1\)

当 \(N=1\) 时,格点就是单位正方形的四个顶点,并且

$$p=\frac12.$$

只有当四个顶点全部被选中时,凸包面积才是 \(1\),其概率为 \(1/16\)。如果恰好选中四个顶点中的三个,则凸包是一个直角三角形,面积为 \(1/2\)。这样的情况共有 \(4\) 种,每种的概率都是 \(1/16\)。其余所有情形的面积都为 \(0\)。

因此

$$\mathbb{E}[\operatorname{Area}(H)]=1\cdot\frac{1}{16}+\frac12\cdot 4\cdot\frac{1}{16}=\frac{3}{16}=0.18750,$$

这与实现中的校验值完全一致。

代码如何实现

C++、Python 和 Java 实现遵循同一条计算流程。首先预处理 \(p\)、\(q\) 以及幂表 \(q^0,q^1,\dots,q^{(N+1)^2}\),因为所有概率因子都由这些数值组合而成。随后枚举正方形内部所有可能的原始方向 \((dx,dy)\)。

对某一个方向,实现会扫描全部格点并计算 \(s(x,y)=dy\,x-dx\,y\)。这一步同时完成两件事:一是统计每条平行直线上的格点数量,二是找出应该从哪个边界点开始沿该方向走,这样每条直线只会被访问一次。再利用这些直线计数的前缀和,就能得到每条直线左侧和右侧分别有多少个格点。

接下来按顺序走过每条直线,把线上点记为 \(P_0,\dots,P_{k-1}\),并让每一对 \(P_i,P_j\) 贡献

$$(P_i\times P_j)\,q^i q^{k-1-j}$$

到该直线的内部求和中。再乘上该直线对应的概率因子

$$p^2 q^R(1-q^L),$$

就得到这条直线的期望有向面积贡献。把所有方向的结果相加,最后除以 \(2\),就是所求的期望面积。带有原生并行能力的实现会把方向集合拆分给多个工作单元,但数学上计算的是同一个总和。

复杂度分析

满足 \(|dx|,|dy|\le N\) 的原始方向个数是 \(O(N^2)\)。对一个固定方向,构建 \((N+1)^2\) 个格点的直线统计需要 \(O(N^2)\) 时间;如果某条直线上有 \(k\) 个点,那么该直线上的点对累加需要 \(O(k^2)\) 时间。对该方向的所有直线求和时,有 \(\sum k=(N+1)^2\),而最坏情况下 \(\sum k^2=O(N^3)\),所以单个方向的最坏时间代价是 \(O(N^3)\)。

因此总工作量为 \(O(N^5)\)。空间复杂度是 \(O(N^2)\):幂表、按方向构造的直线计数数组以及单条直线的临时存储,规模都不超过二次。并行化只能改善实际运行时间,不会改变渐近复杂度。

脚注与参考

  1. 题目页面: https://projecteuler.net/problem=514
  2. 多边形有向面积恒等式,也就是常说的鞋带公式。
  3. 凸包、支撑直线,以及“凸包边由一侧空半平面决定”的几何事实。
  4. 原始格点方向,也就是坐标互素的整数步向量,它们可以无重复地枚举所有格点斜率。

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

Рассматривается квадратная решетка \(\{0,\dots,N\}^2\), содержащая \((N+1)^2\) узловых точек. Каждая точка выбирается независимо с вероятностью

$$p=\frac{1}{N+1},\qquad q=1-p.$$

Пусть \(H\) обозначает выпуклую оболочку выбранных точек. Требуется вычислить \(\mathbb{E}[\operatorname{Area}(H)]\). Перебор всех \(2^{(N+1)^2}\) подмножеств невозможен, поэтому решение переписывает ожидаемую площадь как сумму по возможным ребрам оболочки и затем считает эту сумму по направлениям.

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

Главная идея состоит в том, чтобы разложить площадь на вклады ориентированных ребер. После этого линейность математического ожидания позволяет отдельно рассматривать каждое потенциальное ребро оболочки и затем складывать результаты.

Шаг 1: Записать площадь как сумму по ориентированным ребрам оболочки

Пусть вершины оболочки в порядке против часовой стрелки равны \(v_0,\dots,v_{m-1}\). Тогда формула Гаусса дает

$$2\,\operatorname{Area}(H)=\sum_{t=0}^{m-1} v_t\times v_{t+1},\qquad v_t\times v_{t+1}=x_t y_{t+1}-y_t x_{t+1}.$$

Берем ожидание и используем линейность, после чего можно суммировать по всем упорядоченным парам решеточных точек \((u,v)\):

$$2\,\mathbb{E}[\operatorname{Area}(H)]=\sum_{u,v}(u\times v)\,\mathbb{P}(u\to v\text{ появляется как ребро оболочки в ориентации против часовой стрелки}).$$

Тем самым задача превращается в вероятностную: для фиксированного ориентированного отрезка \(u\to v\) надо понять, когда он действительно является ребром выпуклой оболочки.

Шаг 2: Сгруппировать кандидаты по примитивному направлению и опорной прямой

Любой решеточный отрезок имеет направление, пропорциональное примитивному вектору \(d=(dx,dy)\), где \(\gcd(|dx|,|dy|)=1\). Зафиксируем такое направление. Все точки на прямой, параллельной \(d\), удовлетворяют

$$s(x,y)=dy\,x-dx\,y=\text{constant}.$$

Поэтому реализация сначала группирует кандидаты по примитивному направлению, а затем по значению \(s\), которое и задает конкретную опорную прямую. На фиксированной прямой, если идти по шагам \(d\), точки можно упорядочить так:

$$P_0,P_1,\dots,P_{k-1}.$$

Тогда любое возможное ребро оболочки с этим наклоном является некоторой парой \(P_i\to P_j\), где \(0\le i\lt j\lt k\).

Шаг 3: Вероятность того, что выбранная пара образует ребро оболочки

Зафиксируем опорную прямую и пару \(P_i\to P_j\). Чтобы этот отрезок был ребром оболочки, у которого внутренняя часть оболочки находится слева, должны одновременно выполняться четыре независимых условия:

$$\text{(a) }P_i\text{ и }P_j\text{ выбраны},$$

$$\text{(b) все точки строго справа от прямой не выбраны},$$

$$\text{(c) все коллинеарные точки до }P_i\text{ и после }P_j\text{ не выбраны},$$

$$\text{(d) хотя бы одна точка строго слева от прямой выбрана}.$$

Пусть \(L\) — число точек строго слева, \(R\) — число точек строго справа, а

$$O=i+(k-1-j)$$

— число коллинеарных точек вне самого отрезка. Тогда из независимости получаем

$$\mathbb{P}(P_i\to P_j\text{ является этим ребром оболочки})=p^2 q^{R} q^{O}(1-q^{L}).$$

Точки между \(P_i\) и \(P_j\) на той же прямой могут быть выбраны; они остаются на той же граничной стороне и не меняют вклад в площадь. Множитель \(1-q^L\) также убирает вырожденный случай, когда все выбранные точки лежат на одной прямой и площадь равна нулю.

Шаг 4: Просуммировать все пары на одной опорной прямой

Для фиксированной прямой величины \(L\) и \(R\) зависят только от самой прямой, а \(O=i+(k-1-j)\) зависит от выбранной пары. Поэтому полный ориентированный ожидаемый вклад одной прямой равен

$$p^2 q^{R}(1-q^{L})\sum_{0\le i\lt j\lt k}(P_i\times P_j)\,q^{i}q^{k-1-j}.$$

Именно поэтому реализация перебирает все упорядоченные пары на каждой прямой: произведение \(P_i\times P_j\) дает геометрический член, а степени \(q\) кодируют вероятность того, что эти две точки являются крайними выбранными точками на данной прямой.

Шаг 5: Просуммировать по всем примитивным направлениям

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

Следовательно,

$$\boxed{\mathbb{E}[\operatorname{Area}(H)]=\frac12\sum_{\text{primitive }d}\ \sum_{\ell\parallel d} p^2 q^{R(\ell)}(1-q^{L(\ell)})\sum_{0\le i\lt j\lt k_\ell}(P_i\times P_j)\,q^{i}q^{k_\ell-1-j}.}$$

Именно эту формулу вычисляет программа.

Разобранный пример: \(N=1\)

При \(N=1\) решетка состоит из четырех вершин единичного квадрата, а

$$p=\frac12.$$

Площадь оболочки равна \(1\) только тогда, когда выбраны все четыре вершины; это происходит с вероятностью \(1/16\). Площадь равна \(1/2\), когда выбраны ровно три вершины. Таких конфигураций четыре, и каждая имеет вероятность \(1/16\). Во всех остальных случаях площадь равна \(0\).

Поэтому

$$\mathbb{E}[\operatorname{Area}(H)]=1\cdot\frac{1}{16}+\frac12\cdot 4\cdot\frac{1}{16}=\frac{3}{16}=0.18750,$$

что совпадает с контрольным значением из реализации.

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

Реализации на C++, Python и Java используют один и тот же вычислительный конвейер. Сначала они вычисляют \(p\), \(q\) и таблицу степеней \(q^0,q^1,\dots,q^{(N+1)^2}\), поскольку все вероятностные множители собираются из этих значений. Затем перечисляются все примитивные направления \((dx,dy)\) внутри квадрата.

Для одного направления реализация просматривает все решеточные точки и вычисляет \(s(x,y)=dy\,x-dx\,y\). Это одновременно дает число точек на каждой параллельной прямой и находит граничную точку, с которой следует начинать обход данной прямой, чтобы каждая прямая была посещена ровно один раз. Префиксные суммы по этим численностям позволяют затем узнать, сколько точек лежит строго слева и строго справа от каждой прямой.

После этого каждая прямая обходится по порядку, ее точки сохраняются как \(P_0,\dots,P_{k-1}\), и каждая пара \(P_i,P_j\) добавляет

$$(P_i\times P_j)\,q^i q^{k-1-j}$$

к сумме по прямой. Умножение на линейный множитель

$$p^2 q^R(1-q^L)$$

превращает этот геометрический член в ожидаемый ориентированный вклад в площадь. Затем все направления суммируются, и результат делится на \(2\). Реализации с нативной поддержкой параллелизма распределяют набор направлений между рабочими потоками, но математически вычисляется одна и та же сумма.

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

Число примитивных направлений с \(|dx|,|dy|\le N\) равно \(O(N^2)\). Для фиксированного направления построение статистики прямых по \((N+1)^2\) точкам требует \(O(N^2)\) времени, а накопление вкладов на одной прямой длины \(k\) стоит \(O(k^2)\). Если суммировать по всем прямым этого направления, то \(\sum k=(N+1)^2\), а в худшем случае \(\sum k^2=O(N^3)\), так что одно направление требует \(O(N^3)\) времени в худшем случае.

Следовательно, полный объем работы равен \(O(N^5)\). Память имеет порядок \(O(N^2)\): таблица степеней, массивы численностей прямых для текущего направления и временное хранение одной прямой не превосходят квадратичного размера. Параллельное выполнение уменьшает реальное время работы, но не меняет асимптотических оценок.

Сноски и Ссылки

  1. Страница задачи: https://projecteuler.net/problem=514
  2. Тождество ориентированной площади многоугольника, обычно называемое формулой шнуровки или формулой Гаусса.
  3. Выпуклые оболочки, опорные прямые и критерий ребра оболочки через пустую полуплоскость.
  4. Примитивные направления решетки, то есть взаимно простые целочисленные шаги, перечисляющие все решеточные наклоны без повторений.

ملخص المسألة

ننظر إلى المربع الشبكي \(\{0,\dots,N\}^2\) الذي يحتوي على \((N+1)^2\) نقطة شبكية. كل نقطة تُختار مستقلًا باحتمال

$$p=\frac{1}{N+1},\qquad q=1-p.$$

ولنرمز بـ \(H\) إلى الغلاف المحدب للنقاط المختارة. المطلوب هو حساب \(\mathbb{E}[\operatorname{Area}(H)]\). من المستحيل عمليًا تعداد جميع المجموعات الفرعية وعددها \(2^{(N+1)^2}\)، لذلك يعيد الحل كتابة المساحة المتوقعة على شكل مجموع لمساهمات الحواف الممكنة للغلاف، ثم يحسب هذا المجموع حسب الاتجاهات.

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

الفكرة الأساسية هي تفكيك المساحة إلى مساهمات حواف موجّهة. بعد ذلك تسمح خطية التوقع بتحليل كل حافة مرشحة على حدة ثم جمع النتائج.

الخطوة 1: كتابة المساحة كمجموع على الحواف الموجّهة للغلاف

إذا كانت رؤوس الغلاف المحدب مرتبة بعكس اتجاه عقارب الساعة على الصورة \(v_0,\dots,v_{m-1}\)، فإن صيغة الحذاء تعطي

$$2\,\operatorname{Area}(H)=\sum_{t=0}^{m-1} v_t\times v_{t+1},\qquad v_t\times v_{t+1}=x_t y_{t+1}-y_t x_{t+1}.$$

وبأخذ التوقع واستعمال الخطية، يمكننا الجمع على جميع الأزواج المرتبة \((u,v)\) من نقاط الشبكة. وإذا رمزنا إلى حدث ظهور \(u\to v\) كحافة من حواف الغلاف بالرمز \(E_{u,v}\)، فإن

$$2\,\mathbb{E}[\operatorname{Area}(H)]=\sum_{u,v}(u\times v)\,\mathbb{P}(E_{u,v}).$$

وهكذا تتحول المسألة إلى سؤال احتمالي: متى يكون المقطع الموجّه الثابت \(u\to v\) حافة حقيقية من حواف الغلاف المحدب؟

الخطوة 2: تجميع الحواف المرشحة حسب الاتجاه الأولي والخط الداعم

أي قطعة تصل بين نقطتين شبكيتين يكون اتجاهها متناسبًا مع متجه أولي \(d=(dx,dy)\) بحيث \(\gcd(|dx|,|dy|)=1\). ثبّت مثل هذا الاتجاه. كل النقاط الواقعة على خط موازٍ لـ \(d\) تحقق

$$s(x,y)=dy\,x-dx\,y.$$

لذلك تقوم الخوارزمية أولًا بتجميع المرشحين حسب الاتجاه الأولي، ثم حسب قيمة \(s\) التي تحدد الخط الداعم نفسه. وعلى خط ثابت، وعند السير بخطوات متتالية في اتجاه \(d\)، تصبح نقاط الخط مرتبة على الصورة

$$P_0,P_1,\dots,P_{k-1}.$$

وعليه فإن أي حافة محتملة بهذا الميل هي زوج من الشكل \(P_i\to P_j\) حيث \(0\le i\lt j\lt k\).

الخطوة 3: احتمال أن يكون زوج معين هو حافة الغلاف

لنثبّت خطًا داعمًا وزوجًا معينًا \(P_i\to P_j\). لكي يكون هذا المقطع هو حافة الغلاف التي يقع داخل الغلاف على يسارها، فلا بد من تحقق أربع شروط مستقلة:

\(P_i\) و\(P_j\) نقطتان مختارتان.

جميع النقاط الواقعة تمامًا إلى يمين الخط غير مختارة.

جميع النقاط المسايرة للخط قبل \(P_i\) أو بعد \(P_j\) غير مختارة.

توجد نقطة واحدة على الأقل مختارة تمامًا إلى يسار الخط.

إذا كان \(L\) عدد النقاط الواقعة تمامًا إلى اليسار، و\(R\) عدد النقاط الواقعة تمامًا إلى اليمين، وكان

$$O=i+(k-1-j)$$

هو عدد النقاط المسايرة للخط الواقعة خارج القطعة نفسها. وإذا رمزنا إلى حدث كون \(P_i\to P_j\) هو هذه الحافة بالرمز \(E_{i,j}\)، فإن الاستقلال يعطي

$$\mathbb{P}(E_{i,j})=p^2 q^{R} q^{O}(1-q^{L}).$$

أما النقاط الواقعة بين \(P_i\) و\(P_j\) على الخط نفسه فيمكن أن تُختار أو لا تُختار؛ فهي تبقى على نفس القطعة الحدودية ولا تغيّر مساهمة المساحة. كذلك فإن العامل \(1-q^L\) يستبعد الحالة المتدهورة التي تكون فيها جميع النقاط المختارة على استقامة واحدة، وعندها تكون المساحة صفرًا.

الخطوة 4: جمع كل الأزواج على الخط الداعم الواحد

على خط ثابت، يعتمد كل من \(L\) و\(R\) على الخط نفسه فقط، في حين يعتمد \(O=i+(k-1-j)\) على الزوج المختار. لذلك تكون المساهمة الموجّهة المتوقعة لخط واحد هي

$$p^2 q^{R}(1-q^{L})\sum_{0\le i\lt j\lt k}(P_i\times P_j)\,q^{i}q^{k-1-j}.$$

وهذا يفسر مباشرة الحلقة التي تمر على جميع الأزواج على الخط في التنفيذ: فالجداء الاتجاهي \(P_i\times P_j\) هو الحد الهندسي، في حين تمثل قوى \(q\) احتمال أن يكون هذان الطرفان هما النقطتين المختارتين الأبعد على ذلك الخط.

الخطوة 5: الجمع على جميع الاتجاهات الأولية

بعد ذلك نُجمل التعبير السابق على جميع الخطوط الداعمة لجميع الاتجاهات الأولية. كل حافة حقيقية من حواف الغلاف تُحسب مرة واحدة فقط: ميلها يحدد الاتجاه الأولي، والاتجاه العكسي لعقارب الساعة هو الاتجاه الوحيد الذي يكون فيه داخل الغلاف إلى اليسار والخارج إلى اليمين.

ومن ثم نحصل على

$$\boxed{\mathbb{E}[\operatorname{Area}(H)]=\frac12\sum_{\text{primitive }d}\ \sum_{\ell\parallel d} p^2 q^{R(\ell)}(1-q^{L(\ell)})\sum_{0\le i\lt j\lt k_\ell}(P_i\times P_j)\,q^{i}q^{k_\ell-1-j}.}$$

وهذه هي الصيغة التي يقيّمها البرنامج مباشرة.

مثال محلول: \(N=1\)

عندما \(N=1\)، تكون الشبكة هي رؤوس مربع الوحدة الأربعة، ويصبح

$$p=\frac12.$$

تساوي مساحة الغلاف \(1\) فقط عندما تُختار الزوايا الأربع كلها، واحتمال ذلك هو \(1/16\). وتساوي المساحة \(1/2\) عندما تُختار ثلاث زوايا بالضبط؛ وهناك أربع حالات من هذا النوع، واحتمال كل حالة \(1/16\). أما جميع الحالات الأخرى فمساحتها \(0\).

إذًا

$$\mathbb{E}[\operatorname{Area}(H)]=1\cdot\frac{1}{16}+\frac12\cdot 4\cdot\frac{1}{16}=\frac{3}{16}=0.18750,$$

وهو تمامًا مقدار التحقق المستخدم في التنفيذ.

كيف يعمل الكود

تتبع تطبيقات C++ وPython وJava المسار الحسابي نفسه. فهي تبدأ بحساب \(p\) و\(q\) وجدول للقوى \(q^0,q^1,\dots,q^{(N+1)^2}\)، لأن جميع العوامل الاحتمالية تُبنى من هذه القيم. بعد ذلك تُعدّ جميع الاتجاهات الأولية \((dx,dy)\) داخل المربع.

وبالنسبة إلى اتجاه واحد، يمسح التنفيذ جميع نقاط الشبكة ويحسب \(s(x,y)=dy\,x-dx\,y\). وهذه العملية تعطي في الوقت نفسه عدد النقاط على كل خط موازٍ، كما تحدد نقطة الحد التي يجب أن يبدأ منها السير على ذلك الخط حتى يُزار كل خط مرة واحدة فقط. ثم تُستخدم المجاميع السابقة لأعداد النقاط على الخطوط لمعرفة عدد النقاط الواقعة تمامًا إلى اليسار وعدد الواقعة تمامًا إلى اليمين لكل خط.

بعد ذلك يُسار على كل خط بالترتيب، وتُخزن نقاطه على هيئة \(P_0,\dots,P_{k-1}\)، ويضيف كل زوج \(P_i,P_j\) المقدار

$$(P_i\times P_j)\,q^i q^{k-1-j}$$

إلى مجموع ذلك الخط. وعند ضربه في عامل الخط

$$p^2 q^R(1-q^L)$$

يتحول المجموع الهندسي إلى مساهمة متوقعة في المساحة الموجّهة. ثم تُجمع جميع الاتجاهات ويُقسم الناتج على \(2\) للحصول على المساحة المتوقعة النهائية. وعندما تدعم اللغة التنفيذ المتوازي محليًا، تُقسّم مجموعة الاتجاهات بين عدة عمال، لكن القيمة الرياضية المحسوبة تبقى نفسها.

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

عدد الاتجاهات الأولية التي تحقق \(|dx|,|dy|\le N\) هو \(O(N^2)\). وبالنسبة إلى اتجاه ثابت، فإن بناء إحصاءات الخطوط عبر \((N+1)^2\) نقطة شبكية يحتاج إلى زمن \(O(N^2)\)، بينما يحتاج جمع الأزواج على خط طوله \(k\) إلى \(O(k^2)\). وعند الجمع على جميع خطوط ذلك الاتجاه يكون \(\sum k=(N+1)^2\)، وفي أسوأ الأحوال يكون \(\sum k^2=O(N^3)\)، ولذلك تبلغ كلفة الاتجاه الواحد \(O(N^3)\) في أسوأ حالة.

ومن ثم يكون مجموع العمل الكلي \(O(N^5)\). أما الذاكرة فهي \(O(N^2)\): جدول القوى، ومصفوفات تعداد الخطوط الخاصة بالاتجاه الجاري، والتخزين المؤقت لخط واحد، كلها تبقى ضمن رتبة تربيعية أو أقل. التنفيذ المتوازي يحسن الزمن الفعلي فقط، ولا يغير الحدود التقاربية.

الهوامش والمراجع

  1. صفحة المسألة: https://projecteuler.net/problem=514
  2. متطابقة المساحة الموجّهة للمضلعات، والمعروفة عادة بصيغة الحذاء.
  3. الأغلفة المحدبة والخطوط الداعمة وفكرة أن حافة الغلاف تتحدد بوجود نصف مستوى خالٍ على أحد الجانبين.
  4. الاتجاهات الشبكية الأولية، أي متجهات الخطوة الصحيحة المتباينة أوليًا التي تعدد جميع الميول الشبكية من دون تكرار.