//
// Table 2 of the article
//

Attach("absfact.m");
import "absfact.m": absfact_random_example, multi_factor_hensel_lifting, gao_matrix;

SetVerbose("AbsoluteFactorization",4);
//SetVerbose("PolyFact", 1);
// partial implementations
gao_gaussian:= function(F)
/*
  Returns the time to solve the linear system.
*/
  M := gao_matrix(F);
  ns_time:= Cputime();
  E:= Nullspace(M);
  return Cputime()-ns_time;
end function;

gao_blackbox:=function(F)
/*
  Simulates the computation of the matrix-vector products that would
  occur in the resolution "à la Wiedemann" of Gao's system.
  Returns the corresponding time.
*/

  d:=TotalDegree(F);
  tmp:= [ UnivariatePolynomial(z) : z in Coefficients(F,2)];
  UPx<x>:= Universe(tmp);
  UPy<y>:= PolynomialRing(UPx);
  F:= UPy!tmp;
  delete tmp;
  K:= BaseRing(UPx);

  random_polynomial:= function()
    return UPy![ UPx![ K!Random(10^9) : j in [0..d-i]] : i in [0..d-1]];
  end function;

  G:=random_polynomial();

  begin_time:= Cputime();
  for i:=1 to 3*d do
    _:= Derivative(F)*G;
  end for;
  return (d+1)*(Cputime() - begin_time);
end function;

trager:=function(F)
/*
  Does Trager's stragegy for one factor and lift to precision d/2.
*/
  MP:=Parent(F);
  K:=BaseRing(MP);
  d:=TotalDegree(F);

  UP<t>:= PolynomialRing(K);
  f:=Evaluate(F,[0,t]);

  "Rational factorization";
  time tmp:=Factorization(f);

  e:=tmp[1][1];
  "Degree of the extension";
  Degree(e);
  E:=quo<UP|e>;

  UP_E:= PolynomialRing(E);
  "Rational factorization over E";
  time _:=Factorization(F);
  lf:=[ z[1] : z in Factorization(UP_E!f)];

  tmp:= [ UP_E!UnivariatePolynomial(z) : z in Coefficients(F,2)];
  UPx<x>:= Universe(tmp);
  UPy<y>:= PolynomialRing(UPx);
  F1:= UPy!tmp;
  delete tmp;

  "Multi factor Hensel lifting";
  time _:=multi_factor_hensel_lifting(lf,F1,2*d+1);

//  "Native function";
//  MP_E:= PolynomialRing(E,2);
//  time _:=Factorization(MP_E!F);
  return true;
end function;

K:=GF(754974721);
MP<x,y>:= PolynomialRing(K,2);

range:= [3..7];
TODO:= [ [2^(i - (i div 2)) , 2^(i div 2)] : i in range];

/*
for z in TODO do
  "Building test ", z;
  F:= absfact_random_example(MP,z[1],z[2],10^9);
  "Gao naif...";
  ResetMaximumMemoryUsage();
  "time for linear resolution: ", gao_gaussian(F);
  "MaximumMemoryUsage (MB): ", GetMaximumMemoryUsage()/ 1024.0^2 ;
end for;
*/

/*
for z in TODO do
  "Building test ", z;
  F:= absfact_random_example(MP,z[1],z[2],10^9);

  "Gao blackbox...";
  ResetMaximumMemoryUsage();
  "time: ", gao_blackbox(F);
  "MaximumMemoryUsage (MB): ", GetMaximumMemoryUsage()/ 1024.0^2 ;
end for;
*/

for z in TODO do
  "Building test ", z;
  F:= absfact_random_example(MP,z[1],z[2],10^9);

  "Tager...";
  ResetMaximumMemoryUsage();
  time _:= trager(F);
  "MaximumMemoryUsage (MB): ", GetMaximumMemoryUsage()/ 1024.0^2 ;
end for;

quit;
