generate.py

back to table · edit · history · where entries came from · files · download

21050 bytes, as of the version from 2026-09-18 04:02 (current). Recorded here, not run.

r"""Chern-Simons invariants of the hyperbolic prime knots with at most ten crossings -- numberdb.org/T313

    CS(S^3 \ K),   the Chern-Simons invariant of the complete hyperbolic metric
                   on the complement, normalised as KnotInfo and SnapPy print
                   it: an element of R/(1/2)Z, written as its representative in
                   [-1/4, 1/4),

for every hyperbolic prime knot K = n_k of the Rolfsen table with at most ten
crossings, and for the mirror image of each chiral one. The second
normalisation stored is 2*pi^2 CS, the imaginary part of the complex volume.

Run it with SageMath and SnapPy:

    $ sage -pip install numberdb snappy   # once
    $ sage -python generate.py            # check the table against this code
    $ sage -python generate.py --publish  # send it, with NUMBERDB_API_KEY set

## Where the digits come from

The table first held KnotInfo's nine printed decimals. They are not a limit of
what can be computed: SnapPy's complex volume is available as a *certified
interval* at any working precision, through the extended Bloch group
(Zickert's algorithm) with the shapes bounded by interval Newton. At 400 bits
the interval for a ten-crossing knot has radius under 1e-109, and one knot
costs a tenth of a second, so a hundred digits is not an expense -- it is the
same run.

That matches HREF-free but deliberate company: the sibling table of hyperbolic
volumes (T141) carries a hundred digits of the *real* part of the very same
complex volume, computed from the same braid words. A table holding the
imaginary part to nine decimals beside one holding the real part to a hundred
was an accident of provenance, not a fact about the invariants.

## The 2-torsion, and why the level is `assumed-bound` rather than `proven`

SnapPy's verified complex volume is `verified_modulo_2_torsion`: it determines
Vol + i 2 pi^2 CS modulo i pi^2/2, not modulo the full i pi^2. In the
normalisation used here that means the certified interval pins CS **modulo
1/4**, leaving two candidates in [-1/4, 1/4) that differ by exactly 1/4.

The class is chosen by SnapPy's quad-double value, and then checked twice:

  * the quad-double computation (`Manifold.high_precision().chern_simons()`,
    an independent path through the gluing equations) must agree with the
    chosen candidate to at least 30 digits -- it agrees to about 50;
  * KnotInfo's published nine decimals must agree with it too.

Both would have to be wrong by 0.25 for the choice to be wrong. That is not a
proof, and it is exactly what `assumed-bound` is for: a fixed-precision
computation with a bound asserted for a stated reason, after which everything
is interval arithmetic again. The digits written here follow from the width of
a certified interval; what rests on floating point is *which* number, not how
many of its digits are right.

An amphichiral knot is the one case where the class is proven: CS = -CS in
R/(1/2)Z forces 2 CS = 0, so CS is 0 or 1/4 exactly, and the certified
interval excludes 1/4. Those twenty entries are exactly 0.

## What else is checked, per knot, before an entry is returned

  * the exterior is built from KnotInfo's braid word, not from SnapPy's
    Rolfsen name, which numbers part of the ten-crossing table as Rolfsen did
    before Perko;
  * the triangulation is geometric (all tetrahedra positively oriented);
  * the real part of the same complex volume agrees with SnapPy's quad-double
    volume to 30 digits, which is the check that the manifold is the one the
    volume table (T141) also computed;
  * for a chiral knot, the mirror exterior is computed separately and
    CS(mirror) + CS(K) must contain 0 as an interval -- the congruence
    CS(S^3 \ bar K) = -CS(S^3 \ K) (mod 1/2) at a hundred digits rather than
    at nine;
  * the certified interval must support the hundred digits written, with a
    margin, or the knot raises rather than writing fewer.
"""

import os
import sys
from decimal import Decimal

