reduce-cycle.py

#!/usr/bin/python

"""
Application of Theorem 6.12 to triangles, 4-cycles (Theorem 6.13)
and 5.cycles (Theorem 6.14).
"""

import sys, os
sys.path.append(os.path.dirname(__file__) + "/..")

from graph_tools.base import *
from graph_tools.parameters import CircuitDoubleCover, VertexCount
from graph_tools.misc import sun

from sage.all import MixedIntegerLinearProgram
from itertools import permutations

verbose = False

def print_v(fmt, *args):
  if verbose:
    print(fmt % args)

def solve_gadget(g, alt_gadgets, c2, param=CircuitDoubleCover):
  N = g.size()
  print_v("N: %i", N)

  g_n = VertexCount.finalize(g.eval_gadget(VertexCount))

  P = MixedIntegerLinearProgram(maximization=False)
  V = P.new_variable(real=True, nonnegative=True)
  P_dual = MixedIntegerLinearProgram(maximization=True)
  V_dual = P_dual.new_variable(real=True, nonnegative=True)

  BOUNDARIES = list(param.enumerate_boundaries(N))
  print_v("boundaries: %s", BOUNDARIES)

  def fake_gadget(b):
    return FakeGadget(N, [ BoundaryValue(b, 1) ])

  def join(a, b):
    return Gadget.join([ a, b ], [ ((1, i), (2, i)) for i in range(1,N+1) ], [])

  def gadget_to_multiplicity_vector(gadget):
    return [
      join(fake_gadget(b), gadget).eval(param) for b in BOUNDARIES
    ]

  dual_objective = 0
  dual_constraints = [0] * len(BOUNDARIES)
  def gadget_to_constraint(gadget, var):
    nonlocal dual_objective

    mv = gadget_to_multiplicity_vector(gadget)
    n = VertexCount.finalize(gadget.eval_gadget(VertexCount))

    const = sum( a * V[b] for a, b in zip(mv, BOUNDARIES) ) >= c2**((n - g_n) / 2)
    P.add_constraint(const)
    print_v("constraint: %s", const)
    for i, a in enumerate(mv):
      dual_constraints[i] += a * V_dual[var]
    dual_objective += c2**((n - g_n) / 2) * V_dual[var]

  for i, ag in enumerate(alt_gadgets): gadget_to_constraint(ag, i)
  print_v("dual objective: %s", dual_objective)
  P_dual.set_objective(dual_objective)

  mv = gadget_to_multiplicity_vector(g)
  objective = sum( a * V[b] for a, b in zip(mv, BOUNDARIES) )
  print_v("objective: %s", objective)
  P.set_objective(objective)

  for b, const in zip(mv, dual_constraints):
    print_v("dual constraint: %s", const <= b)
    P_dual.add_constraint(const <= b)

  print_v("%s", P)
  print_v("%s", P_dual)

  ret = ( P.solve(), P_dual.solve() )

  #P_dual.show()
  for i, v in sorted(P_dual.get_values(V_dual).items()):
    print_v('w_%s = %s', i, v)

  return ret


def friend_boundary(g, param):
  N = g.size();
  BOUNDARIES = list(param.enumerate_boundaries(N))

  def test_boundary(b):
    fg = FakeGadget(N, [ BoundaryValue(b, 1) ])
    joins = [ ((1, i), (2, i)) for i in range(1,N+1) ]
    return (b, Gadget.join([ fg, g ], joins, []).eval(param))

  return [ test_boundary(b) for b in BOUNDARIES ]


def rot(l, x):
  return l[x:] + l[:x]


G4 = [
  Gadget.join([CUBIC_VERTEX] * 2, [((1,1), (2,1))], [ (1,2), (1,3), (2,2), (2,3) ]),
  Gadget.join([CUBIC_VERTEX] * 2, [((1,1), (2,1))], [ (1,2), (2,2), (1,3), (2,3) ]),
  Gadget.join([CUBIC_VERTEX] * 2, [((1,1), (2,1))], [ (1,2), (2,2), (2,3), (1,3) ]),
]

G4f = [
  Gadget.join([FREE_EDGE]*2, [], [ (1,1), (1,2), (2,1), (2,2) ]),
  Gadget.join([FREE_EDGE]*2, [], [ (1,1), (2,1), (2,2), (1,2) ]),
  Gadget.join([FREE_EDGE]*2, [], [ (1,1), (2,1), (1,2), (2,2) ]),
]

G5fp = [
  Gadget.join([ CUBIC_VERTEX, FREE_EDGE ], [], rot([ (1,1), (1,2), (1,3), (2,1), (2,2) ], i))
  for i in range(5)
]

G5t = [
  Gadget.join([CUBIC_VERTEX] * 3, [((1,1), (2,1)), ((2,2), (3,2))], list(out))
  for out in permutations([ (1,2), (1,3), (2,3), (3,1), (3,3) ])
]

G5f = [
  Gadget.join([ CUBIC_VERTEX, FREE_EDGE ], [], list(out))
  for out in permutations([ (1,1), (1,2), (1,3), (2,1), (2,2) ])
]

G5f_smart = G5fp + [
  Gadget.join([ CUBIC_VERTEX, FREE_EDGE ], [], rot([ (1,1), (1,2), (2,1), (1,3), (2,2) ], i))
  for i in range(5)
]

G5s = [
  Gadget.join([CUBIC_VERTEX] * 3, [((1,1), (2,1)), ((2,2), (3,2))], rot([ (1,2), (1,3), (2,3), (3,1), (3,3) ], i))
  for i in range(5)
]


if __name__ == "__main__":
  import sys
  verbose = "-v" in sys.argv

  print("3-cycle replaced by vertex: %s" % (
    solve_gadget(sun(3), [ CUBIC_VERTEX ], 2),))

  print("4-cycle replaced by all 3 matchings (c = sqrt(%s)): %s" % (
    2, solve_gadget(sun(4), G4f, 2),))
  print("4-cycle replaced by noncrossing matchings (c = sqrt(%s)): %s" % (
    2, solve_gadget(sun(4), G4f[:2], 2),))

  print("5-cycle replaced by a non-crossing cubic vertex and an edge (c = sqrt(%s)): %s" % (
    (5/2)**0.5, solve_gadget(sun(5), G5fp, (5/2)**0.5),))
  print("5-cycle replaced by a tree (c = sqrt(%s)): %s" % (
    3.75**0.5, solve_gadget(sun(5), G5t, 3.75**0.5),))
  print("5-cycle replaced by a cubic vertex and an edge (c = sqrt(%s)): %s" % (
    3.75**0.5, solve_gadget(sun(5), G5f, 3.75**0.5),))
  print("5-cycle replaced by a cubic vertex and an edge (smart) (c = sqrt(%s)): %s" % (
    3.75**0.5, solve_gadget(sun(5), G5f_smart, 3.75**0.5),))