Problem Summary

The ellipse in the problem is the locus of points \(Q\) such that \(QM+QG=15000\), where the foci are \(M=(-2000,1500)\) and \(G=(8000,1500)\). For an integer lattice point \(P\) outside that ellipse, two tangents can be drawn from \(P\) to the ellipse. We must count the lattice points for which the angle between those two tangent rays is greater than \(45^\circ\).

A naive geometric search would be awkward because every candidate point would seem to require solving for the tangency points. The implementations avoid that completely. They translate the ellipse to its center, derive a closed inequality for the tangent angle, and then count lattice points column by column.

Mathematical Approach

Let \((X,Y)\) be the original coordinates, and shift to centered coordinates

$$x=X-3000,\qquad y=Y-1500.$$

The midpoint of the two foci is \((3000,1500)\), so after translation the foci become \((\pm 5000,0)\). The lattice is unchanged because the shift vector is integral.

Standard Form of the Ellipse

The focal distance is \(2c=10000\), so \(c=5000\). The constant distance sum is \(2a=15000\), so \(a=7500\). Therefore

$$b^2=a^2-c^2=7500^2-5000^2=31\,250\,000.$$

The ellipse is thus

$$\frac{x^2}{a^2}+\frac{y^2}{b^2}=1,\qquad a^2=56\,250\,000,\qquad b^2=31\,250\,000.$$

All later formulas depend only on these centered quantities.

Exterior Points and the First Column Bound

Two real tangents exist only when the point lies strictly outside the ellipse. In centered coordinates this means

$$\frac{x^2}{a^2}+\frac{y^2}{b^2} \gt 1,$$

or equivalently

$$b^2x^2+a^2y^2 \gt a^2b^2.$$

For a fixed nonnegative \(x\), this gives the lower edge of the admissible vertical column. If \(x^2\gt a^2\), then \(y=0\) is already outside the ellipse. If \(x^2\le a^2\), then we need

$$y^2 \gt \frac{a^2b^2-b^2x^2}{a^2},$$

so the first valid integer \(y\) is the smallest integer strictly above that threshold.

Deriving the Tangent-Angle Formula

Take an exterior point \(P=(x_0,y_0)\). A line through \(P\) with slope \(m\) has equation \(y=mx+c\), where \(c=y_0-mx_0\). Such a line is tangent to the ellipse exactly when

$$c^2=a^2m^2+b^2.$$

Substituting \(c=y_0-mx_0\) gives a quadratic equation for the two tangent slopes:

$$\left(x_0^2-a^2\right)m^2-2x_0y_0m+\left(y_0^2-b^2\right)=0.$$

If the roots are \(m_1\) and \(m_2\), then

$$m_1+m_2=\frac{2x_0y_0}{x_0^2-a^2},\qquad m_1m_2=\frac{y_0^2-b^2}{x_0^2-a^2}.$$

The angle \(\theta\) between the two tangent rays satisfies

$$\tan^2\theta=\frac{(m_1-m_2)^2}{(1+m_1m_2)^2}=\frac{4\left(b^2x_0^2+a^2y_0^2-a^2b^2\right)}{\left(x_0^2+y_0^2-a^2-b^2\right)^2}.$$

This identity is the central geometric formula used in all three implementations.

The Director Circle and the \(45^\circ\) Test

The denominator vanishes on

$$x^2+y^2=a^2+b^2,$$

the director circle of the ellipse. On this circle the two tangents are perpendicular. That gives a clean case split.

If a point is outside the ellipse but satisfies \(x^2+y^2\le a^2+b^2\), then the tangent angle is at least \(90^\circ\), so it automatically exceeds \(45^\circ\).

If \(x^2+y^2\gt a^2+b^2\), then \(0\lt\theta\le 90^\circ\), so the problem condition is equivalent to \(\tan\theta\gt 1\). Using the previous formula yields

$$4\left(b^2x^2+a^2y^2-a^2b^2\right)-\left(a^2+b^2-x^2-y^2\right)^2 \gt 0.$$

This polynomial inequality is exactly the test performed by the implementations.

Reducing the Count to One Quadrant

The ellipse, the exterior condition, and the angle inequality all depend only on \(x^2\) and \(y^2\). Therefore the valid set is symmetric under reflections across both coordinate axes. It is enough to count the centered lattice points with \(x\ge 0\) and \(y\ge 0\), then correct for axis overcounting.

Counting a Fixed \(x\)-Column

For one nonnegative integer \(x\), the admissible \(y\)-values form a contiguous interval beginning at the first point outside the ellipse.

The first portion of that interval is automatic: every integer

$$y_{\min}(x)\le y\le \left\lfloor\sqrt{a^2+b^2-x^2}\right\rfloor$$

qualifies whenever it lies between the ellipse and the director circle.

Beyond the director circle, set \(u=y^2\) and write

$$k=a^2+b^2-x^2.$$

Then the angle condition becomes

$$-u^2+(4a^2+2k)u+4\left(b^2x^2-a^2b^2\right)-k^2 \gt 0.$$

This is a downward-opening quadratic in \(u\), so the valid values are exactly those below its larger root. Hence

$$u_{\max}(x)=\frac{4a^2+2k+\sqrt{(4a^2+2k)^2+4\left(4\left(b^2x^2-a^2b^2\right)-k^2\right)}}{2},$$

and the largest admissible integer \(y\) is \(\lfloor\sqrt{u_{\max}(x)}\rfloor\), adjusted by a tiny integer correction because the inequality is strict.

Worked Example: Why the Search Is Finite

Take the \(x\)-axis, so \(y=0\). Being outside the ellipse means \(x^2\gt a^2\). Set

$$z=x^2-a^2.$$

The \(45^\circ\) inequality becomes

$$4b^2z \gt (b^2-z)^2,$$

which simplifies to

$$z^2-6b^2z+b^4 \lt 0.$$

Therefore

$$b^2(3-2\sqrt{2}) \lt z \lt b^2(3+2\sqrt{2}).$$

In particular, every qualifying axis point satisfies \(z\lt 6b^2\), so

$$x^2 \lt a^2+6b^2.$$

The same argument on the \(y\)-axis gives \(y^2\lt b^2+6a^2\). This is why the code only scans finitely many columns and rows, using those rounded-up bounds as safe limits.

How the Code Works

Geometric Precomputation

The C++, Python, and Java implementations first reduce the problem data to \(a^2\), \(b^2\), \(a^2+b^2\), and \(a^2b^2\). Once those are known, every geometric decision becomes an integer comparison in centered coordinates.

Column-by-Column Counting

For each nonnegative integer \(x\), the implementation computes the first integer \(y\) outside the ellipse. It then counts every point up to the director circle immediately, because those points satisfy the angle condition automatically. If the column extends beyond the director circle, it solves the quadratic inequality in \(u=y^2\) to obtain the last admissible \(y\) in that column and adds the remaining points.

No tangent points are ever constructed explicitly. The algorithm works entirely with the derived inequalities, which is why it stays fast and numerically stable.

Exact Integer Arithmetic and Symmetry Repair

The implementations use integer square roots and exact integer comparisons, then make short upward or downward corrections when a strict inequality sits just past the square-root estimate. That avoids off-by-one errors near the boundary.

If \(N_{Q1}\) is the number of qualifying centered lattice points with \(x\ge 0\) and \(y\ge 0\), then the full answer is reconstructed as

$$N=4N_{Q1}-2N_x-2N_y+N_0,$$

where \(N_x\) and \(N_y\) count the qualifying points on the centered axes and \(N_0\) is the origin correction. The C++ and Java implementations parallelize the column sweep, while the Python implementation applies the same mathematics with a serial or process-based split depending on the runtime environment.

Complexity Analysis

Let \(X_{\max}\) be the safe horizontal scan bound, which is on the order of \(\sqrt{a^2+6b^2}\). The main loop processes each integer \(x\) from \(0\) to \(X_{\max}\), and each column requires only constant-time arithmetic plus a few final boundary corrections. The running time is therefore essentially linear in the scanned width.

Memory usage is \(O(1)\) beyond a handful of counters and geometric constants, or \(O(T)\) if one includes thread-local accumulators for \(T\) workers. The important point is that the algorithm never searches over tangent points on the ellipse; it counts lattice points directly from closed formulas.

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=246
  2. Ellipse: Wikipedia - Ellipse
  3. Director circle: Wikipedia - Director circle
  4. Tangent: Wikipedia - Tangent
  5. Quadratic equation: Wikipedia - Quadratic equation

Problemzusammenfassung

Die Ellipse der Aufgabe ist die Punktmenge aller \(Q\) mit \(QM+QG=15000\), wobei die Brennpunkte \(M=(-2000,1500)\) und \(G=(8000,1500)\) sind. Für einen ganzzahligen Gitterpunkt \(P\) außerhalb der Ellipse kann man zwei Tangenten an die Ellipse legen. Gesucht ist die Anzahl der Gitterpunkte, für die der von diesen beiden Tangentenstrahlen eingeschlossene Winkel größer als \(45^\circ\) ist.