import numberdb.sage as numberdb
from numberdb._write import to_text
#The knot table imports sage.functions, which cannot initialise the symbolic
#ring from inside its own import; brought up first, by name, it can.
import sage.symbolic.ring  # noqa: F401
from sage.rings.integer_ring import ZZ
from sage.rings.real_arb import RealBallField
from snappy import Link

from knotinfo_prime_knots import (
    AMPHICHIRAL,
    KNOTINFO_BRAIDS,
    NAMED,
    TORUS,
    is_alternating,
    names,
    subscript,
)


TABLE = os.environ.get("NUMBERDB_TABLE", "T313")

#: Digits written, the corpus's usual hundred and the sibling volume table's.
DIGITS = 100

#: Working precision for the certified computation. Measured: at 400 bits the
#: interval radius is 1e-115 for the figure-eight (2 tetrahedra) and 1e-109 for
#: 10_123 (18 tetrahedra), so a hundred digits has nine to spare on the worst
#: knot in the table. A knot costs about a tenth of a second at this precision.
BITS = 400

#: How far SnapPy's quad-double value must agree with the certified interval
#: before its choice of 2-torsion class is accepted. It agrees to about 50; the
#: wrong class is 0.25 away, so 30 is already an enormous margin, and a value
#: near zero is compared absolutely for the same reason.
CLASS_DIGITS = 30

GEOMETRIC = "all tetrahedra positively oriented"

