Vincent Pilaud | Enseignement | Agrégation

> restart;

1. ALGORITHME DE WIEDEMANN SUR LE CORPS DES RÉELS

A. GÉNÉRATEUR D'UNE SUITE RÉCURRENTE LINÉAIRE

> calcultermes:=proc(p,s)
local i,d,r;
d:=nops(p); r:=s;
for i from 1 to d do
r:=[op(r),add(p[j]*r[i+j-1],j=1..d)];
end do;
return r;
end proc;

> generateur:=proc(s)
local R,R0,R1, V,V0,V1, Q, n;
n:=nops(s)/2;
R0:=X^(2*n); R1:=add(s[i]*X^(2*n-i),i=1..2*n);
V0:=0;V1:=1;
while degree(R1)>=n do
R:=rem(R0,R1,X,'Q'); V:=simplify(V0-Q*V1); V0:=V1; V1:=V; R0:=R1; R1:=R;
end do;
return(V1/lcoeff(V1));
end proc;

> sort(generateur(calcultermes([1,25,3,8],[12,0,4,57])),X);

B. ALGORITHME DE WIEDEMANN

> with(linalg):

Warning, the protected names norm and trace have been redefined and unprotected

> vecteuralea:=n->[seq(rand(1..100)(),i=1..n)];

> Wiedemann:=proc(A,b)
local n,B,d,P,U,Q;
n:=rowdim(A);d:=0;
B:=b;Q:=1;
while norm(B)<>0 do
        U:=vecteuralea(n);
        P:=generateur([seq(evalm(U&*A^i&*B)[1],i=1..2*(n-d))]);
        d:=degree(P)+d;Q:=P*Q;
        B:=evalm(value(Eval(expand(P),X=A))&*B);
end do;
return(Q);
end proc;

> M:=matrix(5,5,[[1,0,4,0,0],[0,0,0,2,0],[0,0,4,1,0],[2,0,0,3,1],[0,5,0,0,2]]);det(M); b:=matrix(5,1,[[1],[5],[0],[4],[2]]); Wiedemann(M,b);

C. SOLUTION AU SYSTÈME LINÉAIRE

> solution:=proc(A,b)
local P,Pmoins;
P:=Wiedemann(A,b);
Pmoins:=quo(P-tcoeff(P),tcoeff(P)*X,X);
return(evalm(-value(Eval(Pmoins,X=A))&*b));
end proc;

> S:=solution(M,b);

> evalm(M&*S);

2. LA MÊME CHOSE DANS LES CORPS FINIS

A. REPRÉSENTATION DES CORPS FINIS

On décrit le corps fini à p^n éléments comme le corps de rupture sur Fp d'un polynôme irréductible P de degrée n. Pour trouver un tel polynôme, on choisit aléatoirement un polynôme de degré n, puis on teste s'il est irréductible ou non.

> irred:=proc(p,n,R)
local d,i,q,L;
q:=p^n; d:=degree(R,Y);
if (Rem(Y^(q^d)-Y,R,Y) mod p)<>0 then return(false); end if;
L:=[op(numtheory[divisors](d))];
for i from 1 to nops(L)-1 do
if (Rem(Y^(q^i)-Y,R,Y) mod p)=0 then return(false); end if;
end do;
return(true);
end proc;

> polyalea:=(p,n)->add(rand(0..p-1)()*Y^i,i=0..n);

> representant:=proc(p,n)
local R;
R:=polyalea(p,n);
while degree(R,Y)<n or not(irred(p,1,R)) do R:=polyalea(p,n); end do;
return R;
end proc;

> representant(2,2);representant(3,2);

> t:=time():R:=representant(3,2);time()-t;

> inv:=proc(P,p,R)
Gcdex(P,R,Y,'s','t') mod p;
return(s);
end proc;

> inv(Y,3,Y^2+Y+2);inv(Y+1,3,Y^2+Y+2);inv(Y+2,3,Y^2+Y+2);

B. GÉNÉRATEUR D'UNE SUITE RÉCURRENTE LINÉAIRE

> calcultermesCF:=proc(l,s,p,P)
local i,d,r;
d:=nops(l);r:=s;
for i from 1 to d do
r:=[op(r),Rem(add(l[j]*r[i+j-1],j=1..d),P,Y) mod p];
end do;
return r;
end proc;