Eine direkte geometrische Suche wäre unhandlich, weil man für jeden Kandidaten scheinbar die Berührpunkte bestimmen müsste. Die Implementierungen umgehen das vollständig: Sie verschieben die Ellipse in den Ursprung, leiten eine geschlossene Ungleichung für den Tangentenwinkel her und zählen danach die Gitterpunkte spaltenweise.

Mathematischer Ansatz

Bezeichne die ursprünglichen Koordinaten mit \((X,Y)\) und gehe zu zentrierten Koordinaten über:

$$x=X-3000,\qquad y=Y-1500.$$

Der Mittelpunkt der beiden Brennpunkte ist \((3000,1500)\), also werden die Brennpunkte nach der Verschiebung zu \((\pm 5000,0)\). Da der Verschiebungsvektor ganzzahlig ist, bleibt das Gitter erhalten.

Standardform der Ellipse

Die Brennpunktdistanz ist \(2c=10000\), also \(c=5000\). Die konstante Abstandssumme ist \(2a=15000\), also \(a=7500\). Damit gilt

$$b^2=a^2-c^2=7500^2-5000^2=31\,250\,000.$$

Die Ellipse hat also die Form

$$\frac{x^2}{a^2}+\frac{y^2}{b^2}=1,\qquad a^2=56\,250\,000,\qquad b^2=31\,250\,000.$$

Alle späteren Formeln hängen nur noch von diesen zentrierten Größen ab.

Äußere Punkte und die erste Spaltenschranke

Zwei reelle Tangenten existieren genau dann, wenn der Punkt strikt außerhalb der Ellipse liegt. In zentrierten Koordinaten heißt das

$$\frac{x^2}{a^2}+\frac{y^2}{b^2} \gt 1,$$

beziehungsweise

$$b^2x^2+a^2y^2 \gt a^2b^2.$$

Für festes \(x\ge 0\) liefert das sofort die Unterkante der zulässigen Spalte. Falls \(x^2\gt a^2\), ist bereits \(y=0\) außerhalb der Ellipse. Falls \(x^2\le a^2\), braucht man

$$y^2 \gt \frac{a^2b^2-b^2x^2}{a^2},$$

also ist das erste gültige ganzzahlige \(y\) die kleinste ganze Zahl strikt oberhalb dieser Schwelle.

Herleitung der Winkelgleichung für die Tangenten

Sei \(P=(x_0,y_0)\) ein äußerer Punkt. Eine Gerade durch \(P\) mit Steigung \(m\) hat die Form \(y=mx+c\), wobei \(c=y_0-mx_0\). Diese Gerade ist genau dann Tangente an die Ellipse, wenn

$$c^2=a^2m^2+b^2$$

gilt. Setzt man \(c=y_0-mx_0\) ein, erhält man für die beiden Tangentensteigungen die quadratische Gleichung

$$\left(x_0^2-a^2\right)m^2-2x_0y_0m+\left(y_0^2-b^2\right)=0.$$

Sind die beiden Wurzeln \(m_1\) und \(m_2\), dann folgt

$$m_1+m_2=\frac{2x_0y_0}{x_0^2-a^2},\qquad m_1m_2=\frac{y_0^2-b^2}{x_0^2-a^2}.$$

Für den Winkel \(\theta\) zwischen den beiden Tangentenstrahlen ergibt sich damit

$$\tan^2\theta=\frac{(m_1-m_2)^2}{(1+m_1m_2)^2}=\frac{4\left(b^2x_0^2+a^2y_0^2-a^2b^2\right)}{\left(x_0^2+y_0^2-a^2-b^2\right)^2}.$$

Diese Identität ist die zentrale geometrische Formel, auf der alle drei Implementierungen beruhen.

Direktorkreis und der \(45^\circ\)-Schwellwert

Der Nenner verschwindet auf

$$x^2+y^2=a^2+b^2,$$

dem Direktorkreis der Ellipse. Auf diesem Kreis stehen die beiden Tangenten senkrecht aufeinander. Das liefert eine saubere Fallunterscheidung.

Liegt ein Punkt außerhalb der Ellipse, aber innerhalb von \(x^2+y^2\le a^2+b^2\), dann ist der Tangentenwinkel mindestens \(90^\circ\) und damit sicher größer als \(45^\circ\).

Für Punkte mit \(x^2+y^2\gt a^2+b^2\) gilt dagegen \(0\lt\theta\le 90^\circ\). Dann ist die Aufgabenbedingung äquivalent zu \(\tan\theta\gt 1\). Mit der vorherigen Formel wird das zu

$$4\left(b^2x^2+a^2y^2-a^2b^2\right)-\left(a^2+b^2-x^2-y^2\right)^2 \gt 0.$$

Genau diese Polynom-Ungleichung wird von den Implementierungen getestet.

Die Zählung auf einen Quadranten reduzieren

Ellipse, Außenbedingung und Winkeltest hängen alle nur von \(x^2\) und \(y^2\) ab. Die gültige Punktmenge ist also spiegelsymmetrisch zu beiden Koordinatenachsen. Deshalb genügt es, die zentrierten Gitterpunkte mit \(x\ge 0\) und \(y\ge 0\) zu zählen und danach die Achsenüberzählung zu korrigieren.

Eine feste \(x\)-Spalte auszählen

Für ein festes ganzzahliges \(x\ge 0\) bilden die zulässigen \(y\)-Werte ein zusammenhängendes Intervall, das beim ersten Punkt außerhalb der Ellipse beginnt.

Der erste Teil dieses Intervalls ist automatisch gültig: Jedes ganzzahlige

$$y_{\min}(x)\le y\le \left\lfloor\sqrt{a^2+b^2-x^2}\right\rfloor$$

erfüllt die Bedingung, solange es zwischen Ellipse und Direktorkreis liegt.

Außerhalb des Direktorkreises setzt man \(u=y^2\) und

$$k=a^2+b^2-x^2.$$

Dann wird die Winkelbedingung zu

$$-u^2+(4a^2+2k)u+4\left(b^2x^2-a^2b^2\right)-k^2 \gt 0.$$

Das ist eine nach unten geöffnete quadratische Ungleichung in \(u\), also sind genau die Werte unterhalb der größeren Nullstelle zulässig. Somit gilt

$$u_{\max}(x)=\frac{4a^2+2k+\sqrt{(4a^2+2k)^2+4\left(4\left(b^2x^2-a^2b^2\right)-k^2\right)}}{2},$$

und das größte zulässige ganzzahlige \(y\) ist \(\lfloor\sqrt{u_{\max}(x)}\rfloor\), ergänzt um eine kleine ganzzahlige Korrektur wegen der strikten Ungleichung.

Durchgerechnetes Beispiel: Warum die Suche endlich ist

Betrachte die \(x\)-Achse, also \(y=0\). Außerhalb der Ellipse bedeutet das \(x^2\gt a^2\). Setze

$$z=x^2-a^2.$$

Die \(45^\circ\)-Bedingung wird dann zu

$$4b^2z \gt (b^2-z)^2,$$

also zu

$$z^2-6b^2z+b^4 \lt 0.$$

Daraus folgt

$$b^2(3-2\sqrt{2}) \lt z \lt b^2(3+2\sqrt{2}).$$

Insbesondere gilt für jeden gültigen Achsenpunkt \(z\lt 6b^2\), also

$$x^2 \lt a^2+6b^2.$$

Dasselbe Argument an der \(y\)-Achse liefert \(y^2\lt b^2+6a^2\). Genau deshalb durchsucht der Code nur endlich viele Spalten und Zeilen und verwendet diese aufgerundeten Schranken als sichere Suchgrenzen.

Wie der Code arbeitet

Geometrische Vorberechnung

Die C++-, Python- und Java-Implementierungen reduzieren die Aufgabendaten zuerst auf \(a^2\), \(b^2\), \(a^2+b^2\) und \(a^2b^2\). Sobald diese Größen feststehen, ist jede geometrische Entscheidung nur noch ein ganzzahliger Vergleich in zentrierten Koordinaten.

Spaltenweise Zählung

Für jedes nichtnegative ganzzahlige \(x\) berechnet die Implementierung zunächst das kleinste \(y\) außerhalb der Ellipse. Anschließend werden alle Punkte bis zum Direktorkreis sofort gezählt, weil sie die Winkelbedingung automatisch erfüllen. Falls die Spalte über den Direktorkreis hinausreicht, wird die quadratische Ungleichung in \(u=y^2\) gelöst, um das letzte zulässige \(y\) dieser Spalte zu bestimmen und den restlichen Beitrag zu addieren.

Berührpunkte auf der Ellipse werden niemals explizit konstruiert. Das gesamte Verfahren arbeitet direkt mit den abgeleiteten Ungleichungen, und genau deshalb bleibt es schnell und numerisch robust.

Exakte Ganzzahlarithmetik und Symmetriereparatur

