unit modulop;

interface

uses algebra;

type  intvect=array[0..max] of longint;
      matrix=array[0..max] of intvect;
      liste=record
             v:intvect;
             next:pointer;
            end;
       pliste=^liste;
     ttt=array[0..2] of longint;
     pttt=^ttt;

var     inv:pttt;
      prime:longint;

Procedure Bezout(a,b:longint;var d,u,v:longint);
function power (a,n,p:longint):longint;
function invnp(n,p:longint):longint;
function getinvtbl(p:longint):pttt;
procedure int_real(u:intvect;var v:vecteur);
Function Int_Degre(p:intvect):longint;
Function Int_Val(a:intvect;n,p:longint):longint;
procedure Div_Poly_mod(N,D:INTVECT;VAR Q,R:INTVECT);
procedure pgcd_poly_modp(a,b:intvect;var r:intvect);
Procedure Set_Prime(p:longint);
Procedure Delete_Prime;
function Ecris_Mod_p(a:intvect):string;

implementation

Procedure Bezout(a,b:longint;var d,u,v:longint);
 var t,r,k,n:longint;
           q:array[0..16] of longint;
 begin
  k:=0;
  repeat
   q[k]:=a div b;
   r:=a mod b;
   a:=b;
   b:=r;
   n:=k;
   inc(k);
  until r=0;
  d:=a;
  if n=0 then begin u:=0; v:=1; end else
  begin
   v:=-q[n-1]; u:=1;
   for k:=n-1 downto 1 do
    begin
     t:=u;
     u:=v;
     v:=t-v*q[k-1]
    end
  end
 end;

function power (a,n,p:longint):longint;
 var r:longint;
 begin
  a:=a mod p;
  r:=1;
  while n<>0 do
   begin
    if odd(n) then r:=(r*a) mod p;
    a:=sqr(a) mod p;
    n:=n shr 1
   end;
  power:=r
 end;

function invnp(n,p:longint):longint;
 var d,u,v:longint;
  begin
   Bezout(n,p,d,u,v);
   invnp:=(u+p) mod p;
  end;

function getinvtbl(p:longint):pttt;
 var  i:longint;
     ad:pttt;
 begin
  getmem(ad,(p+1)*4);
  ad^[0]:=0; ad^[1]:=1;
  for i:=2 to p-1 do
   ad^[i]:=invnp(i,p);
  getinvtbl:=ad;
 end;

procedure int_real(u:intvect;var v:vecteur);
 var i:byte;
 begin
  fillchar(v,sizeof(v),0);
  for i:=0 to max do v[i]:=u[i];
 end;

Function Int_Degre(p:intvect):longint;
var i:longint;
begin
  i:=Max;
  while (i>=0) and (p[i]=0) do i:=i-1;
  if i<0 then i:=-maxint;
  Int_Degre:=i;
end;


Function Int_Val(a:intvect;n,p:longint):longint;
 var r:longint;
     i:byte;
 Begin
  r:=0;
  for i:=Int_Degre(a) downto 0 do r:=(r*n+a[i]) mod p;
  Int_Val:=r;
 End;

PROCEDURE DIV_POLY_mod(N,D:INTVECT;VAR Q,R:INTVECT);
VAR DD,I,A,B:WORD;
    COEFF,s:longint;
 BEGIN
  R:=N;
  DD:=Int_DEGRE(D);
  FILLCHAR(Q,SIZEOF(Q),0);
  WHILE Int_DEGRE(R)>=DD DO
   BEGIN
    A:=Int_DEGRE(R);
    B:=Int_DEGRE(D);
    s:=inv^[( d[b] mod prime +prime) mod prime];
    COEFF:=(R[A]*s) mod prime;
    Q[A-B]:=COEFF;
    FOR I:=0 TO B DO
     R[I+A-B]:=(R[I+A-B]-COEFF*D[I]) mod prime;
   END;
 END;

procedure pgcd_poly_modp(a,b:intvect;var r:intvect);
 var q:intvect;
     i,j,x:longint;

 begin
  while Int_Degre(b)<>-Maxint do
   begin
    Div_poly_mod(a,b,q,r);
    a:=b;
    b:=r;
   end;
   fillchar(r,sizeof(r),0);
   j:=Int_Degre(a);
   x:=inv^[(a[j] mod prime+prime) mod prime];
   for i:=0 to j do r[i]:=((a[i]*x) mod prime+prime) mod prime;
  end;


Procedure Set_Prime(p:longint);
 begin
  prime:=p;
  inv:=getinvtbl(p);
 end;

Procedure Delete_Prime;
 begin
  freemem(inv,(prime+1)*4);
 end;

function Ecris_Mod_p(a:intvect):string;
 var e:vecteur;
 begin
  Int_Real(a,e);
  Ecris_Mod_p:=ecriture_poly(e);
 end;

 end.