#Read("C:\Andrea\\AlgebraI2006.gap");

######################################################
##Support File for the Algebra I course year 2005-2006
######################################################

###########################################
##Complexity of the Euclidean Algorithm
##We count the number of iterations in the euclidean algorithm
###########################################

NumIterEuclidean:=function(a,b)
    local cou,aa,bb,rr;
    cou:=0;
    aa:=a;
    bb:=b;
    repeat
        rr:=aa mod bb;
        aa:=bb;
        bb:=rr;
        cou:=cou+1;
    until rr=0;
    return cou;
end;

##and study the distribution of such values in the range [0..b-1] for a
##given integer b

NumOccurenceList:=function(list)
    local set;
    set:=Set(list);
    return List(set,y->[y,Size(Filtered(list,x->x=y))]);
end;

DistributionEuclidean:=function(int)
    local list, nli;
    list:=List([1..int],a->NumIterEuclidean(a,int));
    nli:=NumOccurenceList(list);
    return  nli;  #List(nli,x->x[2]);
end;


##############################################
##Subgroups of symmetric groups
##############################################

G:=SymmetricGroup(3);
subg:=SubgroupsSolvableGroup(G);
Size(subg);
ConjugacyClassesSubgroups(G);
H:=Group((1,2));
IsSubset(G,H);
List(G,Order);
setg:=Set(G);
List(setg,Order);
List(setg,x->[Order(x),x]);
rcgh:=RightCosets(G,H);
Size(rcgh);
Set(rcgh[2]);
List(rcgh[2],x->List(rcgh[2],y->x*y));
List(rcgh,Set);

################################################
##Factorization Polynomials
################################################

ZZ:=Integers;
x:=Indeterminate(Rationals,"x");

probfact:=function(numatt,deg)
    local t,g,f,counter;
    counter:=0;
    for t in [1..numatt] do
        f:=Sum(List([0..deg-1],i->Random(-1323515432,14326152152)*x^i))+x^deg;
        g:=Factors(f);
        if Length(g)=1 then counter:=counter+1;
        fi;
    od;
    return counter;
end;

f:=x^4+1;
Factors(f);

pri:=Filtered([2..101],IsPrime);

for p in pri do
    x:=Indeterminate(GF(p),"x");
    f:=x^4+1;
    fac:=Factors(f);
    Print("Number of Factors for x^4+1 over GF(",p,") equals ",Length(fac),"\n");
od;

lism1:=[];
lis2:=[];
lism2:=[];

for p in pri do
    K:=GF(p);
    squ:=List(K,x->x^2);
    if -Z(p)^0 in squ then Append(lism1,[p]);fi;
    if 2*Z(p)^0 in squ then Append(lis2,[p]);fi;
    if -2*Z(p)^0 in squ then Append(lism2,[p]);fi;    
od;

pri=Union(lism2,lism1,lis2);

List(lism1,x->x mod 4);
List(lis2,x->x mod 4);
List(lism2,x->x mod 4);

List(lis2,x->x mod 8);
List(lism2,x->x mod 8);

n:=Random(1,1324265);
repeat n :=n+2;until IsPrime(n);
K:=GF(n);
squ:=List(K,x->x^2);

esp:=QuoInt(n-1,2);
(2*Z(n)^0)^esp;
2^esp*Z(n)^0;

###############
##Finite Fields
###############

p:=5;
x:=Indeterminate(GF(p),"x");
deg:=2;  
f:=Sum(List([0..deg-1],i->Random(GF(p))*x^i))+x^deg;                                             
g:=Factors(f);
h:=g[2];
IsIrreducible(h);
dh:=Derivative(h);
Gcd(h,dh);

K:=GF(2^2);
gg:=PrimitiveElement(K);
Order(gg);

NumberPrimitiveElements:=function(q)
    local K;
    if not IsPrimePowerInt(q) then 
        Print(q," is not prime power \n");
    else 
        K:=GF(q);
        return Length(Filtered(K,x->Order(x)=q-1));
    fi;
end;    


MinimalPolynomial(K,gg);
MinimalPolynomial(GF(2),gg);

x:=Indeterminate(GF(2),"x");
List(F,a->IsIrreducible(x^2+a));
List(F,a->IsIrreducible(x^2+a*x+1));
F2:=Cartesian(F,F);
List(F2,a->IsIrreducible(x^2+a[1]*x+a[2]));

#######################
##Generators and cosets
#######################


########################
##Square free reduction
########################


p:=2;
x:=Indeterminate(GF(p),"x");
f:=(x^2+x+1)^2*(x^3+x+1)^3*(x+1);
df:=Derivative(f);
gdf:=Gcd(f,df);
qq:=Quotient(f,gdf);
coe:=CoefficientsOfUnivariatePolynomial(gdf);
Display(coe);
lcoe:=Length(coe);
slcoe:=QuoInt(lcoe,2)+1;
redgdf:=Sum(List([1..slcoe],i->coe[2*i-1]*x^(i-1)));

PartSquareFreeRed:=function(pol)
    df:=Derivative(f);
    gdf:=Gcd(f,df);
    return Quotient(f,gdf);
end;

PolypRootExtraction:=function(pol)
    coe:=CoefficientsOfUnivariatePolynomial(pol);
    bou:=QuoInt(Length(coe),p)+1;
    filcoe:=List([1..bou],i->coe[p*i-1]);
    if filcoe=List([1..bou],0);

#####################
##Berlekamp Algorithm
#####################

Berlekamp:=function(pol)
    qq:=PartSquareFreeRed(pol);
    deg:=Degree(qq);
    mat:=List([0..deg-1],i->CoefficientsOfUnivariatePolynomial(x^(p*i) mod qq));
    id:=IdentityMat(deg,GF(p));
    Q:=mat-id;
    rk:=Rank(Q);
    numirrfacs:=deg-rk;
    ker:=NullspaceMat(Q);
    vv:=Random(ker);
    uu:=Sum(List([0..deg-1],i->vv[i+1]*x^i));
    return List(GF(p),s->Gcd(qq,uu-s));
end;

#################
##Hensel Lifting
#################
    
y:=Indeterminate(Integers,"y");   
pol:=y^4+1;
p:=5;
x:=Indeterminate(GF(p),"x");
polmod:=x^4+1;
Factors(pol);
fac:=Factors(polmod);
f1mod:=fac[1];
g1mod:=fac[2];
Gcd(f1mod,g1mod)=1;
m1:=GcdRepresentation(f1mod,g1mod);
v1mod:=m1[1];
w1mod:=m1[2];

ConvertToIntPol:=function(pol)
    coe:=CoefficientsOfUnivariatePolynomial(pol);
    return Sum(List([0..Length(coe)-1],i->Int(coe[i+1])*y^i));
end;

f1:=ConvertToIntPol(f1mod);
g1:=ConvertToIntPol(g1mod);
a1:=(pol-f1*g1)/p;

v1:=ConvertToIntPol(v1mod);
w1:=ConvertToIntPol(w1mod);

f2:=f1+p*a1*w1;
g2:=g1+p*a1*v1;

a2:=(pol-f2*g2)/p^2;