Die Implementierungen verwenden ganzzahlige Quadratwurzeln und exakte Vergleiche. Danach folgen kurze Korrekturen nach oben oder unten, wenn eine strikte Ungleichung knapp jenseits der Wurzelschätzung liegt. So werden Randfehler um eins vermieden.

Ist \(N_{Q1}\) die Anzahl der gültigen zentrierten Gitterpunkte mit \(x\ge 0\) und \(y\ge 0\), dann wird die Gesamtzahl durch

$$N=4N_{Q1}-2N_x-2N_y+N_0$$

rekonstruiert, wobei \(N_x\) und \(N_y\) die gültigen Punkte auf den zentrierten Achsen zählen und \(N_0\) die Ursprungskorrektur ist. Die C++- und Java-Implementierungen parallelisieren den Spaltendurchlauf; die Python-Implementierung verwendet dieselbe Mathematik seriell oder mit prozessbasierter Aufteilung, je nach Laufzeitumgebung.

Komplexitätsanalyse

Sei \(X_{\max}\) die sichere horizontale Suchschranke, also von der Größenordnung \(\sqrt{a^2+6b^2}\). Die Hauptschleife verarbeitet jedes ganzzahlige \(x\) von \(0\) bis \(X_{\max}\), und pro Spalte fallen nur \(O(1)\) arithmetische Operationen sowie wenige Randkorrekturen an. Die Laufzeit ist daher im Wesentlichen linear in der Breite des abgesuchten Bereichs.

Der Speicherverbrauch ist \(O(1)\) abgesehen von einigen Zählern und geometrischen Konstanten, beziehungsweise \(O(T)\), wenn man thread-lokale Akkumulatoren für \(T\) Arbeiter mitzählt. Der entscheidende Effizienzgewinn ist, dass das Verfahren nicht über Berührpunkte auf der Ellipse sucht, sondern Gitterpunkte direkt über geschlossene Formeln zählt.

Fußnoten und Referenzen

  1. Problemseite: https://projecteuler.net/problem=246
  2. Ellipse: Wikipedia - Ellipse
  3. Direktorkreis: Wikipedia - Director circle
  4. Tangente: Wikipedia - Tangent
  5. Quadratische Gleichung: Wikipedia - Quadratic equation

Problem Özeti

Problemin elipsi, odakları \(M=(-2000,1500)\) ve \(G=(8000,1500)\) olan ve \(QM+QG=15000\) koşulunu sağlayan tüm \(Q\) noktalarından oluşur. Elipsin dışındaki bir tam sayı kafes noktası \(P\) için elipse iki teğet çizilebilir. İstenen şey, bu iki teğet ışını arasındaki açının \(45^\circ\)'den büyük olduğu kafes noktalarının sayısıdır.

Doğrudan geometrik arama verimsiz olurdu; çünkü her aday nokta için teğet noktalarını ayrıca çözmek gerekir gibi görünür. Uygulamalar bunu hiç yapmaz. Elipsi merkezine taşırlar, teğet açısı için kapalı bir eşitsizlik türetirler ve sonra kafes noktalarını sütun sütun sayarlar.

Matematiksel Yaklaşım

Orijinal koordinatları \((X,Y)\) ile gösterelim ve merkezlenmiş koordinatlara geçelim:

$$x=X-3000,\qquad y=Y-1500.$$

İki odağın orta noktası \((3000,1500)\) olduğundan, bu ötelemeden sonra odaklar \((\pm 5000,0)\) olur. Öteleme vektörü tam sayı olduğu için kafes yapısı değişmez.

Elipsin Standart Biçimi

Odaklar arası uzaklık \(2c=10000\) olduğundan \(c=5000\). Sabit uzaklık toplamı \(2a=15000\) olduğundan \(a=7500\). Böylece

$$b^2=a^2-c^2=7500^2-5000^2=31\,250\,000$$

elde edilir. Elips artık

$$\frac{x^2}{a^2}+\frac{y^2}{b^2}=1,\qquad a^2=56\,250\,000,\qquad b^2=31\,250\,000$$

şeklindedir. Sonraki tüm formüller yalnızca bu merkezlenmiş büyüklüklere dayanır.

Dış Noktalar ve İlk Sütun Sınırı

İki gerçek teğetin var olması için noktanın elipsin kesin olarak dışında olması gerekir. Merkezlenmiş koordinatlarda bu koşul

$$\frac{x^2}{a^2}+\frac{y^2}{b^2} \gt 1,$$

ya da eşdeğer olarak

$$b^2x^2+a^2y^2 \gt a^2b^2$$

şeklindedir. Sabit bir \(x\ge 0\) için bu, uygun sütunun alt sınırını verir. Eğer \(x^2\gt a^2\) ise \(y=0\) zaten elipsin dışındadır. Eğer \(x^2\le a^2\) ise

$$y^2 \gt \frac{a^2b^2-b^2x^2}{a^2}$$

olmalıdır; dolayısıyla ilk geçerli tam sayı \(y\), bu eşiğin kesin olarak üstündeki en küçük tam sayıdır.

Teğet Açısı Formülünü Türetmek

Dışarıdaki bir \(P=(x_0,y_0)\) noktasını alalım. \(P\) noktasından geçen ve eğimi \(m\) olan bir doğrunun denklemi \(y=mx+c\) olup \(c=y_0-mx_0\)'dır. Bu doğrunun elipse teğet olması için gerekli ve yeterli koşul

$$c^2=a^2m^2+b^2$$

olur. \(c=y_0-mx_0\) yazınca iki teğetin eğimleri için şu ikinci derece denklem çıkar:

$$\left(x_0^2-a^2\right)m^2-2x_0y_0m+\left(y_0^2-b^2\right)=0.$$

Kökler \(m_1\) ve \(m_2\) ise

$$m_1+m_2=\frac{2x_0y_0}{x_0^2-a^2},\qquad m_1m_2=\frac{y_0^2-b^2}{x_0^2-a^2}$$

elde edilir. İki teğet ışını arasındaki \(\theta\) açısı için

$$\tan^2\theta=\frac{(m_1-m_2)^2}{(1+m_1m_2)^2}=\frac{4\left(b^2x_0^2+a^2y_0^2-a^2b^2\right)}{\left(x_0^2+y_0^2-a^2-b^2\right)^2}$$

olur. Üç uygulamanın da dayandığı temel geometrik eşitlik budur.

Director Circle ve \(45^\circ\) Eşiği

Payda

$$x^2+y^2=a^2+b^2$$

üzerinde sıfır olur; bu eğri elipsin director circle'ıdır. Bu çember üzerinde iki teğet birbirine diktir. Böylece temiz bir durum ayrımı elde ederiz.

Nokta elipsin dışında ama \(x^2+y^2\le a^2+b^2\) bölgesindeyse, teğet açısı en az \(90^\circ\) olur; dolayısıyla otomatik olarak \(45^\circ\)'den büyüktür.

\(x^2+y^2\gt a^2+b^2\) olduğunda ise \(0\lt\theta\le 90^\circ\) ve problem koşulu \(\tan\theta\gt 1\) ile eşdeğerdir. Önceki formülü kullanınca

$$4\left(b^2x^2+a^2y^2-a^2b^2\right)-\left(a^2+b^2-x^2-y^2\right)^2 \gt 0$$

elde edilir. Uygulamaların test ettiği polinom eşitsizliği tam olarak budur.

Sayımı Tek Bölgeye İndirmek

Elips, dışta olma koşulu ve açı eşitsizliği yalnızca \(x^2\) ve \(y^2\)'ye bağlıdır. Bu yüzden geçerli küme her iki eksene göre simetriktir. O halde yalnızca \(x\ge 0\) ve \(y\ge 0\) olan merkezlenmiş kafes noktalarını saymak yeterlidir; sonra eksenlerdeki fazla sayımı düzeltiriz.

Sabit Bir \(x\) Sütununu Saymak

Sabit bir \(x\ge 0\) tam sayısı için geçerli \(y\) değerleri, elipsin dışındaki ilk noktadan başlayan kesintisiz bir aralık oluşturur.

Bu aralığın ilk kısmı otomatiktir: elips ile director circle arasında kalan

$$y_{\min}(x)\le y\le \left\lfloor\sqrt{a^2+b^2-x^2}\right\rfloor$$

şeklindeki tüm tam sayılar doğrudan uygundur.

Director circle'ın ötesinde \(u=y^2\) ve

$$k=a^2+b^2-x^2$$

yazalım. Açı koşulu o zaman

$$-u^2+(4a^2+2k)u+4\left(b^2x^2-a^2b^2\right)-k^2 \gt 0$$

haline gelir. Bu, \(u\) cinsinden aşağı açılan bir ikinci derece eşitsizliktir; dolayısıyla geçerli değerler büyük kökün altındadır. Yani

$$u_{\max}(x)=\frac{4a^2+2k+\sqrt{(4a^2+2k)^2+4\left(4\left(b^2x^2-a^2b^2\right)-k^2\right)}}{2}$$

