1:- module(gf2orth, []).    2:- use_module(library(apply)).    3:- use_module(library(apply_macros)).    4:- use_module(zdd('zdd-array')).    5:- use_module(util(math)).    6:- use_module(util(meta2)).    7:- use_module(pac(basic)).    8:- use_module(pac(meta)).    9:- use_module(util(misc)).   10:- use_module(zdd(zdd)).   11:- use_module(zdd('zdd-gf2')).   12:- use_module(pac(op)).   13
   14% ?- orth_matrix(2, M),
   15%	time((zdd gf2_poly_solve_list(M, Sols), card(Sols, C))).
   16%@ % 1,000 inferences, 0.000 CPU in 0.000 seconds (93% CPU, 8130081 Lips)
   17
   18
   19% ?- orth_matrix(6, M),
   20%	time((zdd gf2_poly_solve_list(M, Sols), card(Sols, C))).
   21
   22% ?- between(1, 7, I), orth_matrix(I, M),
   23%	time((zdd gf2_poly_solve_list(M, Sols), card(Sols, C))),
   24%	format("Num of sols for orth-matrix(~w) = ~w~n~n", [I, C]),
   25%	fail; true.
   26% ?- between(1, 7, I), orth_matrix(I, M),
   27%	time((zdd gf2_poly_solve_list(M, Sols), card(Sols, C))),
   28%	format("Num of sols for orth-matrix(~w) = ~w~n~n", [I, C]),
   29%	fail; true.
   30
   31
   32% Via discourse 1969 paper.
   33% o(2t+1)= q^t PI(i=0, t-1)(q^2t - q^2i)
   34% o(2t+1)= (q^2t - 1)*o(2t)
   35% o(2t) = q^iPi(i=1, t-1)(q^2t-q2i)
 card_orth(+Q, +N, -V) is det
Q: 2,3,5, ... prime numbers N: 2,3,4,5,... V is unified with the number of O(GF(Q), N).
   42% ?- time(forall(between(2, 12, N), (card_orth(2, N, V), writeln(orth(N)= V)))).
   43% ?- time(forall(between(2, 20, N), (card_orth(2, N, V), writeln(orth(N)= V)))).
   44
   45card_orth(Q, N, V):-
   46	(	0 =:=  N mod 2 -> card_even_orth(Q, N, V)
   47	;	card_odd_orth(Q, N, V)
   48	).
   49%
   50card_orth_help(Q, T, I, U, V):- V is U* (Q^(2*T)-Q^(2*I)).
   51%
   52card_even_orth(Q, E, V):- T is E//2,
   53	foldnum(card_orth_help(Q, T), 1-(T-1), 1, U),
   54	V is Q^T* U.
   55%
   56card_odd_orth(Q, O, V):- E is O-1,
   57	card_even_orth(Q, E, U),
   58	V is (Q^E-1)*U.
   59
   60		/***************
   61		*     Notes    *
   62		***************/
   63
   64% For a, b in GF2, 1 and 2 hold.
   65% 1.  a  = 0  iff  1 + a = 1.
   66% 2.  ab = 1  iff  a = 1 and b = 1.
   67% 3.  Orthogonality of matrix A = [ aij ] ( n=< i, j =< n )
   68% is represented as a system of equations Eij = 0,
   69% where Eij is a polynomial with indeterminates aij obtained
   70% from the familiar multiplication of matricies.
   71
   72% It follows from 1, 2, 3 that the set Orth(n, GF2)
   73% and the intersection of sol(Eij=0) have the same cardinality.
   74% where Sol(E) is the set of solutions of E = 0.
   75
   76
   77		/***************************
   78		*     orth(GF2, 2) case    *
   79		***************************/
 gf2_count_orth(+N, -C) is det
N: N>0 C is unified with the number of orthogonal macricies over GF(2)
   85% ?- zdd gf2_count_orth(1, C).
   86% ?- zdd gf2_count_orth(2, C).
   87% ?- zdd gf2_count_orth(3, C).
   88% ?- zdd gf2_count_orth(4, C).
   89% ?- time((zdd gf2_count_orth(5, C))).
   90% ?- time((zdd gf2_count_orth(6, C))).
   91%@ % 25,200,092 inferences, 2.878 CPU in 3.789 seconds (76% CPU, 8756061 Lips)
   92%@ C = 23040.
   93
   94% ?- time((zdd gf2_count_orth(7, C))).
   95%@ % 844,702,506 inferences, 101.117 CPU in 103.467 seconds (98% CPU, 8353719 Lips)
   96%@ C = 1451520.
   97
   98% ?- time(gf2_count_orth(8, C)).
   99%@ % 77,269,642,066 inferences, 8239.463 CPU in 8300.403 seconds (99% CPU, 9377995 Lips)
  100%@ C = 185794560.
  101% 80GB
  102
  103gf2_count_orth(N, C, S):- orth_matrix(N, Ps),
  104					   gf2_poly_solve_list(Ps, X, S),
  105					   card(X, C, S).
  106
  107% ?- orth_matrix(2, X).
  108orth_matrix(N, ListOfEq0):- findall( CJK,
  109	   (	between(1, N, J),
  110			between(1, N, K),
  111			K =< J,
  112			findall(a(J, I) * a(K, I),
  113				   between(1, N, I),
  114				   CJK0),
  115			(	J=K -> CJK = +(CJK0)
  116			;	CJK = 1 + +(CJK0)
  117			)
  118	   ),
  119	   ListOfEq0).
  120
  121% ?- orth_matrix_cond(2, X), maplist(writeln, X).
  122orth_matrix_cond(N, ListOfEq0):- findall( CJK,
  123	   (	between(1, N, J),
  124			between(1, N, K),
  125			K =< J,
  126			findall(a(J, I) * a(K, I),
  127				   between(1, N, I),
  128				   CJK0),
  129			(	J=K -> CJK = 1 + +(CJK0)
  130			;	CJK = +(CJK0)
  131			)
  132	   ),
  133	   ListOfEq0).
  134%
  135orth_matrix_cond_one_poly(N, 1 + *(ListOfEq0)):- findall( CJK,
  136	   (	between(1, N, J),
  137			between(1, N, K),
  138			K =< J,
  139			findall(a(J, I) * a(K, I),
  140				   between(1, N, I),
  141				   CJK0),
  142			(	J=K -> CJK = +(CJK0)
  143			;	CJK = 1+(CJK0)
  144			)
  145	   ),
  146	   ListOfEq0)