(* ::Package:: *)

(*Authors: Erhard Aichinger, Sebastian Kreinecker *)


BeginPackage["Gen`"];


Nearringgenerator::usage="Nearringgenerator[basislist_, var_, rounds_] returns a list of the length of basislist which is sorted by degree and which contains polynomials of the subnearring of (Z[var], +, circ) which is generated by 
basislist. In each round one polynomial is randomly generated and reduced due to the current list";


Begin["`Private`"];


Clear[DeleteZeros];
Clear[Deg];
Clear[Coeff];
Clear[Lcf];
Clear[Comp];
Clear[PolyReduction];
Clear[PolyReductionEv];
Clear[PolyMult];
Clear[ReducedCoeff];
Clear[ExtBasisList];
Clear[ExtBasisAdd];
Clear[NewElem];
Clear[ReducedBasis];
Clear[RedBasis];
Clear[Nearringgenerator];
Clear[NearComp];
Clear[BasisUpdate];
Clear[GetNewPol];


$IterationLimit = Infinity;


DeleteZeros[{L1___, 0, L2___}]:= DeleteZeros[{L1, L2}];
DeleteZeros[{L___}]:={L};


(*Deg returns the degree of f*)
Deg[f_, var_] := Length[CoefficientList[f, var]] - 1;

(*Coeff returns the coefficient at x^i*)
Coeff[f_, i_, var_] /; f===0 := 0;
Coeff[f_, i_, var_] /; Deg[f,var]<i := 0;
Coeff[f_, i_, var_]:= CoefficientList[f, var][[i+1]];

(*Lcf returns the leading coefficient of f*)
Lcf[f_, var_] /; f===0 := 0;
Lcf[f_, var_] := Last[CoefficientList[f, var]];

(*Comp returns the composition f\[Degree]g*)
Comp[f_, g_, var_]/; g===0 || f===0:= Coeff[f, 0, var];
Comp[f_, g_, var_] := Expand[Sum[Coeff[f, i, var]*g^i,{i, 0, Deg[f,var]}]];


(*PolyReduction returns the polynomial which results from p by reducing it with respect to the polynomials of basislist*)
PolyReduction[p_, basislist_, var_] /; p===0 := 0;
PolyReduction[p_, basislist_, var_] := Expand[PolyReductionEv[p, basislist, var, Deg[p, var]]];
PolyReductionEv[p_, basislist_, var_, n_] /; n<0 := p;
PolyReductionEv[p_, basislist_, var_, n_] := PolyReductionEv[p - PolyMult[p, basislist, var, n]* basislist[[n+1]], basislist, var, n-1];
PolyMult[p_, basislist_, var_, n_]/; basislist[[n+1]]===0 := 0;
PolyMult[p_, basislist_, var_, n_] := ReducedCoeff[Coeff[p, n, var],Coeff[basislist[[n+1]], n, var]];
ReducedCoeff[coeffToBeRed_, redCoeff_] := (coeffToBeRed - Mod[coeffToBeRed, redCoeff]) / redCoeff;


(*ExtBasisList returns a list where p is added as reduced polynomial to the list at the right place*)
ExtBasisList[p_, basislist_, var_] := ExtBasisAdd[PolyReduction[p, basislist, var],   basislist[[Deg[PolyReduction[p, basislist, var], var] + 1]], basislist, var];
ExtBasisAdd[pReduced_, q_, basislist_, var_] /; pReduced===0 := basislist;
ExtBasisAdd[pReduced_, q_,  basislist_, var_] := ReplacePart[basislist, Deg[pReduced, var] + 1 -> NewElem[pReduced, q, ExtendedGCD[Lcf[pReduced, var], Lcf[q, var]][[2]] ]];
NewElem[p1_, p2_, cofactors_] := cofactors[[1]]*p1 + cofactors[[2]]*p2;


(*ReducedBasis returns a list where all elements are reduced (degree is the degree of the greatest polynomial which is not checked)*)
ReducedBasis[basislist_, degree_,  var_] := RedBasis[basislist, degree, Length[basislist] - 1, var];
RedBasis[basislist_, m_, n_, var_] /; m > n := basislist;
RedBasis[basislist_, m_, n_, var_] := RedBasis[ReplacePart[basislist, n+1 -> PolyReduction[basislist[[n + 1]], ReplacePart[basislist, n+1 -> 0],var]], m, n - 1 ,var];


Nearringgenerator[basislist_, var_, rounds_] := NearComp[DeleteZeros[Take[basislist, - (Length[basislist] - 2)]], basislist, var, rounds ];
NearComp[original_, newList_, var_, rounds_] /; Mod[rounds, 2000]===0  := If[rounds===0, InputForm[newList], Print[InputForm[newList]]; NearComp[original, newList, var, rounds - 1]];
NearComp[original_, newList_, var_, rounds_] := NearComp[original, BasisUpdate[newList, GetNewPol[original, newList,  var], var], var, rounds - 1];
BasisUpdate[basis_, newPol_, var_] := ReducedBasis[ExtBasisList[newPol, basis, var], Deg[newPol, var], var];

GetNewPol[original_, newList_, var_]:= Module[{maxDeg, newPol, left, right1, right, count, i, j, realList},

        left = RandomChoice[original];
        count = Deg[left, var];
        maxDeg = Floor[(Length[newList] - 1) / count ];

        realList = DeleteZeros[Take[newList, maxDeg + 1 ]];
        right1 = RandomChoice[realList];
    	right = right1 + Sum[ RandomChoice[ {RandomChoice[Table[Length[realList]/j, {j, 1, Length[realList]} ] -> realList], 0 }], {i, 1, count-1} ];
		
        newPol = Comp [left, right, var];

        Expand[newPol]
 ];



End[];


EndPackage[];