KNOTINFO_CS = {
    (4, 1): "0",
    (5, 2): "-0.153204133",
    (6, 1): "0.155977017",
    (6, 2): "-0.202492498",
    (6, 3): "0",
    (7, 2): "-0.055153535",
    (7, 3): "0.187220178",
    (7, 4): "0.021726694",
    (7, 5): "0.120555869",
    (7, 6): "0.182283191",
    (7, 7): "0.132985600",
    (8, 1): "0.222133121",
    (8, 2): "0.212918515",
    (8, 3): "0",
    (8, 4): "-0.181952812",
    (8, 5): "0.182078283",
    (8, 6): "-0.072534896",
    (8, 7): "0.131958556",
    (8, 8): "-0.117099372",
    (8, 9): "0",
    (8, 10): "0.157052315",
    (8, 11): "0.019309343",
    (8, 12): "0",
    (8, 13): "-0.192783664",
    (8, 14): "-0.051993050",
    (8, 15): "-0.018424716",
    (8, 16): "0.151381673",
    (8, 17): "0",
    (8, 18): "0",
    (8, 20): "0.103363447",
    (8, 21): "-0.241804603",
    (9, 2): "-0.007707448",
    (9, 3): "0.116481613",
    (9, 4): "-0.207615097",
    (9, 5): "0.100298381",
    (9, 6): "0.024287061",
    (9, 7): "0.219457098",
    (9, 8): "0.097141668",
    (9, 9): "-0.017764270",
    (9, 10): "-0.089574066",
    (9, 11): "-0.063062742",
    (9, 12): "-0.172870946",
    (9, 13): "-0.197081974",
    (9, 14): "0.006939989",
    (9, 15): "0.197391343",
    (9, 16): "-0.057237415",
    (9, 17): "0.035265127",
    (9, 18): "-0.166289356",
    (9, 19): "0.247414919",
    (9, 20): "0.007626177",
    (9, 21): "0.120109098",
    (9, 22): "-0.010955728",
    (9, 23): "-0.242331920",
    (9, 24): "-0.121749639",
    (9, 25): "-0.219499139",
    (9, 26): "-0.173837793",
    (9, 27): "-0.196169915",
    (9, 28): "0.195408234",
    (9, 29): "-0.029441592",
    (9, 30): "0.158418513",
    (9, 31): "0.134608909",
    (9, 32): "-0.212531774",
    (9, 33): "0.140691430",
    (9, 34): "0.167756984",
    (9, 35): "0.109912886",
    (9, 36): "-0.056326575",
    (9, 37): "-0.199700696",
    (9, 38): "-0.164704133",
    (9, 39): "0.184473921",
    (9, 40): "0.202815317",
    (9, 41): "0.050000000",
    (9, 42): "0.103415652",
    (9, 43): "-0.218003228",
    (9, 44): "0.023833924",
    (9, 45): "-0.101005115",
    (9, 46): "-0.145047960",
    (9, 47): "-0.085168050",
    (9, 48): "-0.102810873",
    (9, 49): "0.066236000",
    (10, 1): "-0.242222319",
    (10, 2): "0.169666840",
    (10, 3): "0.072170674",
    (10, 4): "-0.126209321",
    (10, 5): "0.193205657",
    (10, 6): "-0.173623346",
    (10, 7): "0.116340907",
    (10, 8): "0.106309910",
    (10, 9): "-0.088973251",
    (10, 10): "0.222320933",
    (10, 11): "-0.219389223",
    (10, 12): "0.031490227",
    (10, 13): "0.139343184",
    (10, 14): "-0.194044237",
    (10, 15): "0.229298725",
    (10, 16): "-0.101238697",
    (10, 17): "0",
    (10, 18): "-0.181673094",
    (10, 19): "0.183855815",
    (10, 20): "0.000141354",
    (10, 21): "-0.064632260",
    (10, 22): "0.096051141",
    (10, 23): "-0.090530220",
    (10, 24): "0.134075911",
    (10, 25): "-0.190643197",
    (10, 26): "0.214239473",
    (10, 27): "0.023328145",
    (10, 28): "0.219341461",
    (10, 29): "-0.179148898",
    (10, 30): "0.123862719",
    (10, 31): "-0.080975515",
    (10, 32): "0.136232538",
    (10, 33): "0",
    (10, 34): "-0.179187486",
    (10, 35): "-0.091648377",
    (10, 36): "0.033820488",
    (10, 37): "0",
    (10, 38): "0.053142588",
    (10, 39): "-0.152342696",
    (10, 40): "-0.009943898",
    (10, 41): "-0.111474841",
    (10, 42): "-0.054167761",
    (10, 43): "0",
    (10, 44): "-0.166712665",
    (10, 45): "0",
    (10, 46): "0.117601502",
    (10, 47): "0.233824120",
    (10, 48): "-0.030368616",
    (10, 49): "-0.142813363",
    (10, 50): "-0.110994882",
    (10, 51): "-0.052808895",
    (10, 52): "0.169102971",
    (10, 53): "0.164530375",
    (10, 54): "0.226255147",
    (10, 55): "0.070960601",
    (10, 56): "-0.206716780",
    (10, 57): "0.034925503",
    (10, 58): "0.069592985",
    (10, 59): "-0.135970795",
    (10, 60): "-0.199702174",
    (10, 61): "0.107287730",
    (10, 62): "-0.243566962",
    (10, 63): "0.097962685",
    (10, 64): "-0.142858798",
    (10, 65): "-0.022363657",
    (10, 66): "-0.196684628",
    (10, 67): "0.091195441",
    (10, 68): "0.185277564",
    (10, 69): "-0.013618542",
    (10, 70): "0.181592128",
    (10, 71): "-0.002720688",
    (10, 72): "-0.241280314",
    (10, 73): "-0.060830498",
    (10, 74): "0.203041242",
    (10, 75): "-0.155784059",
    (10, 76): "0.223370687",
    (10, 77): "0.094005380",
    (10, 78): "-0.088623060",
    (10, 79): "0",
    (10, 80): "-0.164435070",
    (10, 81): "0",
    (10, 82): "-0.104004668",
    (10, 83): "-0.039387087",
    (10, 84): "-0.064360610",
    (10, 85): "0.239779044",
    (10, 86): "0.183147375",
    (10, 87): "0.076759058",
    (10, 88): "0",
    (10, 89): "-0.073918755",
    (10, 90): "0.146644036",
    (10, 91): "-0.010769866",
    (10, 92): "-0.179635306",
    (10, 93): "-0.227496792",
    (10, 94): "-0.122288498",
    (10, 95): "-0.044751139",
    (10, 96): "-0.188546503",
    (10, 97): "0.109616979",
    (10, 98): "-0.144138872",
    (10, 99): "0",
    (10, 100): "-0.246476660",
    (10, 101): "0.095028752",
    (10, 102): "0.117069293",
    (10, 103): "0.002156507",
    (10, 104): "-0.009383138",
    (10, 105): "-0.171725261",
    (10, 106): "-0.126045750",
    (10, 107): "-0.035682646",
    (10, 108): "-0.241341556",
    (10, 109): "0",
    (10, 110): "-0.138072510",
    (10, 111): "-0.117089254",
    (10, 112): "0.117723670",
    (10, 113): "-0.046145975",
    (10, 114): "-0.142409299",
    (10, 115): "0",
    (10, 116): "0.126086904",
    (10, 117): "0.018788434",
    (10, 118): "0",
    (10, 119): "-0.149213537",
    (10, 120): "0.159383820",
    (10, 121): "0.030415577",
    (10, 122): "-0.096862024",
    (10, 123): "0",
    (10, 125): "-0.072814255",
    (10, 126): "0.230410569",
    (10, 127): "0.124073851",
    (10, 128): "0.220275737",
    (10, 129): "0.191263345",
    (10, 130): "-0.021262717",
    (10, 131): "-0.084977695",
    (10, 132): "0.186748986",
    (10, 133): "-0.211510648",
    (10, 134): "0.067235354",
    (10, 135): "-0.200777540",
    (10, 136): "0.247290630",
    (10, 137): "-0.147855661",
    (10, 138): "0.060735926",
    (10, 139): "-0.228927561",
    (10, 140): "-0.103360013",
    (10, 141): "-0.175623698",
    (10, 142): "0.161733527",
    (10, 143): "-0.224020127",
    (10, 144): "0.071725956",
    (10, 145): "0.228447210",
    (10, 146): "-0.173581855",
    (10, 147): "0.119232019",
    (10, 148): "-0.214138314",
    (10, 149): "0.076427941",
    (10, 150): "0.140030728",
    (10, 151): "0.241018368",
    (10, 152): "0.240454021",
    (10, 153): "0.226641032",
    (10, 154): "0.240380161",
    (10, 155): "0.177006746",
    (10, 156): "0.186162912",
    (10, 157): "-0.062900125",
    (10, 158): "-0.082871522",
    (10, 159): "0.198702064",
    (10, 160): "0.212926331",
    (10, 161): "-0.190989852",
    (10, 162): "-0.124019328",
    (10, 163): "-0.230670258",
    (10, 164): "0.231024036",
    (10, 165): "0.098889515",
}

