Mirrors from conic sections

The conic sections focus. A ray from one focus of an ellipse, reflected by the ellipse, passes through the other focus, and every such path has the same length \(2a\). In an elliptical room, a whisper at one focus is heard at the other. A parabola reflects every ray parallel to its axis into its focus, and a hyperbola reflects every ray aimed at one of its foci toward the other. The Cassegrain telescope combines the two, a parabolic primary mirror and a small convex hyperbolic secondary in front of its focus, which sends the light back through a hole in the primary to a focus behind it.

All three curves are exact rational B-splines in planex2d. The example traces rays through a whispering gallery and a Cassegrain telescope. Every reflected ray must reach its focus, every path must be equally long, and since planex2d intersects a spline on its approximation, the misses must shrink with the spline tolerance.

import math
import numpy as np
import planex as px

The conics as B-splines

A rational quadratic B-spline traces a conic exactly. Between two points, with its middle pole where the tangents there meet, the middle weight sets the kind. Below 1 it is an ellipse, at 1 a parabola, and above 1 a hyperbola. A quarter of the ellipse with the semi-axes \(a\) and \(b\) has the poles \((a, 0)\), \((a, b)\), \((0, b)\) and the middle weight \(\cos 45^\circ\). Four quarters make the whole ellipse. A parabola needs no weights. The branch of the hyperbola \((y - y_0)^2/a^2 - x^2/b^2 = 1\) between \(x = \pm b \sinh u\) has its middle pole at \(y_0 + a/\cosh u\) and the middle weight \(\cosh u\):

def ellipse(a, b):
    """the whole ellipse around the origin, four quarters with the middle weights cos 45 degrees"""
    poles = [(a, 0), (a, b), (0, b), (-a, b), (-a, 0), (-a, -b), (0, -b), (a, -b), (a, 0)]
    w = math.cos(math.radians(45))
    return px.BSpline(poles, [0, 0, 0, 1, 1, 2, 2, 3, 3, 4, 4, 4], 2, weights=[1, w, 1, w, 1, w, 1, w, 1])


def parabola(f, x):
    """y = h^2 / (4 f) for |h| <= x, its focus at (0, f)"""
    return px.BSpline([(-x, x * x / (4 * f)), (0, -x * x / (4 * f)), (x, x * x / (4 * f))], [0, 0, 0, 1, 1, 1], 2)


def hyperbola(y0, a, b, x):
    """the upper branch of (v - y0)^2 / a^2 - h^2 / b^2 = 1 for |h| <= x"""
    u = math.asinh(x / b)
    ends = [(-x, y0 + a * math.cosh(u)), (x, y0 + a * math.cosh(u))]
    return px.BSpline([ends[0], (0, y0 + a / math.cosh(u)), ends[1]], [0, 0, 0, 1, 1, 1], 2, weights=[1, math.cosh(u), 1])


room = ellipse(5, 3)
print("ellipse: one edge, area", room.Face().mass, "exact pi a b", math.pi * 5 * 3)

ellipse: one edge, area 47.12388980384688 exact pi a b 47.12388980384689

Rays

trace follows a ray from a point along a direction to the nearest point where it meets a mirror, and reflects its direction about the mirror’s tangent there:

def trace(mirror, p, d):
    """the point where the ray from p along d first meets the mirror, and the reflected direction"""
    ray = px.Line(p, (p[0] + 100 * d[0], p[1] + 100 * d[1]))
    q = min(((h.x, h.y) for h in mirror.Intersections(ray) if math.dist((h.x, h.y), p) > 1e-9),
            key=lambda h: math.dist(h, p))
    t = mirror.edges[0].edge.Tangent(q)
    t = (t.x / math.hypot(t.x, t.y), t.y / math.hypot(t.x, t.y))
    dot = d[0] * t[0] + d[1] * t[1]
    return q, (2 * dot * t[0] - d[0], 2 * dot * t[1] - d[1])