> generateurCF:=proc(s,p,P)
local R,R0,R1,Rpartiel, Q,Qpartiel, V,V0,V1, n;
n:=nops(s)/2;
R0:=X^(2*n);R1:=add(s[i]*X^(2*n-i),i=1..2*n);
V0:=0;V1:=1;
while degree(R1,X)>=n do
        Rpartiel:=sprem(R0,R1,X,'m','q')*inv(m,p,P); R:=Rem(Rpartiel,P,Y) mod p;
        Qpartiel:=q*inv(m,p,P); Q:=Rem(Qpartiel,P,Y) mod p;
        V:=Rem(V0-Q*V1,P,Y) mod p; V0:=V1; V1:=V; R0:=R1; R1:=R;
end do;
return(Rem(inv(lcoeff(V1,X),p,P)*V1,P,Y) mod p);
end proc;

> generateurCF(calcultermesCF([Y,2*Y,Y+1,2],[Y+2,1,Y,2*Y+1],3,Y^2+Y+2),3,Y^2+Y+2);

C. ALGORITHME DE WIEDEMANN

> elementaleaCF:=proc(p,P)
local d,i,u;
d:=degree(P,Y);u:=0;
for i from 0 to d-1 do u:=u+rand(0..p-1)()*Y^i;end do;
return u;
end proc;

> vecteuraleaCF:=proc(n,p,P)
local i,U;
U:=[];
for i from 1 to n do U:=[op(U),elementaleaCF(p,P)];end do;
return U;
end proc;

> WiedemannCF:=proc(A,b,p,P)
local n,B,d,R,U,Q;
n:=rowdim(A); d:=0;
B:=b; Q:=1;
while norm(B)<>0 do
        U:=vecteuraleaCF(n,p,P);
        R:=generateurCF([seq(Rem(evalm(U&*A^i&*B)[1],P,Y) mod p,i=1..2*(n-d))],p,P);
        d:=degree(R,X)+d;Q:=Rem(R*Q,P,Y) mod p;
        B:=map(z->Rem(z,P,Y) mod p,evalm(value(Eval(expand(R),X=A))&*B));
end do;
return(Q);
end proc;

> M:=matrix(5,5,[[Y,0,0,0,0],[0,0,2*Y,1,0],[0,0,2*Y+1,1,0],[2,0,0,Y+2,1],[0,2*Y+2,0,0,2]]);Rem(det(M),Y^2+Y+2,Y) mod 3; b:=matrix(5,1,[[1],[2*Y+2],[0],[Y+1],[2]]); WiedemannCF(M,b,3,Y^2+Y+2);

D. SOLUTION AU SYSTÈME LINÉAIRE

> solutionCF:=proc(A,b,p,P)
local R,Rmoins;
R:=WiedemannCF(A,b,p,P);
Rmoins:=Rem(quo(R-tcoeff(R,X),tcoeff(R,X)*X,X),P,Y) mod p;
return(map(z->Rem(z,P,Y) mod p,evalm(-value(Eval(Rmoins,X=A))&*b)));
end proc;

> S:=solutionCF(M,b,3,Y^2+Y+2);

> map(z->Rem(z,Y^2+Y+2,Y) mod 3,evalm(M&*S));

E. VITESSE DE CONVERGENCE

> WiedemannCFcompteur:=proc(A,b,p,P)
local n,B,d,R,U,Q,compteur;
n:=rowdim(A);d:=0;
B:=b;Q:=1;compteur:=0;
while norm(B)<>0 do
        U:=vecteuraleaCF(n,p,P);
        R:=BerlekampMasseyCF([seq(Rem(evalm(U&*A^i&*B)[1],P,Y) mod p,i=1..2*(n-d))],p,P);
        if degree(R,X)>0 then d:=degree(R,X)+d;Q:=Rem(R*Q,P,Y) mod p;
        B:=map(z->Rem(z,P,Y) mod p,evalm(value(Eval(expand(R),X=A))&*B));end if;
        compteur:=compteur+1;
end do;
return(compteur);
end proc;

> matricealeaCF:=proc(n,p,P)
local i,j,A;
A:=Matrix(n,n,0);
for i from 1 to n do
        for j from 1 to n do
                A[i,j]:=elementaleaCF(p,P);
        end do;
end do;
return A;
end proc;

> vitesseconv:=proc(N,n,p,P)
local A,b,m,i;
m:=0;
for i from 1 to N do
A:=matricealeaCF(n,p,P);b:=transpose([vecteuraleaCF(n,p,P)]);
if Rem(det(A),P,Y) mod p<>0 then m:=m+WiedemannCFcompteur(A,b,p,P); else i:=i-1; end if;
end do;
return (floor(m/N*100));
end proc;

> E:={seq([2^i,vitesseconv(20,5,2,representant(2,i))],i=1..4),seq([3^i,vitesseconv(20,5,3,representant(3,i))],i=1..3),seq([5^i,vitesseconv(20,5,5,representant(5,i))],i=1..2)};

> with(plots):pointplot(E);