Solvable Evaporation Models

Solvable two-dimensional evaporation models isolate the ingredients behind an island transition: a nongravitating radiation algebra, a semiclassical gravitational entropy cost, matter correlations between radiation and a candidate island, and a competition with the empty island. This page supplies a standard-library Python reference calculation with an analytic QES and Page time. It tests saddle selection and numerical convergence; it is not a simulation of dilaton-gravity field equations.

Required background. Quantum Information supplies joint-region entropy, purity, and the meaning of a Page transition.

Helpful background. Semiclassical Gravity supplies dilaton generalized entropy and backreaction; Quantum Fields in Curved Spacetime supplies the bath entropy; Black-Hole Information supplies the QES and replica-wormhole construction.

A dimensionless evaporation benchmark

Use a dimensionless bath time $t$ and island endpoint coordinate $q\in[0,2]$. The empty-island and nonempty-island candidates are

$$S_{\varnothing}(t)=s_0+r t,$$ $$S_I(t,q) =B_0+\delta-bt +\frac{k}{2}(q-q_*)^2.$$

$B_0+\delta-bt$ represents an already-renormalized dilaton term plus the joint matter entropy of $R\cup I$ after the leading correlations are combined. The positive curvature $k$ makes $q_*$ a stable minimum in this one spatial coordinate. A physical Lorentzian QES must also be stationary in the independent time-normal direction, which the surrogate omits. The defaults are

$$(s_0,r,B_0,\delta,b,k,q_*) =(0.2,1.2,10,0.4,0.3,6,0.733).$$

Analytically,

$$q_{\mathrm{QES}}=q_*, \qquad t_{\mathrm P} =\frac{B_0+\delta-s_0}{r+b} =6.8.$$

The numerical method samples $q$ on a uniform grid, chooses the least island value, and locates its crossing with $S_{\varnothing}$ by bisection to a time-bracket tolerance of $10^{-12}$. The QES residual is $\lvert\partial_qS_I\rvert=k\lvert q-q_*\rvert$. The calculation is deterministic and uses no random seed.

The model preserves the selection structure of explicit JT-gravity bath calculations but replaces their conformal-field-theory entropy functions with a quadratic local benchmark. The physical derivations are in Almheiri, Mahajan, Maldacena, and Zhao 2020, §§2–4 and Almheiri et al. 2020, §§2–5.

Executable reference program

The retained environment is CPython 3.11 using only the Python standard library; no third-party package is imported. IEEE-754 binary64 arithmetic is sufficient because discretization error is much larger than roundoff.

Only deterministic scalar arithmetic, a finite grid search, and bisection enter the result. A rerun on another conforming interpreter should reproduce the printed digits; disagreement above roundoff is therefore a useful signal of a transcription, parameter, or stopping-rule change rather than package drift. The program exits nonzero if an external analytic reference disagrees, the default refinement sequence fails its declared convergence thresholds, or the admissible interval excludes the interior QES.

Save the following program as island_reference.py. The program emits the analytic benchmark and platform-independent columns for grid size, mesh spacing, QES coordinate, stationarity residual, Page time, and absolute Page-time error.

import argparse
from dataclasses import dataclass


@dataclass(frozen=True)
class Model:
    s0: float = 0.2
    rate_r: float = 1.2
    black0: float = 10.0
    state_offset: float = 0.4
    rate_b: float = 0.3
    curvature: float = 6.0
    q_star: float = 0.733


def no_island(t, model):
    return model.s0 + model.rate_r * t


def island(t, q, model):
    displacement = q - model.q_star
    return (model.black0 + model.state_offset - model.rate_b * t
            + 0.5 * model.curvature * displacement * displacement)


def discrete_qes(t, size, q_min, q_max, model):
    step = (q_max - q_min) / (size - 1)
    points = (q_min + index * step for index in range(size))
    q_num = min(points, key=lambda q: island(t, q, model))
    entropy = island(t, q_num, model)
    residual = abs(model.curvature * (q_num - model.q_star))
    return q_num, entropy, residual, step


