The Uranus Problem

Le Verrier’s three memoirs and the prediction of Neptune’s position.

A Planet That Wouldn’t Behave

William Herschel discovered Uranus on 13 March 1781. Within decades, astronomers noticed that the planet refused to follow the orbit predicted by Newtonian theory — even after accounting for Jupiter and Saturn’s gravitational pull. By the 1840s, the discrepancy had grown to over 100 arcseconds, far beyond observational error.

In 1845, François Arago — director of the Paris Observatory — pointed Le Verrier toward this problem. Could the anomaly be explained by an unknown planet beyond Uranus?

The First Memoir (10 November 1845)

Le Verrier’s first memoir established the negative result: no adjustment to the known planets’ masses or orbits could eliminate the Uranus residuals. Something else was pulling on Uranus.

We can reproduce this finding with our simulation. Using JPL state vectors for Jupiter, Saturn, and Uranus — a modernization of Le Verrier’s starting point — we integrate without Neptune and compare against observed longitudes:

NoteA note on methodology

Our reconstruction uses JPL state vectors and a numerical integrator where Le Verrier used hand-built planetary tables and analytical perturbation theory. His “baseline” was noisier than ours: he had to manually subtract Jupiter and Saturn’s perturbations using his own imperfect tables before isolating the Uranus anomaly. Our simulation does this “for free.” The results are qualitatively the same, but our residuals are cleaner than what Le Verrier had to work with.

from pathlib import Path
import numpy as np
import pandas as pd

from discoverneptune.simulation import (
    StateVector,
    simulate_uranus_longitudes_from_vectors,
    simulate_uranus_longitudes_no_neptune,
)

# Load bundled observation data
from discoverneptune.data import find_bundled_data_dir

data_dir = find_bundled_data_dir()

obs = pd.read_csv(data_dir / "uranus_observations_1781_1846.csv")
svs = pd.read_csv(data_dir / "outer_planet_state_vectors_1781_01_01.csv")
state = {
    row.planet: StateVector(row.x, row.y, row.z, row.vx, row.vy, row.vz)
    for row in svs.itertuples(index=False)
}

jd_times = obs["datetime_jd"].to_numpy()
truth_lon = obs["longitude_rad"].to_numpy()
start_jd = float(jd_times[0])

# Simulate WITHOUT Neptune
model_no_neptune = simulate_uranus_longitudes_no_neptune(
    jd_times, start_jd, state["jupiter"], state["saturn"], state["uranus"]
)
residuals_no_neptune = np.degrees(np.unwrap(truth_lon - model_no_neptune)) * 3600

print(f"Residual RMS without Neptune: {np.sqrt(np.mean(residuals_no_neptune**2)):.1f} arcsec")
print(f"Peak-to-peak: {np.ptp(residuals_no_neptune):.1f} arcsec")

Residual RMS without Neptune: 31.5 arcsec
Peak-to-peak: 131.4 arcsec

This is the signal Le Verrier saw — over 100 arcseconds of unexplained drift. No tweak to the known planets could make it go away.

The Second Memoir (1 June 1846)

Six months later, Le Verrier presented his mathematical framework for the inverse problem: given the pattern of residuals, what mass, distance, and position could an unknown planet have? He assumed the unknown planet’s orbit was roughly circular and coplanar with the ecliptic — reasonable starting assumptions that reduced the search space.

His approach was a 19th-century version of what we now call optimization: systematically adjusting the unknown planet’s parameters to minimize the residuals.

from discoverneptune.historical_values import (
    LEVERRIER_MEMOIR_1846_JUN,
    LEVERRIER_MEMOIR_1846_AUG,
    MODERN_NEPTUNE,
    GALLE_OBSERVATION,
)

# Le Verrier's June 1846 prediction
print("Le Verrier's June 1846 prediction:")
print(f"  Semi-major axis: {LEVERRIER_MEMOIR_1846_JUN.semi_major_axis_au:.3f} AU")
print(f"  Mass: {LEVERRIER_MEMOIR_1846_JUN.mass_solar:.6f} solar masses")
print(f"  Longitude: {LEVERRIER_MEMOIR_1846_JUN.longitude_deg:.1f}\u00b0")

Le Verrier's June 1846 prediction:
  Semi-major axis: 36.154 AU
  Mass: 0.000108 solar masses
  Longitude: 325.0°

The Third Memoir and the Letter to Galle (31 August 1846)

Le Verrier’s third memoir refined the prediction. He wrote to Johann Galle at the Berlin Observatory on 18 September 1846 with specific coordinates. On the night of 23 September, Galle and his student Heinrich d’Arrest found the new planet within 1° of Le Verrier’s predicted position.