def _key_from_stdin():
    if os.environ.get("NUMBERDB_KEY_FROM_STDIN") != "1":
        return
    token = sys.stdin.read().strip()
    if "=" in token and token.split("=", 1)[0].isupper():
        token = token.split("=", 1)[1].strip().strip("'\"")
    if token:
        os.environ["NUMBERDB_API_KEY"] = token


def hyperbolic_names():
    keys = [key for key in names() if key in KNOTINFO_CS]
    if len(keys) != 243:
        raise ArithmeticError("expected 243 hyperbolic knots, got %d" % len(keys))
    if set(keys) != set(KNOTINFO_BRAIDS) - set(TORUS):
        raise ArithmeticError("the Chern-Simons source and braid table name different hyperbolic knots")
    return keys


def _exterior(n, k, mirror=False):
    """The knot exterior from KnotInfo's braid word, geometrically triangulated."""
    strands, word = KNOTINFO_BRAIDS[(n, k)]
    link = Link(braid_closure=[int(letter) for letter in word])
    if mirror:
        link = link.mirror()
    manifold = link.exterior()
    if manifold.num_cusps() != 1:
        raise ArithmeticError("%d_%d: the exterior has %d cusps"
                              % (n, k, manifold.num_cusps()))
    for attempt in range(50):
        if manifold.solution_type() == GEOMETRIC:
            break
        manifold.randomize()
    else:
        raise ArithmeticError("%d_%d: no geometric triangulation found" % (n, k))
    return manifold


