#Author := Ayoola Jinadu
# If you find a bug in the code please email me at ajinadu@sfu.ca

x1 := 0: y1 := 0: z1 := 0: y2 := 0: z2 := 0: z3 := 0:
Sys := [ (x-x1)^2+ (y-y1)^2 +(z-z1)^2-t1^2,
         (x-x2)^2+ (y-y2)^2 +(z-z2)^2-t2^2,
         (x-x3)^2+ (y-y3)^2 +(z-z3)^2-t3^2,
         (x-x4)^2+ (y-y4)^2 +(z-z4)^2-t4^2,
         (x1-x2)^2+ (y1-y2)^2+(z1-z2)^2-d1^2,
         (x1-x3)^2+ (y1-y3)^2+(z1-z3)^2-d2^2,
         (x1-x4)^2+ (y1-y4)^2+(z1-z4)^2-d3^2,
         (x2-x3)^2+ (y2-y3)^2+(z2-z3)^2-d4^2,
         (x2-x4)^2+ (y2-y4)^2+(z2-z4)^2-d5^2,
         (x3-x4)^2+ (y3-y4)^2+(z3-z4)^2-d6^2]:        
 #Vaariables := [x, x2, x3, x4, y, y3, y4, z, z4];
 #Parameters :=  [d1, d2, d3, d4, d5, d6, t1, t2, t3, t4];
 #Eliminate :=  [x, x2, x3, x4, y, y3, y4, z, z4]:
 #Main variable := [t1]:

 elim := [x, x2, x3, x4, y, y3, y4, z, z4]:


kernelopts(gcmaxthreads=1);
kernelopts(threadlocalheapmultiplier=16);

read det; read newdeg; read dixres; read Bmbot;
read dixon;
read minor;
with(LinearAlgebra):

  M := dixonmatrix( Sys, elim ):
  rank,rows,cols := minor(M):
  printf(" rank = %d\n", rank ):
  B := M[rows,cols]: # extract minor
  L := StronglyConnectedBlocks(B):
  BlockStructure = map(RowDimension,L);

    X := [t1, d2, d3, d4, d5, d6, d1, t2, t3, t4]:
    R := convert(L[1],listlist):
    rt := [R,X,RowDimension(R)]:
    n :=rt[3];
    deg_M := [seq( max(seq(seq(degree(rt[1][i][j],var),i=1..n),j=1..n)), var in X )];
    GlobalCArray := Array(0..n-1,0..n-1,order=C_order,datatype=integer[8]):
    E :=  DixonRes(BB,X);
    nops( expand(E) );