def miss(point, p, d):
    """how far point lies from the line through p along d"""
    return abs((point[0] - p[0]) * d[1] - (point[1] - p[1]) * d[0]) / math.hypot(*d)

The whispering gallery

The room is an ellipse with the semi-axes 5 and 3, its foci at \((\pm 4, 0)\). Thirty-six rays leave the left focus every 10°. Every reflected ray has to pass through the right focus, and every path from focus to focus has to be \(2a = 10\) long:

F1, F2 = (-4.0, 0.0), (4.0, 0.0)
gallery = []
for k in range(36):
    angle = math.radians(5 + 10 * k)
    q, out = trace(room, F1, (math.cos(angle), math.sin(angle)))
    gallery.append((q, out))
print("largest miss of the right focus", max(miss(F2, q, out) for q, out in gallery))
print("largest |F1 P + P F2 - 2a|", max(abs(math.dist(F1, q) + math.dist(q, F2) - 10) for q, out in gallery))

largest miss of the right focus 1.636441753938289e-05
largest |F1 P + P F2 - 2a| 2.201990458416958e-05

The Cassegrain telescope

The primary is the parabola \(y = x^2/40\), 12 wide, its focus at (0, 10). The secondary is a branch of the hyperbola with the foci (0, 10) and (0, −2), its center at (0, 4), the semi-axis \(a = 4\) and \(b = \sqrt{6^2 - 4^2}\). It is 3 wide, with its vertex at (0, 8), in front of the primary’s focus. Light comes down parallel to the axis, from the height 20, and the secondary shades the middle, so the rays start between 2 and 5.8 from the axis. Every ray has to reach the focus (0, −2) behind the primary, and every path from the height 20 to that focus has to be equally long:

primary = parabola(10, 6)
secondary = hyperbola(4, 4, math.sqrt(6**2 - 4**2), 1.5)
focus = (0.0, -2.0)
telescope = []
for x in np.linspace(2, 5.8, 12):
    p0 = (float(x), 20.0)
    q1, d1 = trace(primary, p0, (0.0, -1.0))
    q2, d2 = trace(secondary, q1, d1)
    telescope.append((p0, q1, q2, d2))
print("largest miss of the focus", max(miss(focus, q2, d2) for p0, q1, q2, d2 in telescope))
paths = [math.dist(p0, q1) + math.dist(q1, q2) + math.dist(q2, focus) for p0, q1, q2, d2 in telescope]
print("paths from", min(paths), "to", max(paths))

largest miss of the focus 4.932002495294797e-06
paths from 37.999993609433 to 38.00000413898293

The next figure draws the rays of both mirrors, the telescope’s on both sides of its axis:

The spline tolerance

The misses are small but not zero. The tangent is that of the exact curve, but planex2d finds the point where a ray meets a spline on its arc approximation, which stays within px.config.spline_tol of the exact curve (see Tolerances). As the tolerance tightens, the reflected rays come closer to the focus and the paths closer to \(2a\):

default = px.config.spline_tol
for tol in (1e-4, 1e-6, 1e-8):
    px.config.spline_tol = tol
    room = ellipse(5, 3)
    traced = [trace(room, F1, (math.cos(math.radians(5 + 10 * k)), math.sin(math.radians(5 + 10 * k)))) for k in range(36)]
    print("spline_tol", tol, "| largest miss", max(miss(F2, q, out) for q, out in traced),
          "| largest path error", max(abs(math.dist(F1, q) + math.dist(q, F2) - 10) for q, out in traced))
px.config.spline_tol = default

spline_tol 0.0001 | largest miss 0.0003318395171638677 | largest path error 0.0005573403676706334
spline_tol 1e-06 | largest miss 1.636441753938289e-05 | largest path error 2.201990458416958e-05
spline_tol 1e-08 | largest miss 1.0664983163266585e-07 | largest path error 7.244408628537258e-08