//G28.sing
//
//  Singular script to compute numbers of real solutions to a problem in Schubert calculus,
//
//   Compute the number of 2-planes in C^8 that meet 4 four planes, nontrivially.
//
//////////////////////////////////////////////////////////////////////////////////////////
option(redSB);     // computes reduces Groebner bases
LIB "solve.lib";   // One way to get at the solutions
LIB "matrix.lib"; 
LIB "linalg.lib";
/////////////////////////////////////////////////////////////////////////////////
//
//   This is one way to compute an eliminant.
//
proc univarpol (ideal G, int i)
{
  ideal B = kbase(G);
  return(charpoly(coeffs(reduce(var(i)*B,G),B), varstr(i)));
}
/////////////////////////////////////////////////////////////////////////////////
// random number
int rand = 100;

//  Set up a ring
ring R = 0, (x(1..6)(1..2)), dp;

matrix H[8][2];   
matrix K[8][4];

ideal I;

H = x(1)(1), x(1)(2),
    x(2)(1), x(2)(2),
    x(3)(1), x(3)(2),
    x(4)(1), x(4)(2),
    x(5)(1), x(5)(2),
    x(6)(1), x(6)(2),
       1   ,    0   ,
       0   ,    1  ;
       

I = 0;
for (int i=1; i<=4; i++) {
 K = random(rand,8,4);
 I = I + minor( concat(H,K), 6 );
}
I = std(I);

//  This has dimension 0 and consists of four points
dim(I);
mult(I);

//  Eliminates to the second variable
//
//ideal J = univarpol(I,2);
//  Need to map this to a univariate ring, and then could apply symbolic methods, or just
//    do what comes next
//

// Let us solve it
def AC=solve(I,6,0,"nodisplay");
setring AC;

//  Number of solutions
size(SOL);
//
//  Print out some coordinates for each solution.
//
for (i=1; i<=size(SOL); i++){
 print([SOL[i][1],SOL[i][2],SOL[i][3],SOL[i][4]]);
 }
quit;