print("Le Verrier's final prediction (August 1846):")
print(
    "  Published longitude: "
    f"{LEVERRIER_MEMOIR_1846_AUG.longitude_deg:.2f}\u00b0 "
    "(quoted for 1 Jan 1847)"
)
print(f"  Predicted semi-major axis: {LEVERRIER_MEMOIR_1846_AUG.semi_major_axis_au:.3f} AU")
print()
print("Galle's observation (23 September 1846):")
print(f"  Observed longitude: {GALLE_OBSERVATION.longitude_deg:.1f}\u00b0")
print()
error_deg = abs(LEVERRIER_MEMOIR_1846_AUG.longitude_deg - GALLE_OBSERVATION.longitude_deg)
print(f"Raw longitude offset: {error_deg:.2f}\u00b0")
print("  Note: the memoir longitude is quoted for 1 Jan 1847,")
print("  so this is only a rough historical comparison.")

Le Verrier's final prediction (August 1846):
  Published longitude: 326.00° (quoted for 1 Jan 1847)
  Predicted semi-major axis: 36.154 AU

Galle's observation (23 September 1846):
  Observed longitude: 328.0°

Raw longitude offset: 2.00°
  Note: the memoir longitude is quoted for 1 Jan 1847,
  so this is only a rough historical comparison.

Our Modernized Search

We can run the same kind of search — a modernized analogue, not a literal recreation — using differential evolution to find Neptune’s parameters from the Uranus residuals. The key difference: we hold Jupiter, Saturn, and Uranus at their JPL truth positions and optimize only Neptune’s 4 parameters (mass, semi-major axis, eccentricity, mean longitude).

from discoverneptune.solver import solve

result = solve(
    jd_times, truth_lon, start_jd,
    state["jupiter"], state["saturn"], state["uranus"],
    de_maxiter=30, de_popsize=5, polish=True,
)
rms_arcsec = np.degrees(np.sqrt(result.sse / len(jd_times))) * 3600
print(f"Optimized Neptune parameters:")
print(f"  Semi-major axis: {result.neptune.a:.2f} AU (modern: {MODERN_NEPTUNE.semi_major_axis_au:.3f} AU)")
print(f"  Mass: {result.neptune.mass_solar:.6f} solar masses")
print(f"  Eccentricity: {result.neptune.e:.4f}")
print(f"  RMS residual: {rms_arcsec:.2f} arcsec")

Optimized Neptune parameters:
  Semi-major axis: 30.14 AU (modern: 30.069 AU)
  Mass: 0.000044 solar masses
  Eccentricity: 0.0645
  RMS residual: 2.87 arcsec

The discovery in one figure

Plotting the residuals before and after the Neptune fit tells the story in a single panel. Without Neptune, even with Jupiter, Saturn, and Uranus pinned at their JPL truth positions, Uranus’s longitude drifts away from observation by ~130 arcsec across the historical record. Adding one unseen planet — fit by four parameters — collapses that drift to a flat ~2.5-arcsec band that sits at the forward-model noise floor.

import matplotlib.pyplot as plt

sim_with_neptune = simulate_uranus_longitudes_from_vectors(
    jd_times, start_jd,
    state["jupiter"], state["saturn"], state["uranus"],
    result.neptune,
)

ARCSEC_PER_RAD = 3600.0 * 180.0 / np.pi
wrap = lambda diff: ((diff + np.pi) % (2 * np.pi)) - np.pi
resid_no = wrap(truth_lon - model_no_neptune) * ARCSEC_PER_RAD
resid_with = wrap(truth_lon - sim_with_neptune) * ARCSEC_PER_RAD
years = 1781.0 + (jd_times - start_jd) / 365.25


def _build(theme="light"):
    apply_style(theme)
    pal = palette(theme)
    fig, ax = plt.subplots(figsize=(9.5, 4))
    ax.plot(years, resid_no, "o-", lw=0.8, ms=3, color=pal.before,
            label="No Neptune — the anomaly")
    ax.plot(years, resid_with, "o-", lw=0.8, ms=3, color=pal.after,
            label=f"With fitted Neptune — RMS {rms_arcsec:.1f}\u2033")
    ax.axhline(0, color=pal.spine, lw=0.5, ls=":")
    ax.set_xlabel("Year")
    ax.set_ylabel("Uranus longitude residual (arcsec)")
    ax.set_title("Le Verrier's problem, solved")
    ax.legend(loc="lower left")
    return fig


display(dual_render(_build, alt="Uranus longitude residuals before and after Neptune fit"))

Figure 1: Uranus longitude residuals versus time under the 4-body no-Neptune model (red, the anomaly Le Verrier saw) and under the 5-body model with the fitted Neptune (black, flat at the ~2.5-arcsec forward-model noise floor). Residuals are angle-wrapped so no detrending or unwrapping is needed.

The Mass-Distance Degeneracy

Le Verrier got the longitude almost exactly right but overestimated the distance (36.2 AU predicted vs 30.1 AU actual). This isn’t a failure — it reflects a fundamental degeneracy: a more massive planet farther away produces nearly the same perturbation pattern as a less massive planet closer in.

The longitude is better constrained because it determines where in the sky the perturbation signal points. Distance and mass trade off against each other along a valley of nearly-equal-fit solutions.

In the figure below, note that our optimizer’s propagated longitude sits a few degrees farther from Galle than Le Verrier’s — that gap is not a worse fit but the ω = Ω = 0 convention bias described in the caption, while the fitted distance (≈30.1 AU) is nearly exact where Le Verrier’s was 20% high: the two methods land on different points of the same degeneracy valley.

