"""Simple Hurwitz numbers H_{g,mu} -- numberdb.org/T463.

Run with SageMath (or a compatible passagemath installation):
    sage -pip install numberdb
    sage -python generate.py                 # verify every stored entry
    sage -python checks.py                   # independent offline checks

Connected covers, labelled preimages of infinity, fixed distinct simple
branch points and inverse-automorphism weights. Only stable (g,len(mu))
are stored. The recurrence also uses the genuine unstable Hurwitz numbers.
It computes A(g,mu)=H(g,mu)/(2g-2+sum(mu)+len(mu))! exactly over QQ.
See Eynard, Mulase and Safnuk, Theorem 3.2, equation (3.5).
No floating point or approximate arithmetic is used.

The default enumeration includes every partition of degrees 1 through 10
and genera 0 through 5 satisfying 2*g-2+len(mu)>0. This is a complete
initial range by degree and genus; the maximum entry has 40 characters.

This file verifies by default. Filling the draft is performed separately
by fill.py. Neither this file nor that script publishes or offers a draft.
For a private draft, set NUMBERDB_API_KEY in the environment, or set
NUMBERDB_KEY_FROM_STDIN=1 and pipe the key on standard input.
"""

import numberdb.sage as numberdb
from sage.rings.rational_field import QQ
from functools import lru_cache
from math import factorial
import os
import sys

MAX_GENUS = 5
MAX_DEGREE = 10


def partitions(total, largest=None):
    """Non-increasing partitions, with larger parts first."""
    if total == 0:
        yield ()
        return
    largest = total if largest is None else min(largest, total)
    for first in range(largest, 0, -1):
        for rest in partitions(total-first, first):
            yield (first,) + rest


def canonical(parts):
    return tuple(sorted(parts, reverse=True))


@lru_cache(None)
def normalized(g, mu):
    """H/r! by the connected cut-and-join recurrence with labelled cycles."""
    if g < 0 or not mu:
        return QQ(0)
    r = 2*g - 2 + sum(mu) + len(mu)
    if r == 0:
        assert g == 0 and mu == (1,)
        return QQ(1)
    if r < 0:
        return QQ(0)
    result = QQ(0)
    n = len(mu)
    for i in range(n):
        for j in range(i+1, n):
            rest = tuple(mu[k] for k in range(n) if k not in (i,j))
            result += (mu[i]+mu[j])*normalized(g, canonical(rest+(mu[i]+mu[j],)))
    for i, m in enumerate(mu):
        rest = mu[:i]+mu[i+1:]
        for a in range(1, m):
            b = m-a
            contribution = normalized(g-1, canonical(rest+(a,b)))
            for mask in range(1 << len(rest)):
                left = tuple(rest[j] for j in range(len(rest)) if mask & (1 << j))
                right = tuple(rest[j] for j in range(len(rest)) if not mask & (1 << j))
                for ga in range(g+1):
                    contribution += normalized(ga, canonical(left+(a,))) * normalized(g-ga, canonical(right+(b,)))
            result += QQ(a*b)*contribution/QQ(2)
    return result/QQ(r)


def hurwitz(g, mu):
    """Exact labelled connected Hurwitz number, including unstable helpers."""
    if not isinstance(g, int) or g < 0 or not mu or any(int(x) != x or x <= 0 for x in mu):
        raise ValueError('A nonnegative genus and a nonempty positive profile are required')
    mu = canonical(tuple(int(x) for x in mu))
    r = 2*g-2+sum(mu)+len(mu)
    return normalized(g, mu)*factorial(r)


class SimpleHurwitz(numberdb.Generator):
    table = 'T463'
    parameters = ('g','mu')
    type = 'Q'
    rigour = 'exact'

    def enumerate(self, max_genus=MAX_GENUS, max_degree=MAX_DEGREE):
        for g in range(max_genus+1):
            for degree in range(1,max_degree+1):
                for mu in partitions(degree):
                    if 2*g-2+len(mu)>0:
                        yield {'g':g, 'mu':','.join(map(str,mu))}

    def value(self, params, digits=100):
        g = int(params['g'])
        mu = tuple(map(int, str(params['mu']).split(',')))
        assert mu == canonical(mu) and 2*g-2+len(mu)>0
        return hurwitz(g,mu)


if __name__ == '__main__':
    if os.environ.get('NUMBERDB_KEY_FROM_STDIN') == '1':
        os.environ['NUMBERDB_API_KEY'] = sys.stdin.read().strip()
    report = SimpleHurwitz().verify(sample=None)
    print(report)
    sys.exit(0 if report.ok else 1)