olur ve en büyük geçerli tam sayı \(y\), \(\lfloor\sqrt{u_{\max}(x)}\rfloor\) değerinin sıkı eşitsizlik için yapılan küçük bir tam sayı düzeltmesiyle bulunur.

Çalışılmış Örnek: Aramanın Neden Sonlu Olduğu

\(x\)-ekseni üzerinde, yani \(y=0\) durumunu ele alalım. Elipsin dışında olmak \(x^2\gt a^2\) demektir. Şimdi

$$z=x^2-a^2$$

tanımlayalım. \(45^\circ\) eşitsizliği

$$4b^2z \gt (b^2-z)^2$$

olur; bu da

$$z^2-6b^2z+b^4 \lt 0$$

şeklinde yazılır. Buradan

$$b^2(3-2\sqrt{2}) \lt z \lt b^2(3+2\sqrt{2})$$

elde edilir. Özellikle her geçerli eksen noktası için \(z\lt 6b^2\) olduğundan

$$x^2 \lt a^2+6b^2$$

gerekir. Aynı düşünce \(y\)-ekseni için \(y^2\lt b^2+6a^2\) verir. Kodun yalnızca sonlu sayıda sütun ve satır taramasının nedeni tam olarak budur; bu yuvarlatılmış sınırlar güvenli arama üst sınırı olarak kullanılır.

Kod Nasıl Çalışır

Geometrik Ön Hesaplama

C++, Python ve Java uygulamaları önce problem verisini \(a^2\), \(b^2\), \(a^2+b^2\) ve \(a^2b^2\) sabitlerine indirger. Bu noktadan sonra her geometrik karar, merkezlenmiş koordinatlarda yapılan bir tam sayı karşılaştırmasına dönüşür.

Sütun Sütun Sayım

Her \(x\ge 0\) tam sayısı için uygulama önce elipsin dışındaki ilk \(y\) değerini bulur. Ardından director circle'a kadar olan bütün noktaları doğrudan sayar; çünkü onlar açı koşulunu otomatik olarak sağlar. Sütun director circle'ın ötesine uzanıyorsa, \(u=y^2\) cinsinden ikinci derece eşitsizlik çözülür, o sütundaki son geçerli \(y\) belirlenir ve kalan katkı eklenir.

Teğet noktaları hiçbir zaman açıkça inşa edilmez. Algoritma tamamen türetilmiş eşitsizliklerle çalışır; bu yüzden hem hızlıdır hem de sayısal olarak kararlıdır.

Kesin Tam Sayı Aritmetiği ve Simetri Düzeltmesi

Uygulamalar tam sayı karekökleri ve kesin karşılaştırmalar kullanır. Sonra sıkı eşitsizlik sınırda kaldığında, tahmini kökün biraz yukarısına ya da aşağısına giderek kısa düzeltmeler yaparlar. Böylece sınırda oluşabilecek birer birimlik hatalar engellenir.

\(x\ge 0\) ve \(y\ge 0\) için bulunan geçerli merkezlenmiş kafes noktası sayısı \(N_{Q1}\) ise, toplam sonuç

$$N=4N_{Q1}-2N_x-2N_y+N_0$$

ile elde edilir. Burada \(N_x\) ve \(N_y\), merkezlenmiş eksenlerdeki geçerli noktaları; \(N_0\) ise orijin düzeltmesini gösterir. C++ ve Java uygulamaları sütun taramasını paralelleştirir; Python uygulaması aynı matematiği çalışma ortamına göre tek süreçte ya da süreçlere bölerek uygular.

Karmaşıklık Analizi

\(X_{\max}\), güvenli yatay tarama sınırı olsun; büyüklük olarak \(\sqrt{a^2+6b^2}\) mertebesindedir. Ana döngü \(0\) ile \(X_{\max}\) arasındaki her tam sayı \(x\) için çalışır ve her sütunda sabit sayıda aritmetik işlem ile birkaç son sınır düzeltmesi yapılır. Dolayısıyla çalışma süresi, taranan genişlik bakımından özünde lineerdir.

Bellek kullanımı birkaç sayaç ve geometrik sabit dışında \(O(1)\)'dir; iş parçacığı yerel sayaçları da dahil edilirse \(T\) çalışan için \(O(T)\) denebilir. Verimliliğin esas nedeni, algoritmanın elips üzerindeki teğet noktalarını aramaması ve kafes noktalarını doğrudan kapalı form eşitsizliklerden saymasıdır.

Dipnotlar ve Referanslar

  1. Problem sayfası: https://projecteuler.net/problem=246
  2. Ellipse: Wikipedia - Ellipse
  3. Director circle: Wikipedia - Director circle
  4. Tangent: Wikipedia - Tangent
  5. Quadratic equation: Wikipedia - Quadratic equation

Resumen del Problema

La elipse del problema es el conjunto de puntos \(Q\) tales que \(QM+QG=15000\), con focos \(M=(-2000,1500)\) y \(G=(8000,1500)\). Para un punto de la retícula \(P\) situado fuera de la elipse, se pueden trazar dos tangentes desde \(P\) hasta la elipse. Debemos contar los puntos enteros para los que el ángulo entre esos dos rayos tangentes es mayor que \(45^\circ\).

Una búsqueda geométrica directa sería incómoda, porque parecería necesario calcular los puntos de tangencia para cada candidato. Las implementaciones evitan por completo ese paso. Trasladan la elipse a su centro, obtienen una desigualdad cerrada para el ángulo de las tangentes y luego cuentan puntos de la retícula columna por columna.

Enfoque Matemático

Denotemos por \((X,Y)\) las coordenadas originales y pasemos a coordenadas centradas:

$$x=X-3000,\qquad y=Y-1500.$$

El punto medio de los dos focos es \((3000,1500)\), así que después de la traslación los focos pasan a ser \((\pm 5000,0)\). Como el vector de traslación es entero, la red de puntos enteros no cambia.

Forma Estándar de la Elipse

La distancia focal es \(2c=10000\), por lo que \(c=5000\). La suma constante de distancias es \(2a=15000\), de modo que \(a=7500\). Entonces

$$b^2=a^2-c^2=7500^2-5000^2=31\,250\,000.$$

La elipse queda

$$\frac{x^2}{a^2}+\frac{y^2}{b^2}=1,\qquad a^2=56\,250\,000,\qquad b^2=31\,250\,000.$$

Todas las fórmulas posteriores dependen solo de estas cantidades centradas.

Puntos Exteriores y la Primera Cota Vertical

Existen dos tangentes reales únicamente si el punto está estrictamente fuera de la elipse. En coordenadas centradas eso equivale a

$$\frac{x^2}{a^2}+\frac{y^2}{b^2} \gt 1,$$

o, de forma equivalente,

$$b^2x^2+a^2y^2 \gt a^2b^2.$$

Para un \(x\ge 0\) fijo, esto da el borde inferior de la columna admisible. Si \(x^2\gt a^2\), entonces \(y=0\) ya está fuera de la elipse. Si \(x^2\le a^2\), necesitamos

$$y^2 \gt \frac{a^2b^2-b^2x^2}{a^2},$$

y por tanto el primer entero válido es el menor entero estrictamente por encima de ese umbral.

Deducción de la Fórmula del Ángulo de las Tangentes

Tomemos un punto exterior \(P=(x_0,y_0)\). Una recta que pasa por \(P\) con pendiente \(m\) tiene ecuación \(y=mx+c\), donde \(c=y_0-mx_0\). Esa recta es tangente a la elipse exactamente cuando

$$c^2=a^2m^2+b^2.$$

Al sustituir \(c=y_0-mx_0\), obtenemos una ecuación cuadrática para las dos pendientes tangentes:

$$\left(x_0^2-a^2\right)m^2-2x_0y_0m+\left(y_0^2-b^2\right)=0.$$

Si las raíces son \(m_1\) y \(m_2\), entonces

$$m_1+m_2=\frac{2x_0y_0}{x_0^2-a^2},\qquad m_1m_2=\frac{y_0^2-b^2}{x_0^2-a^2}.$$

El ángulo \(\theta\) entre los dos rayos tangentes satisface

$$\tan^2\theta=\frac{(m_1-m_2)^2}{(1+m_1m_2)^2}=\frac{4\left(b^2x_0^2+a^2y_0^2-a^2b^2\right)}{\left(x_0^2+y_0^2-a^2-b^2\right)^2}.$$

Esta identidad es la fórmula geométrica central que usan las tres implementaciones.

El Círculo Director y el Umbral de \(45^\circ\)

El denominador se anula sobre

$$x^2+y^2=a^2+b^2,$$

que es el círculo director de la elipse. Sobre esa circunferencia las dos tangentes son perpendiculares. Esto produce una división muy útil del problema.

Si un punto está fuera de la elipse pero satisface \(x^2+y^2\le a^2+b^2\), entonces el ángulo tangente es al menos \(90^\circ\), así que automáticamente es mayor que \(45^\circ\).

