#
# Test the Kontsevich conjecture on the spanning tree generating
# function of a graph, and related conjectures.
# Should work with Maple V.2,3,4,5. Tested thoroughly with V.3.
# Requires the reduce() procedure (a separate file).
#
# There are two main data structures:
# A GRAPH is defined to be a set of pairs (2-sets) of positive integers.
#
#   Ex.: {{1,2},{2,3},{3,4},{1,4}}      (a 4-cycle).
#
# A BIGRAPH is defined to be a set of ordered pairs (lists) [-i,j],
# where i and j are positive integers. 
#
#   Ex.: {[-1,1],[-1,2],[-2,1],[-2,2]}  (also a 4-cycle).
#
# A bigraph is simply a particular way to describe a bipartite graph.
#
##########################################
#
# Main procedures defined in this file:
#
# Q         compute the spanning tree generating function of a graph
# AQ        compute the determinant of a generic matrix (symmetric with
#           indeterminates on the diagonal in the graph case, non-symmetric
#           in the bigraph case) supported on the edges of a (bi)graph
# kontpoly  use reduce() to compute the Kontsevich "polynomial" for the
#           graph or bigraph G; i.e., the number of *non-zero* evaluations
#           of Q(G) or AQ(G) over F_q, as a function of q. 
# test      apply kontpoly() to each member of a list of graphs or bigraphs,
#           returning the indices of those for which reduce() fails or the
#           result is not transparently polynomial.
#
# Read the documentation below to learn about various options that
# control the behavior of kontpoly() and test().
#
# Note: Z, q, and x are used as global variables.
#
# Additional procedures:
#
# K         produce a complete or complete bipartite graph (not bigraph)
# C         produce an n-cycle
# stree     produce a spanning tree/forest for a graph or bigraph
#
# See the files 'kgraphs', 'agraphs', 'bigraphs', and 'blocks' for
# lists of graphs and bigraphs that are useful for testing purposes.
#
# Read the documentation for reduce() to learn about the flags that can
# be set for verbosity and saving intermediate results to files.
#
##########################################
#
# Examples:
#
# G:=C(4); Q(G); kontpoly(G);
#
# read kgraphs; verbosity:=2;
# test(gr[6,9]);
#
# verbosity:=1; G:=K(3,2); kontpoly(G,apex);
# verbosity:=0; kontpoly(G,apex,tree);
# 
# read agraphs; verbosity:=2;
# test(gr[7,8],apex,tree);
#
# read bigraphs; G:=gr[5,15][1];
# verbosity:=1; test(gr[5,15],tree);
# verbosity:=2; kontpoly(G);
#
###########################################################################

if not assigned(reduce) then read reduce fi;

# The spanning tree generating function of graph G, using variables x[i,j].

Q:=proc(G) local A,i,e,n;
  if not type(G,set(set)) then ERROR(`use only with graphs`) fi;
  n:=max(op(map(op,G)));
  A:=array('symmetric','sparse',1..n,1..n);
  for i to n do A[i,i]:=0 od;
  for e in G do A[op(e)]:=-x[op(e)];
    A[e[1],e[1]]:=A[e[1],e[1]]+x[op(e)];
    A[e[2],e[2]]:=A[e[2],e[2]]+x[op(e)];
  od;
  linalg[det](linalg[submatrix](A,1..n-1,1..n-1));
end;

# The determinant of a generic matrix (symmetric in the graph case,
# non-symmetric in the bigraph case) supported on the edges of G.
#
# For a graph, the diagonal entries are x[i], and the off-diagonals are
#    x[i,j]=x[j,i] if {i,j} is an edge of G
#          0       if {i,j} is not an edge of G
#
# For a bigraph, the entries are
#    x[i,j]  if [-i,j] is an edge of G
#       0    if [-i,j] is not an edge of G

AQ:=proc(G) local A,e,n,i;
  n:=max(op(map(op,G))); A:=NULL;
  if type(G,set(set)) then A:='symmetric' fi; # we have a graph
  A:=array(A,'sparse',1..n,1..n);
  for e in G do i:=abs(e[1]),e[2]; A[i]:=x[i] od;
  if type(G,set(set)) then for i to n do A[i,i]:=x[i] od fi;
  linalg[det](A);
end;

# Find a spanning tree/forest for G. Works with graphs and bigraphs. 

stree:=proc(G) local C,T,e,X,GX,i,n;
  X:=map(op,G); n:=nops(X);
  C:=[$1..n]; T:=NULL;
  GX:=subs({seq(X[i]=i,i=1..n)},G);
  for e in GX while nops([T])<n-1 do
    if C[e[1]]<>C[e[2]] then C:=subs(C[e[2]]=C[e[1]],C); T:=T,e fi
  od;
  subs({seq(i=X[i],i=1..n)},{T});
end;

# Compute the Kontsevich "polynomial" for graph or bigraph G; i.e.,
# the number of NON-ZERO evaluations of Q(G) or AQ(G) over F_q.
#   kontpoly(G)      uses Q(G)  (for graphs) or AQ(G) (for bigraphs)
#   kontpoly(G,apex) uses AQ(G) (for graphs)
# If there is an additional argument 'tree', then a spanning tree/forest
# T is chosen for G, and 1's are assigned to the variables of T.
# Ex:   kontpoly(G,tree);  kontpoly(G,apex,tree);

kontpoly:=proc(G) local f,vars,i,T,e;
  if type(G,set(list)) or member('apex',[args])
    then f:=AQ(G) else f:=Q(G) fi;
  if member('tree',[args]) then T:=stree(G);
    f:=subs({seq(x[abs(e[1]),e[2]]=1,e=T)},f);
  fi;
  vars:=indets(f);
  f:=subs({seq(vars[i]=cat('x',i),i=1..nops(vars))},f);
  sort(expand(q^nops(vars)*subs(Z[]=1,1-reduce(Z[f]))))
end;

# Given a list of graphs lg, display the result of applying kontpoly()
# to each member, and return a list of indices for the graphs for which
# kontpoly() fails, or produces a result that is not transparently a
# polynomial in q. Additional arguments are passed to kontpoly().

test:=proc(lg) local g,i,sa,flags;
  if not (type(lg,list(set(set))) or type(lg,list(set(list)))) 
    then ERROR(`input must be a list of graphs or bigraphs`) fi;
  sa:=NULL; flags:=args[2..nargs];
  for i to nops(lg) do
    print(i,lg[i]); g:=traperror(kontpoly(lg[i],flags)); print(g);
    if g=lasterror or not type(g,polynom(integer,q))
      then sa:=sa,i; lprint(g) fi;
  od; [sa]
end;

# K(n)   returns the complete graph on n vertices.
# K(m,n) returns a complete bipartite graph on m+n vertices. (Not bigraph.)
# C(n)   returns a circuit on n vertices.

C:=proc(n) local i; {seq({i,i+1},i=1..n-1),{n,1}} end;
K:=proc(n) local i,j;
  if nargs=1 then {seq(seq({i,j},j=1..i-1),i=1..n)}
    else {seq(seq({i,j},i=1..n),j=n+1..n+args[2])} fi
end;