def page_time(size, q_min, q_max, tolerance, model):
    def difference(t):
        _, entropy, _, _ = discrete_qes(t, size, q_min, q_max, model)
        return entropy - no_island(t, model)

    lower, upper = 0.0, 12.0
    f_lower, f_upper = difference(lower), difference(upper)
    if f_lower * f_upper > 0.0:
        raise ValueError("time bracket contains no Page crossing")
    for _ in range(200):
        middle = 0.5 * (lower + upper)
        f_middle = difference(middle)
        if upper - lower <= tolerance:
            break
        if f_lower * f_middle <= 0.0:
            upper, f_upper = middle, f_middle
        else:
            lower, f_lower = middle, f_middle
    return 0.5 * (lower + upper)


def main():
    parser = argparse.ArgumentParser()
    parser.add_argument("--state-offset", type=float, default=0.4)
    parser.add_argument("--q-min", type=float, default=0.0)
    parser.add_argument("--q-max", type=float, default=2.0)
    parser.add_argument("--expected-page", type=float)
    args = parser.parse_args()
    if not args.q_min < args.q_max:
        raise ValueError("q-min must be smaller than q-max")

    model = Model(state_offset=args.state_offset)
    if not args.q_min < model.q_star < args.q_max:
        raise ValueError("admissible interval excludes the interior QES")

    exact = ((model.black0 + model.state_offset - model.s0)
             / (model.rate_r + model.rate_b))
    reference_tolerance = 1.0e-12
    if (args.expected_page is not None
            and abs(args.expected_page - exact) > reference_tolerance):
        raise RuntimeError("analytic Page benchmark mismatch")

    rows = []
    for size in (41, 81, 161, 321):
        q_num, _, residual, step = discrete_qes(
            exact, size, args.q_min, args.q_max, model
        )
        crossing = page_time(
            size, args.q_min, args.q_max, reference_tolerance, model
        )
        page_error = abs(crossing - exact)
        expected_crossing = (exact
                             + 0.5 * model.curvature
                             * (q_num - model.q_star) ** 2
                             / (model.rate_r + model.rate_b))
        if abs(crossing - expected_crossing) > 2.0e-12:
            raise RuntimeError("numerical and analytic grid benchmarks disagree")
        if abs(q_num - model.q_star) > 0.5 * step + 1.0e-14:
            raise RuntimeError("grid search did not select the nearest QES point")
        rows.append((size, step, q_num, residual, crossing, page_error))

    if args.q_min == 0.0 and args.q_max == 2.0:
        residuals = [row[3] for row in rows]
        page_errors = [row[5] for row in rows]
        residuals_shrink = all(
            later < earlier for earlier, later in zip(residuals, residuals[1:])
        )
        errors_shrink = all(
            later < earlier
            for earlier, later in zip(page_errors, page_errors[1:])
        )
        if (not residuals_shrink or not errors_shrink
                or residuals[-1] > 1.1e-2
                or page_errors[-1] > 7.0e-6):
            raise RuntimeError("default refinement failed convergence thresholds")

    print(f"analytic_qes={model.q_star:.9f} analytic_page={exact:.9f}")
    print("N h qes residual page_time abs_error")
    for size, step, q_num, residual, crossing, page_error in rows:
        print(f"{size:d} {step:.8f} {q_num:.8f} {residual:.8e} "
              f"{crossing:.9f} {page_error:.8e}")


if __name__ == "__main__":
    main()

A clean rerun consists of saving exactly that source and executing these commands:

python3 --version
python3 island_reference.py
python3 island_reference.py --expected-page 6.8
python3 island_reference.py --state-offset 1.0 --expected-page 7.2

The first line should report a CPython 3.11 release. Both default-parameter runs print analytic_qes=0.733000000 analytic_page=6.800000000; the second also checks that value against an independently supplied reference. The final physics run changes only the state-preparation offset and should move the analytic Page time to $7.2$ without changing $q_*$. Both referenced runs must exit zero.

Convergence and independent benchmark

The following values are from a clean default run. The analytic reference is derived above, independently of the grid and bisection implementation.

Grid $N$	Mesh spacing	QES result	Residual	Page time result	Absolute error
41	0.05000000	0.75000000	1.02000000e-1	6.800578000	5.78000000e-4
81	0.02500000	0.72500000	4.80000000e-2	6.800128000	1.28000000e-4
161	0.01250000	0.73750000	2.70000000e-2	6.800040500	4.04999997e-5
321	0.00625000	0.73125000	1.05000000e-2	6.800006125	6.12500032e-6