Si \(x^2+y^2\gt a^2+b^2\), entonces \(0\lt\theta\le 90^\circ\), y la condición del problema equivale a \(\tan\theta\gt 1\). Con la fórmula anterior esto se convierte en

$$4\left(b^2x^2+a^2y^2-a^2b^2\right)-\left(a^2+b^2-x^2-y^2\right)^2 \gt 0.$$

Esa es exactamente la desigualdad polinómica que evalúa la implementación.

Reducir el Conteo a un Solo Cuadrante

La elipse, la condición de estar fuera y la desigualdad angular dependen solo de \(x^2\) y \(y^2\). Por tanto, el conjunto válido es simétrico respecto de ambos ejes coordenados. Basta contar los puntos enteros centrados con \(x\ge 0\) y \(y\ge 0\), y luego corregir la sobrecuenta en los ejes.

Contar una Columna Fija

Para un entero fijo \(x\ge 0\), los valores admisibles de \(y\) forman un intervalo continuo que empieza en el primer punto situado fuera de la elipse.

La primera parte de ese intervalo es automática: todo entero

$$y_{\min}(x)\le y\le \left\lfloor\sqrt{a^2+b^2-x^2}\right\rfloor$$

sirve siempre que esté entre la elipse y el círculo director.

Más allá del círculo director, escribimos \(u=y^2\) y

$$k=a^2+b^2-x^2.$$

Entonces la condición angular se transforma en

$$-u^2+(4a^2+2k)u+4\left(b^2x^2-a^2b^2\right)-k^2 \gt 0.$$

Es una cuadrática cóncava en \(u\), así que los valores válidos son exactamente los que quedan por debajo de la raíz mayor. Por tanto,

$$u_{\max}(x)=\frac{4a^2+2k+\sqrt{(4a^2+2k)^2+4\left(4\left(b^2x^2-a^2b^2\right)-k^2\right)}}{2},$$

y el mayor entero admisible \(y\) es \(\lfloor\sqrt{u_{\max}(x)}\rfloor\), con una pequeña corrección entera porque la desigualdad es estricta.

Ejemplo Trabajado: Por Qué la Búsqueda es Finita

Tomemos el eje \(x\), es decir, \(y=0\). Estar fuera de la elipse significa \(x^2\gt a^2\). Definamos

$$z=x^2-a^2.$$

La desigualdad de \(45^\circ\) se vuelve

$$4b^2z \gt (b^2-z)^2,$$

o de forma equivalente

$$z^2-6b^2z+b^4 \lt 0.$$

De ahí se deduce

$$b^2(3-2\sqrt{2}) \lt z \lt b^2(3+2\sqrt{2}).$$

En particular, todo punto válido del eje satisface \(z\lt 6b^2\), así que

$$x^2 \lt a^2+6b^2.$$

El mismo razonamiento en el eje \(y\) produce \(y^2\lt b^2+6a^2\). Por eso el código solo necesita recorrer un número finito de columnas y filas y usa esas cotas redondeadas hacia arriba como límites seguros.

Cómo Funciona el Código

Precomputación Geométrica

Las implementaciones en C++, Python y Java reducen primero los datos del problema a \(a^2\), \(b^2\), \(a^2+b^2\) y \(a^2b^2\). A partir de ahí, toda decisión geométrica se convierte en una comparación entera en coordenadas centradas.

Barrido Columna por Columna

Para cada entero \(x\ge 0\), la implementación calcula primero el menor \(y\) que ya queda fuera de la elipse. Después cuenta de inmediato todos los puntos hasta el círculo director, porque cumplen automáticamente la condición angular. Si la columna continúa más allá del círculo director, resuelve la desigualdad cuadrática en \(u=y^2\), obtiene el último \(y\) admisible de esa columna y suma el resto.

Nunca se construyen explícitamente los puntos de tangencia. Todo el algoritmo trabaja con las desigualdades derivadas, y eso explica tanto su rapidez como su estabilidad numérica.

Aritmética Entera Exacta y Corrección por Simetría

Las implementaciones usan raíces cuadradas enteras y comparaciones exactas. Después aplican pequeñas correcciones hacia arriba o hacia abajo cuando una desigualdad estricta queda justo más allá de la estimación inicial. Así se evitan errores de una unidad en la frontera.

Si \(N_{Q1}\) es el número de puntos enteros centrados válidos con \(x\ge 0\) y \(y\ge 0\), el total se reconstruye mediante

$$N=4N_{Q1}-2N_x-2N_y+N_0,$$

donde \(N_x\) y \(N_y\) cuentan los puntos válidos sobre los ejes centrados y \(N_0\) es la corrección del origen. Las versiones en C++ y Java paralelizan el barrido de columnas; la versión en Python aplica la misma matemática de forma secuencial o con división por procesos según el entorno de ejecución.

Análisis de Complejidad

Sea \(X_{\max}\) la cota horizontal segura de exploración, del orden de \(\sqrt{a^2+6b^2}\). El bucle principal recorre todos los enteros \(x\) entre \(0\) y \(X_{\max}\), y cada columna requiere solo una cantidad constante de aritmética, más unas pocas correcciones finales de frontera. Por tanto, el tiempo de ejecución es esencialmente lineal en el ancho escaneado.

El uso de memoria es \(O(1)\) aparte de unos pocos contadores y constantes geométricas, u \(O(T)\) si se incluyen acumuladores locales para \(T\) trabajadores. La clave de la eficiencia es que el algoritmo no busca puntos de tangencia sobre la elipse; cuenta directamente puntos de la retícula a partir de fórmulas cerradas.

Notas y Referencias

  1. Página del problema: https://projecteuler.net/problem=246
  2. Ellipse: Wikipedia - Ellipse
  3. Director circle: Wikipedia - Director circle
  4. Tangent: Wikipedia - Tangent
  5. Quadratic equation: Wikipedia - Quadratic equation

问题概述

题目中的椭圆由条件 \(QM+QG=15000\) 定义,其中两个焦点分别是 \(M=(-2000,1500)\) 与 \(G=(8000,1500)\)。对每一个位于椭圆外部的整点 \(P\),都可以向椭圆作两条切线。我们要统计的是:有多少个整点 \(P\),使得这两条切线射线之间的夹角大于 \(45^\circ\)。

如果直接从几何图形出发逐点尝试,看起来就必须对每个候选点显式求出切点,这样既麻烦又低效。三份实现都没有这样做。它们先把椭圆平移到中心坐标系,推导出一个只含 \(x\)、\(y\) 的封闭不等式,然后按纵列逐列统计整点。

数学方法

设原始坐标为 \((X,Y)\),改用以椭圆中心为原点的坐标:

$$x=X-3000,\qquad y=Y-1500.$$

两个焦点的中点正好是 \((3000,1500)\),因此平移后焦点变成 \((\pm 5000,0)\)。由于平移向量本身是整向量,整点格不会发生变化,所以在新坐标系中统计整点与原题完全等价。

把椭圆写成标准方程

焦距的一半为 \(c=5000\),因为两焦点相距 \(10000\)。题目给出的常和是 \(15000\),所以长半轴 \(a=7500\)。于是

$$b^2=a^2-c^2=7500^2-5000^2=31\,250\,000.$$

椭圆的标准方程就是

$$\frac{x^2}{a^2}+\frac{y^2}{b^2}=1,\qquad a^2=56\,250\,000,\qquad b^2=31\,250\,000.$$

后面的所有判断都只依赖这些中心化之后的量。

外点条件与每一列的起点

只有当点严格位于椭圆外部时,才能向椭圆作出两条实切线。用中心坐标表示,这个条件是

$$\frac{x^2}{a^2}+\frac{y^2}{b^2} \gt 1,$$

也就是

$$b^2x^2+a^2y^2 \gt a^2b^2.$$

对固定的非负整数 \(x\) 来说,这个条件立刻给出了该列允许的最低 \(y\)。如果 \(x^2\gt a^2\),那么 \(y=0\) 已经在椭圆外部;如果 \(x^2\le a^2\),就必须满足

$$y^2 \gt \frac{a^2b^2-b^2x^2}{a^2}.$$

因此该列的第一个合法整点,就是严格超过这个阈值的最小整数 \(y\)。

推导切线夹角公式

取一个外点 \(P=(x_0,y_0)\)。过 \(P\) 且斜率为 \(m\) 的直线可写成 \(y=mx+c\),其中 \(c=y_0-mx_0\)。这条直线与椭圆相切,当且仅当

$$c^2=a^2m^2+b^2.$$

把 \(c=y_0-mx_0\) 代入,就得到两条切线斜率满足的二次方程:

$$\left(x_0^2-a^2\right)m^2-2x_0y_0m+\left(y_0^2-b^2\right)=0.$$

若两个根是 \(m_1,m_2\),则

