#43.maple
#
# Frank Sottile
# 29 October 1998
# Berkeley, CA
#
#   We work to show instances of Shapiro's total positivity conjecture
# for the case of (135)^4 = 8 in G(3,4).
#
# Number   Seed                    Time       Size   
#  1000     1     Make systems    1110      4828684 \
#                 Eliminants       436       986913  > Schur
#                 Real roots       131              /
#
#  3000   123     Make systems     3319     14471954 \
#                 Eliminants       1308      2956498  > Schur
#                 Real roots       1905              /
#
interface(quiet=true):

with(linalg):
_seed:=123:
die:=rand(1..10):

ELM:=proc(X,Y,Z,W,S,T)
matrix([[1,X,0,0,0,0,0],[0,1,Y,0,0,0,0],[0,0,1,Z,0,0,0],
[0,0,0,1,W,0,0],[0,0,0,0,1,S,0],[0,0,0,0,0,1,T],[0,0,0,0,0,0,1]])
end:

Pos:=proc(a1,b1,b2,c1,c2,c3,d1,d2,d3,d4,e1,e2,e3,e4,e5,f1,f2,f3,f4,f5,f6)
multiply(ELM(a1,0,0,0,0,0),ELM(b1,b2,0,0,0,0),ELM(c1,c2,c3,0,0,0),
	ELM(d1,d2,d3,d4,0,0),ELM(e1,e2,e3,e4,e5,0),ELM(f1,f2,f3,f4,f5,f6))
end:

Coords:=matrix([
	[1,a,b,0,0,0,0],
	[0,0,1,c,d,0,0],
	[0,0,0,0,1,e,f]]):

COLS:=[[1,2,3,4,5,6],[1,2,3,4,5,7],[1,2,3,4,6,7],
[1,2,3,5,6,7],[1,2,4,5,6,7],[1,3,4,5,6,7],[2,3,4,5,6,7]]:

ROWS:=[[1,2,3,4,5,6,7],[1,2,3,4,5,6,8],[1,2,3,4,5,7,8],[1,2,3,4,6,7,8],
[1,2,3,5,6,7,8],[1,2,4,5,6,7,8],[1,3,4,5,6,7,8],[2,3,4,5,6,7,8]]:

M.0:=band([1],7):

ST_Nreal:="lprint(`realroots = `,NREAL,`seed = `,":
ST_end:=");":

#############################################################################
#
#
#lprint(`timer = 1;`);
lprint(`option(redSB);`);
lprint(`int t=timer;`);
lprint(`ring R= 0, (a,b,c,d,e,f), (dp(5),dp(1));`);
lprint(`ideal I;`);
lprint(`ideal G;`);

lprint(`print("interface(quiet=true):");`);
lprint(`print("readlib(realroot):");`);


for ntrials from 1 to 3000 do

SEED:=_seed:

 
 for j from 1 to 2 do
  Params:=[]:
  for i from 1 to 21 do Params:=[Params[],die()]: od:
  M.j:=multiply(Pos(Params[]),M.(j-1));
 od:

 lprint(`I=`);

 for j from 1 to 2 do
   BM:=stackmatrix(submatrix(M.j,[1,2,3,4,5],[1,2,3,4,5,6,7]),Coords):
   SM:=stackmatrix(submatrix(M.j,[1,2,3],[1,2,3,4,5,6,7]),Coords):

 for Row in ROWS do
  lprint(det(submatrix(BM,Row,[1,2,3,4,5,6,7])),`,`);:
 od:

 for Col in COLS do
  lprint(det(submatrix(SM,[1,2,3,4,5,6],Col)),`,`);:
 od: od:
 lprint(`0;`);

 lprint(`G = std(I);`);
 lprint(`print("POLY:=");`);
 lprint(`short=0;`);
 lprint(`G[1];`);
 lprint(`print(":");`);
 lprint(`print("NREAL:= nops(realroot(POLY, 1/100)): ");`);
 lprint(`print(" if NREAL<6 then");`);
 lprint(`print("`,convert(ST_Nreal,symbol));
 lprint(convert(SEED,symbol));
 lprint(convert(ST_end,symbol),` ");`);
 lprint(`print("fi: ");`);
 lprint(`print(" ");`);
 lprint(`print(" ");`);
 lprint(`print("#######################################################");`);
 lprint(`print("#                                                     #");`);
 lprint(`print("#          Next One                                   #");`);
 lprint(`print("#                                                     #");`);
 lprint(`print("#######################################################");`);
 lprint(`print(" ");`);
od:

lprint(`print("time(); ");`);
lprint(`print("quit; ");`);
lprint(`timer-t;`);
lprint(`quit;`);

time();

quit;












BigMat := s -> matrix([
		[ 1 , s ,s^2, s^3 , s^4 ,  s^5 ,   s^6],
		[ 0 , 1 ,2*s,3*s^2,4*s^3, 5*s^4, 6*s^5],
		[ 0 , 0 , 1 , 3*s ,6*s^2,10*s^3,15*s^4],
		[ 0 , 0 , 0 ,  1  , 4*s ,10*s^2,20*s^3],
		[ 0 , 0 , 0 ,  0  ,  1  , 5*s  ,15*s^2],
		[ 1 , a , b ,  0  ,  0  ,  0   ,  0   ],
		[ 0 , 0 , 1 ,  c  ,  d  ,  0   ,  0   ],
		[ 0 , 0 , 0 ,  0  ,  1  ,  e   ,   f  ]]):

#		[ 1 ,x12,x13,  0  ,  0  ,  0   ,  0   ],
#		[ 0 , 0 , 1 , x24 , x25 ,  0   ,  0   ],
#		[ 0 , 0 , 0 ,  0  ,  1  , x36  ,  x37 ]]):

SmMat  := s -> submatrix(BigMat(s),[1,2,3,6,7,8],[1,2,3,4,5,6,7]):

DD:=det(submatrix(SmMat(s),[1,2,3,4,5,6],[1,2,3,4,5,6])):
EE:=det(submatrix(SmMat(s),[1,2,3,4,5,6],[1,2,3,4,5,7])):
FF:=det(submatrix(SmMat(s),[1,2,3,4,5,6],[1,2,3,4,6,7])):
GG:=det(submatrix(SmMat(s),[1,2,3,4,5,6],[1,2,3,5,6,7])):

AA:=det(submatrix(BigMat(s),[1,2,3,4,5,6,7],[1,2,3,4,5,6,7])):
BB:=det(submatrix(BigMat(s),[1,2,3,4,5,6,8],[1,2,3,4,5,6,7])):
CC:=det(submatrix(BigMat(s),[1,2,3,4,5,7,8],[1,2,3,4,5,6,7])):

K:=gbasis([DD,EE,FF,GG,AA,BB,CC],wdeg([1,2,1,2,1,2],[a,b,c,d,e,f])):


lprint(`//Execute this with:`);
lprint(`//Singular --ticks-per-sec=100<37.sing>37.eqs`);
lprint(`//`);
lprint(`timer=1;`);
lprint(`int T = timer;`);
lprint(`short=0;`);
lprint(`ring R= 0, (s,t,a,b,c,d,e,f), wp(100,100,10,10,20,20,30,30);`);
lprint(`ideal I=`);
for i from 1 to 6 do   lprint(K[i],`,`);  od:
for i from 1 to 6 do   lprint(subs(s=t,K[i]),`,`);  od:
lprint(`0;`);

lprint(`ideal G=std(I);`);
lprint(`timer-T;`);
lprint(`G;`);
lprint(`quit;`);

quit;