Both the stationarity residual and Page-time error decrease across four default refinements. The retained program requires strict decrease for this declared sequence, a final residual no larger than $1.1\times10^{-2}$, and a final Page-time error no larger than $7.0\times10^{-6}$; violation gives a nonzero exit. These are tolerances for the declared benchmark, not a claim about physical evaporation accuracy.

There is also an analytic error bound. With mesh spacing $h$ and an interior $q_*$, nearest-grid minimization gives $\lvert q_N-q_*\rvert\le h/2$, hence

$$\begin{aligned} \lvert\partial_q S_I(q_N)\rvert&\le kh/2,\\ 0\le t_{\mathrm P,N}-t_{\mathrm P} &=\frac{k(q_N-q_*)^2}{2(r+b)} \le\frac{kh^2}{8(r+b)}. \end{aligned}$$

The second line excludes the negligible bisection error. Thus the residual has an $O(h)$ upper bound and the Page-time error an $O(h^2)$ upper bound. Actual refinement ratios oscillate as $q_*$ moves relative to the grid; strict decrease at every refinement is not guaranteed for other parameters. This bound explains convergence without inferring an order from four fortunate samples.

The time bisection error is below $10^{-12}$ in the bracket width and is negligible here. The reported Page-time error is instead dominated by the spatial grid’s positive bias in the island minimum. This separation is useful: tightening bisection alone cannot repair an under-resolved QES.

Ingredients, sensitivities, and failure modes

The nongravitating bath makes $S(R)$ a definite subsystem entropy. The nonempty admissible island and its gravitational cost create the second branch. Correlations encoded in the lower late-time island entropy permit the branch exchange. The quadratic profile, uniform grid, and linear time dependence merely simplify evaluation; none by itself explains why an island saddle is allowed.

Changing --state-offset models a different finite state-preparation contribution and shifts the crossing. Changing the allowed $q$ interval models a boundary-condition change. If $q_*$ lies outside that interval, the grid would select an endpoint with a nonzero derivative; that is a constrained boundary minimum, not a QES. The retained program therefore rejects the run before labeling or tabulating the endpoint.

Two adverse commands exercise distinct failure paths:

python3 island_reference.py --q-max 0.70
python3 island_reference.py --expected-page 6.7

The first must exit nonzero because the allowed interval excludes $q_*=0.733$. The second must exit nonzero because the supplied reference disagrees with the independently evaluated analytic value $6.8$ by much more than $10^{-12}$.

Other failure modes are an unbracketed crossing, negative curvature $k$ producing an unstable stationary point, a Page time beyond the model’s weak-backreaction interval, a third omitted saddle, or a Planckian dilaton where semiclassical entropy fails. The reference calculation does not evolve a metric, solve a matter stress tensor, determine a replica contour, or establish a microscopic S-matrix. A maintained interactive gravitational simulation would require a separately reviewed executable artifact and is not implied here.

Exercises

Predict the outputs affected by --state-offset 1.0, and diagnose a run with --q-max 0.70. Which change represents state preparation, and which invalidates the interior QES benchmark?

Solution

Changing the offset from $0.4$ to $1.0$ gives $t_{\mathrm P}=(10+1-0.2)/1.5=7.2$. It leaves the analytic QES $q_*=0.733$, mesh, and residual sequence unchanged because the offset has no $q$ dependence. Setting $q_{\max}=0.70$ excludes the stationary point. A constrained calculation would drive the minimizer to the upper boundary and give residual $6(0.733-0.70)=0.198$; the retained program instead exits nonzero before presenting that endpoint as an interior QES.

References

• Almheiri et al. 2020a, Ahmed Almheiri, Raghu Mahajan, Juan Maldacena, and Ying Zhao, “The Page Curve of Hawking Radiation from Semiclassical Geometry,” Journal of High Energy Physics 03 (2020) 149, DOI: 10.1007/JHEP03(2020)149.
• Almheiri et al. 2020b, Ahmed Almheiri, Thomas Hartman, Juan Maldacena, Edgar Shaghoulian, and Amirhossein Tajdini, “Replica Wormholes and the Entropy of Hawking Radiation,” Journal of High Energy Physics 05 (2020) 013, DOI: 10.1007/JHEP05(2020)013.