def _certified(n, k, mirror=False):
    """(CS, 2 pi^2 CS) as certified intervals, with the 2-torsion class fixed.

    The interval from SnapPy pins 2 pi^2 CS modulo pi^2/2. The quad-double
    value says which multiple of pi^2/2 to add, and is then required to agree
    with the result -- so a wrong class, which is 0.25 away in CS, cannot pass.
    """
    manifold = _exterior(n, k, mirror)
    complex_volume = manifold.complex_volume(verified_modulo_2_torsion=True,
                                             bits_prec=BITS)
    unnormalised = complex_volume.imag()
    field = unnormalised.parent()
    step = field.pi() ** 2 / 2
    scale = 2 * field.pi() ** 2

    quad = manifold.high_precision()
    approximate, accuracy = quad.chern_simons(accuracy=True)
    if accuracy < CLASS_DIGITS:
        raise ArithmeticError(
            "%d_%d%s: SnapPy estimates only %d digits of Chern-Simons accuracy"
            % (n, k, " mirror" if mirror else "", accuracy))
    target = field(repr(approximate)) * scale
    steps = ZZ(((target - unnormalised) / step).center().round())
    unnormalised = unnormalised + steps * step

    #The class, checked rather than trusted: the two candidates are 0.25 apart
    #in CS, so agreeing to thirty digits settles it.
    tolerance = (field(10) ** (-CLASS_DIGITS)).center()
    if abs((unnormalised / scale).center() - field(repr(approximate)).center()) \
            > tolerance:
        raise ArithmeticError(
            "%d_%d%s: the certified value %s does not agree with SnapPy's %s"
            % (n, k, " mirror" if mirror else "",
               (unnormalised / scale).center(), approximate))

    #The manifold, checked: the real part is the volume, and the volume table
    #computed the same one from the same braid word.
    if abs(complex_volume.real().center() - field(repr(quad.volume())).center()) \
            > tolerance:
        raise ArithmeticError(
            "%d_%d%s: volume %s does not agree with SnapPy's %s"
            % (n, k, " mirror" if mirror else "",
               complex_volume.real().center(), quad.volume()))

    return unnormalised / scale, unnormalised


def _supports(value, digits):
    """Does this interval support that many significant digits?"""
    centre = abs(Decimal(str(value.center())))
    radius = Decimal(str(value.absolute_diameter())) / 2
    if not centre:
        return True
    return radius / centre <= Decimal(10) ** (-digits - 1)


def _text(value, digits):
    """The written digits of a certified interval."""
    ball = RealBallField(BITS)(value)
    return to_text(ball, digits=digits)


def _entry_comment(n, k, image):
    name = subscript(n, k)
    if (n, k) in NAMED:
        name += ", " + NAMED[(n, k)]
    if image == "mirror":
        name = "mirror image of " + name
    parts = [name]
    parts.append("alternating" if is_alternating(n, k) else "non-alternating")
    if (n, k) in AMPHICHIRAL:
        parts.append("amphichiral")
    return "; ".join(parts)