import matplotlib.pyplot as plt

from astropy.time import Time

from discoverneptune.simulation import build_simulation_from_vectors

# Put every marker on the discovery night, 23 Sep 1846.
#
# (a) Our optimizer: result.neptune.l_rad is the mean longitude at the 1781
#     epoch; integrate the fitted orbit forward to the discovery date.
sim_opt = build_simulation_from_vectors(
    state["jupiter"], state["saturn"], state["uranus"], result.neptune
)
t_discovery_yr = (Time("1846-09-23").jd - start_jd) / 365.25
sim_opt.integrate(t_discovery_yr)
_sun = sim_opt.particles[0]
_nep = sim_opt.particles[-1]
opt_lon_1846_deg = float(
    np.degrees(np.arctan2(_nep.y - _sun.y, _nep.x - _sun.x)) % 360.0
)

# (b) Le Verrier's memoir longitudes are quoted for 1 Jan 1847
#     (see src/discoverneptune/historical_values.py). Walk them back to the
#     discovery night with each prediction's own Keplerian mean motion.
_dt_yr = (Time("1846-09-23").jd - Time("1847-01-01").jd) / 365.25  # ≈ −0.274
def _to_discovery_night(memoir):
    n_deg_per_yr = 360.0 / memoir.semi_major_axis_au**1.5
    return (memoir.longitude_deg + n_deg_per_yr * _dt_yr) % 360.0

lv_jun_1846_deg = _to_discovery_night(LEVERRIER_MEMOIR_1846_JUN)
lv_aug_1846_deg = _to_discovery_night(LEVERRIER_MEMOIR_1846_AUG)

labels = ["Le Verrier\n(Jun 1846)", "Le Verrier\n(Aug 1846)", "Galle\n(observed)", "Our optimizer", "Modern\n(JPL)"]
a_values = [
    LEVERRIER_MEMOIR_1846_JUN.semi_major_axis_au,
    LEVERRIER_MEMOIR_1846_AUG.semi_major_axis_au,
    GALLE_OBSERVATION.semi_major_axis_au,
    result.neptune.a,
    MODERN_NEPTUNE.semi_major_axis_au,
]
lon_values = [
    lv_jun_1846_deg,
    lv_aug_1846_deg,
    GALLE_OBSERVATION.longitude_deg,
    opt_lon_1846_deg,
]
lon_labels = ["Le Verrier\n(Jun 1846)", "Le Verrier\n(Aug 1846)", "Galle\n(observed)", "Our optimizer\n(propagated)"]


def _build(theme="light"):
    apply_style(theme)
    cycle = category_cycle(theme)
    # Map: Jun/Aug predictions → cycle[0], Galle observation → cycle[1],
    # our optimizer → cycle[2], modern JPL → cycle[3]
    bar_colors_a = [cycle[0], cycle[0], cycle[1], cycle[2], cycle[3]]
    bar_colors_lon = [cycle[0], cycle[0], cycle[1], cycle[2]]
    fig, axes = plt.subplots(1, 2, figsize=(12, 4))

    ax = axes[0]
    ax.barh(labels, a_values, color=bar_colors_a)
    ax.set_xlabel("Semi-major axis (AU)")
    ax.set_title("Distance estimates")
    truth_line(ax, MODERN_NEPTUNE.semi_major_axis_au, label="Modern value",
               axis="x", theme=theme)

    ax = axes[1]
    for y_pos, (lon, col) in enumerate(
        zip(lon_values, bar_colors_lon, strict=True)
    ):
        ax.plot([lon], [y_pos], "D", ms=11, color=col)
    ax.set_yticks(range(len(lon_labels)))
    ax.set_yticklabels(lon_labels)
    ax.set_xlim(min(lon_values) - 3, max(lon_values) + 3)
    ax.set_xlabel("Heliocentric longitude on 23 Sep 1846 (degrees)")
    ax.set_title("Position estimates")

    return fig


display(dual_render(_build, alt="Comparison of Le Verrier predictions, Galle observation, and modern Neptune"))

Figure 2: Position estimates propagated to the discovery night (23 Sep 1846): Le Verrier’s memoir longitudes walked back from their 1 Jan 1847 epoch via their own mean motions, our fitted orbit integrated forward from 1781, and Galle’s observed position. The optimizer’s residual offset (≈9°) is the longitude bias of the simplified model: with ω = Ω = 0 fixed, part of the orbit geometry is absorbed into the fitted mean longitude (see the JPL-comparison caveat on the landing page).

References

• Le Verrier, U. “Première mémoire sur la théorie d’Uranus.” Comptes Rendus 21 (10 Nov 1845).
• Le Verrier, U. “Recherches sur les mouvements d’Uranus.” Comptes Rendus 22 (1 Jun 1846).
• Le Verrier, U. “Sur la planète qui produit les anomalies observées dans le mouvement d’Uranus.” Comptes Rendus 23 (31 Aug 1846).
• Grosser, M. The Discovery of Neptune (1962), chapters 5–7.
• Standage, T. The Neptune File (2000).