pjt33

def extend_closure(omega, included, extension, excluded):
	# We start with a closure which we want to extend
	closure = set(included)

	# Use a queue
	q = set([extension])
	while q:
		x = next(iter(q))
		q.remove(x)
		closure.add(x)

pjt33

from collections import defaultdict
from itertools import combinations
from sage.rings.polynomial.real_roots import real_roots

R.<alpha, x, y> = PolynomialRing(QQ, order='lex')
S.<z> = PolynomialRing(QQ)

def P(n, p, q):
    return ((1 - x^p) + alpha * (x^n - x^q)) // (1 - x^gcd(p, abs(n-q)))

pjt33

#!/usr/bin/pypy3

from heapq import heappush, heappop

def mod_inv(a: int, m: int) -> int:
	lastremainder, remainder = a % m, m
	x, lastx, y, lasty = 0, 1, 1, 0
	while remainder:
		lastremainder, (quotient, remainder) = remainder, divmod(lastremainder, remainder)
		x, lastx = lastx - quotient*x, x

pjt33

1, 2, 3, 4, 5, 6
1, 2, 3, 4, 6, 5
1, 2, 3, 5, 4, 6
1, 2, 3, 5, 6, 4
1, 2, 3, 6, 4, 5
1, 2, 3, 6, 5, 4
1, 2, 4, 3, 5, 6
1, 2, 4, 5, 3, 6
1, 2, 5, 3, 4, 6
1, 2, 5, 4, 3, 6

pjt33

using System;
using System.Collections.Generic;
using System.Linq;
using BITMASK = System.UInt64;

namespace Sandbox
{
    using STATE = System.ValueTuple<BITMASK, System.Byte, System.Byte>;

    // One way to formalise the snake constraint is to treat internal vertices as deleted.

pjt33

def swap_sequences_inlined(n, swap_vars):
  f = factorial(n)

  def extend_swap_sequence(seq, equivalence_classes, distrib, swap_vars):
    if len(distrib) * (2^len(swap_vars)) < f:
      # We can't hit all permutations, so fail fast
      return

    if not swap_vars:
      yield (seq, [weight - 1/f for weight in distrib.values()])

pjt33

def integer_relation(vals_to_relate, precision = 200):
    n = len(vals_to_relate)
    m = matrix(n, n+1)
    for i in range(n): m[i,i] = 1
    for i in range(n):
        m[i,n] = int(vals_to_relate[i]*RealField(precision)(2)^(precision-20))//2^10

    L = m.LLL()
    return L.row(0)

pjt33

In principle, we should derive essentially the same solutions for $(k,l,m)$ as for $(l,k,m)$ by symmetry in the upper arguments,
and for $(k,l,m)$ as for $(-k,-l,-m)$ by effectively negating the sign of $t$.
I haven't implemented negative $l$ or positive $m$ so can't check the latter.
For the former, I see the symmetry correctly observed in many cases but there is the odd troubling exception.

(1,0,-1)
 Ideal (x - 2, b0, a0 - b0 + c0 - 1, A - 1)
(0,1,-1)
 Ideal (x, a0)
 Ideal (a0, A - 1)

pjt33

from collections import Counter
from fractions import Fraction


def _gcd(a, b):
	while a:
		a, b = b % a, a
	return b

pjt33

#!/usr/bin/python3

# https://math.stackexchange.com/q/3444274

from fractions import Fraction


def build_states(target_dragon):
    states = {}
    def toBinary(x, width):