class ChernSimonsKnotComplements(numberdb.Generator):

    table = TABLE
    parameters = ("n", "k", "knot", "normalisation")
    type = "R"
    digits = DIGITS
    rigour = "assumed-bound"
    files = ("generate.py", "knotinfo_prime_knots.py")

    #: (n, k) -> {"K": (cs, unnormalised), "mirror": ...}, so that the two
    #: normalisations of one knot cost one computation rather than two, and the
    #: mirror congruence can be checked across them.
    _computed = {}

    def enumerate(self):
        for n, k in hyperbolic_names():
            for image in ("K", "mirror"):
                if image == "mirror" and (n, k) in AMPHICHIRAL:
                    continue
                for normalisation in ("cs", "two-pi-squared-cs"):
                    yield {"n": int(n), "k": int(k), "knot": image,
                           "normalisation": normalisation}

    def _knot(self, n, k):
        if (n, k) in self._computed:
            return self._computed[(n, k)]

        found = {"K": _certified(n, k)}

        #KnotInfo's nine decimals, which is what this table held before: an
        #independent source agreeing with the class as well as the digits.
        cs = found["K"][0]
        printed = Decimal(KNOTINFO_CS[(n, k)])
        if abs(Decimal(str(cs.center())) - printed) > Decimal("5e-9"):
            raise ArithmeticError("%d_%d: computed %s, KnotInfo prints %s"
                                  % (n, k, cs.center(), printed))

        if (n, k) in AMPHICHIRAL:
            #CS = -CS in R/(1/2)Z, so CS is 0 or 1/4 exactly. The interval
            #decides which, and the twenty knots of this table are all 0.
            if not (cs.contains_zero() and abs(cs.center()) < 0.1):
                raise ArithmeticError(
                    "%d_%d is amphichiral but its invariant is %s, not 0"
                    % (n, k, cs.center()))
            found["K"] = (None, None)          # exactly zero, written as "0"
        else:
            found["mirror"] = _certified(n, k, mirror=True)
            #CS(mirror) = -CS(K) in R/(1/2)Z, and both representatives lie in
            #(-1/4, 1/4), so the sum is exactly zero.
            total = found["mirror"][0] + found["K"][0]
            if not total.contains_zero():
                raise ArithmeticError(
                    "%d_%d: CS(mirror) + CS(K) = %s, which does not contain 0"
                    % (n, k, total))

        self._computed[(n, k)] = found
        return found

    def value(self, params, digits):
        n, k = int(params["n"]), int(params["k"])
        image, normalisation = params["knot"], params["normalisation"]
        if (n, k) not in KNOTINFO_CS:
            raise ValueError("%d_%d is not a hyperbolic prime knot with at most ten crossings" % (n, k))
        if image not in ("K", "mirror"):
            raise ValueError("knot is 'K' or 'mirror', not %r" % image)
        if image == "mirror" and (n, k) in AMPHICHIRAL:
            raise ValueError("%d_%d is amphichiral and has one entry" % (n, k))
        if normalisation not in ("cs", "two-pi-squared-cs"):
            raise ValueError("normalisation is not recognised: %r" % normalisation)

        found = self._knot(n, k)[image]
        comment = _entry_comment(n, k, image)
        if found[0] is None:
            return {"number": 0, "comment": comment}

        value = found[0] if normalisation == "cs" else found[1]
        if not _supports(value, digits):
            raise ArithmeticError(
                "%d_%d %s %s: the certified interval supports fewer than %d "
                "digits (radius %s)"
                % (n, k, image, normalisation, digits,
                   value.absolute_diameter() / 2))
        return {"number": _text(value, digits), "comment": comment,
                "digits": digits}


if __name__ == "__main__":
    _key_from_stdin()
    generator = ChernSimonsKnotComplements()
    publishing = (os.environ.get("NUMBERDB_PUBLISH") == "1"
                  or "--publish" in sys.argv)
    if publishing:
        print(generator.publish(
            overwrite=True,
            message="a hundred certified digits from SnapPy's verified complex "
                    "volume, replacing KnotInfo's nine; the 2-torsion class is "
                    "fixed by the quad-double value and by KnotInfo"))
    else:
        report = generator.verify(sample=None)
        print(report)
        sys.exit(0 if report.ok else 1)