$$m_1+m_2=\frac{2x_0y_0}{x_0^2-a^2},\qquad m_1m_2=\frac{y_0^2-b^2}{x_0^2-a^2}.$$

两条切线射线之间的夹角记为 \(\theta\),则有

$$\tan^2\theta=\frac{(m_1-m_2)^2}{(1+m_1m_2)^2}=\frac{4\left(b^2x_0^2+a^2y_0^2-a^2b^2\right)}{\left(x_0^2+y_0^2-a^2-b^2\right)^2}.$$

这就是三份实现共同使用的核心几何恒等式。

导圆与 \(45^\circ\) 判定

上式分母在

$$x^2+y^2=a^2+b^2$$

上为零,这条曲线正是椭圆的导圆。位于导圆上的点看到的两条切线互相垂直,因此这里自然形成了一个非常有用的分区。

如果某点在椭圆外部,同时满足 \(x^2+y^2\le a^2+b^2\),那么它对应的切线夹角至少是 \(90^\circ\),当然也就一定大于 \(45^\circ\)。

如果 \(x^2+y^2\gt a^2+b^2\),那么 \(0\lt\theta\le 90^\circ\),题目的要求就等价于 \(\tan\theta\gt 1\)。把前面的公式代入,得到

$$4\left(b^2x^2+a^2y^2-a^2b^2\right)-\left(a^2+b^2-x^2-y^2\right)^2 \gt 0.$$

这正是实现中真正用于筛选点的多项式不等式。

把平面计数压缩到第一象限

椭圆方程、外点条件以及角度不等式都只含 \(x^2\) 与 \(y^2\),所以合法点集关于两条坐标轴完全对称。因此只统计中心坐标系中满足 \(x\ge 0\)、\(y\ge 0\) 的整点就够了,最后再对坐标轴上的重复计数做修正。

固定一列 \(x\) 时如何计数

对某个固定的非负整数 \(x\),合法的 \(y\) 值会形成一个从“刚离开椭圆”开始的连续区间。

这个区间的第一部分不需要额外判断:凡是位于椭圆与导圆之间的整数

$$y_{\min}(x)\le y\le \left\lfloor\sqrt{a^2+b^2-x^2}\right\rfloor$$

都自动满足角度条件。

超过导圆以后,记 \(u=y^2\),并设

$$k=a^2+b^2-x^2.$$

那么角度不等式会变成

$$-u^2+(4a^2+2k)u+4\left(b^2x^2-a^2b^2\right)-k^2 \gt 0.$$

这是一个开口向下的关于 \(u\) 的二次不等式,因此合法范围恰好落在较大实根之下。于是

$$u_{\max}(x)=\frac{4a^2+2k+\sqrt{(4a^2+2k)^2+4\left(4\left(b^2x^2-a^2b^2\right)-k^2\right)}}{2},$$

而该列最大的合法整数 \(y\) 就是 \(\lfloor\sqrt{u_{\max}(x)}\rfloor\),再加上一点点整数修正,以确保严格不等式处理正确。

具体例子:为什么搜索范围一定有限

先看 \(x\) 轴,即 \(y=0\) 的情况。位于椭圆外意味着 \(x^2\gt a^2\)。令

$$z=x^2-a^2.$$

则 \(45^\circ\) 条件化为

$$4b^2z \gt (b^2-z)^2,$$

也就是

$$z^2-6b^2z+b^4 \lt 0.$$

因此有

$$b^2(3-2\sqrt{2}) \lt z \lt b^2(3+2\sqrt{2}).$$

特别地,必有 \(z\lt 6b^2\),所以任何合法的轴上点都满足

$$x^2 \lt a^2+6b^2.$$

对 \(y\) 轴做同样分析,可得 \(y^2\lt b^2+6a^2\)。这正是代码只需要扫描有限多列和有限多行的原因;实现把这些略向上取整后的界当作安全搜索范围。

代码如何工作

几何常量预处理

C++、Python 与 Java 实现都先把题目数据压缩成 \(a^2\)、\(b^2\)、\(a^2+b^2\) 和 \(a^2b^2\) 这几个不变量。做到这一步之后,所有几何判断都能写成中心坐标系下的整数比较。

逐列扫描第一象限

对每个非负整数 \(x\),实现先求出该列中第一个落在椭圆外部的整数 \(y\)。接着,直到导圆为止的点可以直接计入,因为那一整段都会自动满足夹角条件。如果这一列在导圆之外还有合法部分,就求解关于 \(u=y^2\) 的二次不等式,得到该列最后一个可取的 \(y\),并把剩余整点数量加入总和。

整个过程中从不显式构造椭圆上的切点。算法完全依赖前面推导出的不等式,因此既快又稳定。

精确整数运算与对称修正

实现使用整数平方根和精确整数比较,然后在边界附近做很短的向上或向下修正,以处理严格不等式带来的边界偏差,从而避免一格之差的错误。

如果 \(N_{Q1}\) 表示中心坐标系中 \(x\ge 0\)、\(y\ge 0\) 的合法整点个数,那么总答案通过

$$N=4N_{Q1}-2N_x-2N_y+N_0$$

恢复,其中 \(N_x\) 和 \(N_y\) 分别是中心坐标轴上的合法点数,\(N_0\) 是原点修正项。C++ 与 Java 版本把列扫描并行化;Python 版本则在运行环境允许时使用顺序或多进程分块,但数学逻辑完全相同。

复杂度分析

设 \(X_{\max}\) 为安全的水平扫描上界,它的量级是 \(\sqrt{a^2+6b^2}\)。主循环依次处理 \(0\) 到 \(X_{\max}\) 的每个整数 \(x\),每一列只需要常数次算术操作,再加上极少量边界修正,所以总运行时间本质上与扫描宽度成线性关系。

除少数计数器和几何常量之外,空间使用量是 \(O(1)\);如果把 \(T\) 个工作单元各自的局部累加器也算进去,可以写成 \(O(T)\)。算法高效的关键在于:它从不搜索切点,而是直接用闭式不等式统计整点。

注释与参考资料

  1. 题目页面:https://projecteuler.net/problem=246
  2. Ellipse:Wikipedia - Ellipse
  3. Director circle:Wikipedia - Director circle
  4. Tangent:Wikipedia - Tangent
  5. Quadratic equation:Wikipedia - Quadratic equation

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

Эллипс в задаче задан условием \(QM+QG=15000\), где фокусы равны \(M=(-2000,1500)\) и \(G=(8000,1500)\). Для любой целочисленной точки решетки \(P\), лежащей вне эллипса, к эллипсу можно провести две касательные. Нужно посчитать, для скольких целых точек угол между этими двумя касательными лучами больше \(45^\circ\).

Если атаковать задачу напрямую, кажется, что для каждой точки пришлось бы отдельно находить точки касания. Все три реализации обходят это полностью. Они переносят эллипс в систему координат с центром в начале, выводят замкнутое неравенство для угла между касательными, а затем считают решетчатые точки по вертикальным столбцам.

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

Обозначим исходные координаты через \((X,Y)\) и перейдем к центрированным координатам:

$$x=X-3000,\qquad y=Y-1500.$$

Середина отрезка между фокусами равна \((3000,1500)\), поэтому после сдвига фокусы становятся \((\pm 5000,0)\). Так как вектор сдвига целочисленный, решетка целых точек остается той же самой.

Стандартное уравнение эллипса

Половина расстояния между фокусами равна \(c=5000\), потому что сами фокусы удалены на \(10000\). Постоянная сумма расстояний равна \(15000\), значит большая полуось \(a=7500\). Тогда

$$b^2=a^2-c^2=7500^2-5000^2=31\,250\,000.$$

Стандартное уравнение эллипса принимает вид

$$\frac{x^2}{a^2}+\frac{y^2}{b^2}=1,\qquad a^2=56\,250\,000,\qquad b^2=31\,250\,000.$$

Дальнейшая математика зависит только от этих центрированных параметров.

Внешние точки и нижняя граница столбца

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

$$\frac{x^2}{a^2}+\frac{y^2}{b^2} \gt 1,$$

или, что то же самое,

$$b^2x^2+a^2y^2 \gt a^2b^2.$$

Для фиксированного \(x\ge 0\) это сразу дает нижнюю границу допустимого столбца. Если \(x^2\gt a^2\), то уже \(y=0\) находится вне эллипса. Если \(x^2\le a^2\), требуется

$$y^2 \gt \frac{a^2b^2-b^2x^2}{a^2},$$

поэтому первым допустимым целым \(y\) будет наименьшее целое число, строго превышающее этот порог.

Вывод формулы для угла между касательными

Возьмем внешнюю точку \(P=(x_0,y_0)\). Прямая через \(P\) с угловым коэффициентом \(m\) имеет вид \(y=mx+c\), где \(c=y_0-mx_0\). Она касается эллипса тогда и только тогда, когда

$$c^2=a^2m^2+b^2.$$

Подставляя \(c=y_0-mx_0\), получаем квадратное уравнение для двух наклонов касательных:

$$\left(x_0^2-a^2\right)m^2-2x_0y_0m+\left(y_0^2-b^2\right)=0.$$

Если его корни равны \(m_1\) и \(m_2\), то

$$m_1+m_2=\frac{2x_0y_0}{x_0^2-a^2},\qquad m_1m_2=\frac{y_0^2-b^2}{x_0^2-a^2}.$$

Для угла \(\theta\) между двумя касательными лучами отсюда следует

$$\tan^2\theta=\frac{(m_1-m_2)^2}{(1+m_1m_2)^2}=\frac{4\left(b^2x_0^2+a^2y_0^2-a^2b^2\right)}{\left(x_0^2+y_0^2-a^2-b^2\right)^2}.$$

Именно эта формула является ключевой геометрической инвариантой во всех трех реализациях.

Директорная окружность и порог \(45^\circ\)

Знаменатель обращается в нуль на кривой

$$x^2+y^2=a^2+b^2,$$

которая называется директорной окружностью эллипса. На ней касательные взаимно перпендикулярны. Это дает удобное разбиение на случаи.

Если точка лежит вне эллипса, но при этом удовлетворяет \(x^2+y^2\le a^2+b^2\), то угол между касательными не меньше \(90^\circ\), а значит автоматически больше \(45^\circ\).

Если же \(x^2+y^2\gt a^2+b^2\), то \(0\lt\theta\le 90^\circ\), и условие задачи эквивалентно \(\tan\theta\gt 1\). Подстановка предыдущей формулы дает

$$4\left(b^2x^2+a^2y^2-a^2b^2\right)-\left(a^2+b^2-x^2-y^2\right)^2 \gt 0.$$

Это в точности то полиномиальное неравенство, которое проверяют реализации.

Сведение подсчета к одному квадранту

Уравнение эллипса, условие внешности и угловая проверка зависят только от \(x^2\) и \(y^2\). Следовательно, множество подходящих точек симметрично относительно обеих координатных осей. Достаточно посчитать центрированные целые точки с \(x\ge 0\) и \(y\ge 0\), а затем исправить переучет на осях.

Подсчет в фиксированном столбце \(x\)

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

Первая часть этого интервала определяется автоматически: любой целый

$$y_{\min}(x)\le y\le \left\lfloor\sqrt{a^2+b^2-x^2}\right\rfloor$$

подходит, если находится между эллипсом и директорной окружностью.

Дальше, за пределами директорной окружности, вводим \(u=y^2\) и

$$k=a^2+b^2-x^2.$$

Тогда угловое условие принимает вид

$$-u^2+(4a^2+2k)u+4\left(b^2x^2-a^2b^2\right)-k^2 \gt 0.$$

Это квадратное неравенство по \(u\) с ветвями вниз, поэтому допустимы ровно значения ниже большего корня. Следовательно,

$$u_{\max}(x)=\frac{4a^2+2k+\sqrt{(4a^2+2k)^2+4\left(4\left(b^2x^2-a^2b^2\right)-k^2\right)}}{2},$$

а максимальное допустимое целое \(y\) равно \(\lfloor\sqrt{u_{\max}(x)}\rfloor\) с небольшой целочисленной поправкой из-за строгого неравенства.

Разобранный пример: почему перебор конечен

Рассмотрим ось \(x\), то есть случай \(y=0\). Быть вне эллипса означает \(x^2\gt a^2\). Положим

$$z=x^2-a^2.$$

Тогда условие \(45^\circ\) превращается в

$$4b^2z \gt (b^2-z)^2,$$

то есть в

$$z^2-6b^2z+b^4 \lt 0.$$

Отсюда получаем

$$b^2(3-2\sqrt{2}) \lt z \lt b^2(3+2\sqrt{2}).$$

В частности, для любой подходящей осевой точки выполняется \(z\lt 6b^2\), а значит

$$x^2 \lt a^2+6b^2.$$

Тот же аргумент на оси \(y\) дает \(y^2\lt b^2+6a^2\). Именно поэтому код просматривает только конечное число столбцов и строк и берет эти округленные вверх оценки как безопасные пределы перебора.

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

Предварительное вычисление геометрических констант

Реализации на C++, Python и Java сначала сводят данные задачи к четырем инвариантам: \(a^2\), \(b^2\), \(a^2+b^2\) и \(a^2b^2\). После этого любое геометрическое решение записывается как целочисленное сравнение в центрированной системе координат.

Поштучный просмотр столбцов

Для каждого неотрицательного целого \(x\) реализация находит наименьшее \(y\), уже лежащее вне эллипса. Затем она сразу засчитывает все точки вплоть до директорной окружности, потому что там угловое условие выполняется автоматически. Если столбец продолжается за директорной окружностью, решается квадратное неравенство по \(u=y^2\), находится последнее допустимое \(y\) в этом столбце и добавляется оставшийся вклад.

Точки касания на эллипсе никогда не строятся явно. Алгоритм целиком работает с выведенными неравенствами, что и обеспечивает скорость и устойчивость.

Точная целочисленная арифметика и восстановление симметрией

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

Если \(N_{Q1}\) обозначает число подходящих центрированных точек с \(x\ge 0\) и \(y\ge 0\), то полный ответ восстанавливается по формуле

$$N=4N_{Q1}-2N_x-2N_y+N_0,$$

где \(N_x\) и \(N_y\) подсчитывают подходящие точки на центрированных осях, а \(N_0\) отвечает за поправку в начале координат. Версии на C++ и Java распараллеливают просмотр столбцов, а версия на Python применяет ту же самую математику последовательно или с разбиением по процессам, в зависимости от среды выполнения.

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

Пусть \(X_{\max}\) обозначает безопасную горизонтальную границу просмотра; она имеет порядок \(\sqrt{a^2+6b^2}\). Основной цикл обрабатывает все целые \(x\) от \(0\) до \(X_{\max}\), причем на каждый столбец приходится лишь константное число арифметических операций и несколько коротких поправок на границе. Поэтому время работы по существу линейно по ширине просматриваемого диапазона.

Память равна \(O(1)\), если не считать нескольких счетчиков и геометрических констант, либо \(O(T)\), если включать локальные накопители для \(T\) потоков или процессов. Главный выигрыш состоит в том, что алгоритм не перебирает точки касания на эллипсе, а считает решетчатые точки напрямую по закрытым формулам.

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

  1. Страница задачи: https://projecteuler.net/problem=246
  2. Ellipse: Wikipedia - Ellipse
  3. Director circle: Wikipedia - Director circle
  4. Tangent: Wikipedia - Tangent
  5. Quadratic equation: Wikipedia - Quadratic equation

ملخص المسألة

القطع الناقص في هذه المسألة هو مجموعة النقاط \(Q\) التي تحقق \(QM+QG=15000\)، حيث البؤرتان هما \(M=(-2000,1500)\) و\(G=(8000,1500)\). إذا كانت \(P\) نقطة شبكية صحيحة خارج القطع الناقص، فيمكن رسم مماسين من \(P\) إلى القطع الناقص. المطلوب هو عدّ نقاط الشبكة الصحيحة التي تكون الزاوية بين شعاعي هذين المماسين فيها أكبر من \(45^\circ\).

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

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

لنرمز إلى الإحداثيات الأصلية بـ \((X,Y)\)، ثم ننتقل إلى إحداثيات متمركزة:

$$x=X-3000,\qquad y=Y-1500.$$

منتصف القطعة الواصلة بين البؤرتين هو \((3000,1500)\)، ولذلك تصبح البؤرتان بعد الإزاحة \((\pm 5000,0)\). وبما أن متجه الإزاحة ذو إحداثيات صحيحة، فإن شبكة النقاط الصحيحة تبقى كما هي.

الصيغة القياسية للقطع الناقص

نصف البعد البؤري هو \(c=5000\) لأن المسافة بين البؤرتين تساوي \(10000\). ومجموع المسافتين الثابت هو \(15000\)، لذا فإن نصف المحور الأكبر \(a=7500\). ومن ثم

$$b^2=a^2-c^2=7500^2-5000^2=31\,250\,000.$$

فتصبح معادلة القطع الناقص

$$\frac{x^2}{a^2}+\frac{y^2}{b^2}=1,\qquad a^2=56\,250\,000,\qquad b^2=31\,250\,000.$$

وجميع الصيغ اللاحقة تعتمد فقط على هذه القيم بعد التمركز.

شرط النقطة الخارجية وبداية كل عمود

وجود مماسين حقيقيين يتطلب أن تكون النقطة خارج القطع الناقص بصورة صارمة. في الإحداثيات المتمركزة يصبح الشرط

$$\frac{x^2}{a^2}+\frac{y^2}{b^2} \gt 1,$$

أو بصورة مكافئة

$$b^2x^2+a^2y^2 \gt a^2b^2.$$

ولكل \(x\ge 0\) ثابت، يعطينا هذا الحد الأدنى للعمود المسموح. إذا كان \(x^2\gt a^2\)، فإن \(y=0\) يقع أصلًا خارج القطع الناقص. وإذا كان \(x^2\le a^2\)، فلا بد من تحقق

$$y^2 \gt \frac{a^2b^2-b^2x^2}{a^2},$$

وبالتالي فإن أول قيمة صحيحة صالحة لـ \(y\) هي أصغر عدد صحيح يقع فوق هذا الحد مباشرة.

اشتقاق صيغة زاوية المماسين

لنأخذ نقطة خارجية \(P=(x_0,y_0)\). المستقيم المارّ بـ \(P\) والميل فيه \(m\) يكتب على الصورة \(y=mx+c\)، حيث \(c=y_0-mx_0\). ويكون هذا المستقيم مماسًا للقطع الناقص إذا وفقط إذا تحقق

$$c^2=a^2m^2+b^2.$$

وبالتعويض عن \(c\) نحصل على معادلة تربيعية لميلَي المماسين:

$$\left(x_0^2-a^2\right)m^2-2x_0y_0m+\left(y_0^2-b^2\right)=0.$$

إذا كان الجذران هما \(m_1\) و\(m_2\)، فإن

$$m_1+m_2=\frac{2x_0y_0}{x_0^2-a^2},\qquad m_1m_2=\frac{y_0^2-b^2}{x_0^2-a^2}.$$

ومن ثم فإن الزاوية \(\theta\) بين شعاعي المماسين تحقق

$$\tan^2\theta=\frac{(m_1-m_2)^2}{(1+m_1m_2)^2}=\frac{4\left(b^2x_0^2+a^2y_0^2-a^2b^2\right)}{\left(x_0^2+y_0^2-a^2-b^2\right)^2}.$$

هذه هي المتطابقة الهندسية الأساسية التي تعتمد عليها جميع التطبيقات.

دائرة الدليل وحدّ \(45^\circ\)

يصبح المقام صفرًا على المنحنى

$$x^2+y^2=a^2+b^2,$$

وهو دائرة الدليل الخاصة بالقطع الناقص. على هذه الدائرة يكون المماسان متعامدين، وهذا يمنحنا تقسيمًا واضحًا للحالات.

إذا كانت النقطة خارج القطع الناقص ولكنها تحقق \(x^2+y^2\le a^2+b^2\)، فإن زاوية المماسين لا تقل عن \(90^\circ\)، وبالتالي فهي بالتأكيد أكبر من \(45^\circ\).

أما إذا كان \(x^2+y^2\gt a^2+b^2\)، فحينئذٍ \(0\lt\theta\le 90^\circ\)، ويصبح شرط المسألة مكافئًا لـ \(\tan\theta\gt 1\). وباستخدام الصيغة السابقة نحصل على

$$4\left(b^2x^2+a^2y^2-a^2b^2\right)-\left(a^2+b^2-x^2-y^2\right)^2 \gt 0.$$

وهذه هي المتراجحة متعددة الحدود التي تختبرها التطبيقات فعليًا.

اختزال العد إلى ربع واحد

معادلة القطع الناقص، وشرط الخروج عنه، واختبار الزاوية كلها تعتمد فقط على \(x^2\) و\(y^2\). لذا فإن مجموعة النقاط المقبولة متناظرة حول المحورين. ولذلك يكفي أن نعدّ النقاط الصحيحة المتمركزة التي تحقق \(x\ge 0\) و\(y\ge 0\)، ثم نصحح فرط العد على المحاور.

كيف نعدّ عمودًا ثابتًا

عند تثبيت عدد صحيح غير سالب \(x\)، فإن قيم \(y\) المقبولة تشكل فترة متصلة تبدأ من أول نقطة تقع خارج القطع الناقص.

الجزء الأول من هذه الفترة سهل مباشرة: كل عدد صحيح يحقق

$$y_{\min}(x)\le y\le \left\lfloor\sqrt{a^2+b^2-x^2}\right\rfloor$$

يكون صالحًا ما دام يقع بين القطع الناقص ودائرة الدليل.

بعد دائرة الدليل نكتب \(u=y^2\)، ونعرّف

$$k=a^2+b^2-x^2.$$

فتتحول متراجحة الزاوية إلى

$$-u^2+(4a^2+2k)u+4\left(b^2x^2-a^2b^2\right)-k^2 \gt 0.$$

وهذه متراجحة تربيعية مقعّرة في \(u\)، لذلك تكون القيم المقبولة هي تمامًا القيم الواقعة دون الجذر الأكبر. وعليه فإن

$$u_{\max}(x)=\frac{4a^2+2k+\sqrt{(4a^2+2k)^2+4\left(4\left(b^2x^2-a^2b^2\right)-k^2\right)}}{2},$$

وأكبر قيمة صحيحة مقبولة لـ \(y\) هي \(\lfloor\sqrt{u_{\max}(x)}\rfloor\)، مع تصحيح صحيح صغير بسبب كون المتراجحة صارمة.

مثال محسوب: لماذا مجال البحث منتهٍ

لننظر إلى محور \(x\)، أي الحالة \(y=0\). الخروج عن القطع الناقص يعني \(x^2\gt a^2\). لنضع

$$z=x^2-a^2.$$

عندئذ تصبح متراجحة \(45^\circ\)

$$4b^2z \gt (b^2-z)^2,$$

أي

$$z^2-6b^2z+b^4 \lt 0.$$

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

$$b^2(3-2\sqrt{2}) \lt z \lt b^2(3+2\sqrt{2}).$$

وبوجه خاص، كل نقطة محورية صالحة تحقق \(z\lt 6b^2\)، وبالتالي

$$x^2 \lt a^2+6b^2.$$

والحجة نفسها على محور \(y\) تعطي \(y^2\lt b^2+6a^2\). ولهذا السبب لا يفحص البرنامج إلا عددًا منتهيًا من الأعمدة والصفوف، ويستعمل هذه الحدود المرفوعة قليلًا كحدود بحث آمنة.

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

التهيئة الهندسية

تبدأ تطبيقات C++ وPython وJava باختزال معطيات المسألة إلى الثوابت \(a^2\) و\(b^2\) و\(a^2+b^2\) و\(a^2b^2\). بعد ذلك تصبح كل مفاضلة هندسية مجرد مقارنة صحيحة في الإحداثيات المتمركزة.

المسح عمودًا عمودًا

لكل عدد صحيح غير سالب \(x\)، تحسب الشيفرة أول قيمة \(y\) تقع خارج القطع الناقص. ثم تضيف مباشرة جميع النقاط حتى دائرة الدليل، لأن تلك النقاط تحقق شرط الزاوية تلقائيًا. وإذا امتد العمود بعد دائرة الدليل، فإنها تحل المتراجحة التربيعية في \(u=y^2\)، وتحدد آخر قيمة مسموحة لـ \(y\) في ذلك العمود، ثم تضيف بقية النقاط.

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

حساب صحيح دقيق وتصحيح التناظر

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

إذا رمزنا بعدد النقاط الصحيحة المقبولة في الربع المتمركز \(x\ge 0\) و\(y\ge 0\) بالرمز \(N_{Q1}\)، فإن العدد الكلي يُستعاد من الصيغة

$$N=4N_{Q1}-2N_x-2N_y+N_0,$$

حيث يعدّ \(N_x\) و\(N_y\) النقاط المقبولة على المحورين المتمركزين، بينما \(N_0\) هو تصحيح الأصل. تطبق نسختا C++ وJava مسح الأعمدة على نحو متوازٍ، أما نسخة Python فتستخدم الرياضيات نفسها إما تسلسليًا أو بتقسيم قائم على العمليات بحسب بيئة التشغيل.

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

لنرمز إلى حد المسح الأفقي الآمن بـ \(X_{\max}\)، وهو من رتبة \(\sqrt{a^2+6b^2}\). تعالج الحلقة الرئيسية جميع القيم الصحيحة لـ \(x\) من \(0\) حتى \(X_{\max}\)، وكل عمود يحتاج فقط إلى عدد ثابت من العمليات الحسابية مع بضعة تصحيحات حدودية قصيرة. لذلك فإن زمن التنفيذ خطي تقريبًا في عرض المجال الممسوح.

استهلاك الذاكرة هو \(O(1)\) عدا عدد قليل من العدادات والثوابت الهندسية، أو \(O(T)\) إذا احتسبنا المجاميع المحلية الخاصة بـ \(T\) من العمال أو الخيوط. ومصدر الكفاءة الحقيقي هو أن الخوارزمية لا تبحث عن نقاط التماس نفسها، بل تعدّ نقاط الشبكة مباشرة من صيغ مغلقة.

هوامش ومراجع

  1. صفحة المسألة: https://projecteuler.net/problem=246
  2. Ellipse: Wikipedia - Ellipse
  3. Director circle: Wikipedia - Director circle
  4. Tangent: Wikipedia - Tangent
  5. Quadratic equation: Wikipedia - Quadratic equation