/////////////////////////////////////////////////////////////////////////////// version="$Id$"; category="Teaching"; info=" LIBRARY: crypto.lib Procedures for teaching cryptography AUTHOR: Gerhard Pfister, pfister@mathematik.uni-kl.de OVERVIEW: The library contains procedures to compute the discrete logarithm, primality-tests, factorization included elliptic curves. The library is intended to be used for teaching purposes but not for serious computations. Sufficiently high printlevel allows to control each step, thus illustrating the algorithms at work. PROCEDURES: decimal(s); number corresponding to the hexadecimal number s exgcdN(a,n) compute s,t,d such that d=gcd(a,n)=s*a+t*n eexgcdN(L) T with sum_i L[i]*T[i]=T[n+1]=gcd(L[1],...,L[n]) gcdN(a,b) compute gcd(a,b) lcmN(a,b) compute lcm(a,b) powerN(m,d,n) compute m^d mod n chineseRem(T,L) compute x such that x = T[i] mod L[i] Jacobi(a,n) the generalized Legendre symbol of a and n primList(n) the list of all primes <=n primL(q) first primes p_1,...,p_r such that q=2;i--) { if((d mod i)==0){return(powerN(m,d/i,n)^i);} } return(m*powerN(m,d-1,n)); } for(i=12;i>=2;i--) { if((d mod i)==0){return(powerN(m,d/i,n)^i mod n);} } return(m*powerN(m,d-1,n) mod n); } example { "EXAMPLE:"; echo = 2; ring R = 0,x,dp; powerN(24,15,7); } proc chineseRem(list T,list L) "USAGE: chineseRem(T,L); RETURN: x such that x = T[i] mod L[i] NOTE: chinese remainder theorem EXAMPLE:example chineseRem; shows an example " { int i; int n=size(L); number N=1; for(i=1;i<=n;i++) { N=N*L[i]; } list M; for(i=1;i<=n;i++) { M[i]=N/L[i]; } list S=eexgcdN(M); number x; for(i=1;i<=n;i++) { x=x+S[i]*M[i]*T[i]; } x=x mod N; return(x); } example { "EXAMPLE:"; echo = 2; ring R = 0,x,dp; chineseRem(list(24,15,7),list(2,3,5)); } proc Jacobi(number a, number n) "USAGE: Jacobi(a,n); RETURN: the generalized Legendre symbol NOTE: if n is an odd prime then Jacobi(a,n)=0,1,-1 if n|a, a=x^2 mod n,else EXAMPLE:example Jacobi; shows an example " { int i; int z=1; number t=1; number k; if((((n-1)/2) mod 2)!=0){z=-1;} if(a<0){return(z*Jacobi(-a,n));} a=a mod n; if(n==1){return(1);} if(a==0){return(0);} while(a!=0) { while((a mod 2)==0) { a=a/2; if(((n mod 8)==3)||((n mod 8)==5)){t=-t;} } k=a;a=n;n=k; if(((a mod 4)==3)&&((n mod 4)==3)){t=-t;} a=a mod n; } if (n==1){return(t);} return(0); } example { "EXAMPLE:"; echo = 2; ring R = 0,x,dp; Jacobi(13580555397810650806,5792543); } proc primList(int n) "USAGE: primList(n); RETURN: the list of all primes <=n EXAMPLE:example primList; shows an example " { int i,j; list re; re[1]=2; re[2]=3; for(i=5;i<=n;i=i+2) { j=1; while(j<=size(re)) { if((i mod re[j])==0){break;} j++; } if(j==size(re)+1){re[size(re)+1]=i;} } return(re); } example { "EXAMPLE:"; echo = 2; list L=primList(100); size(L); L[size(L)]; } proc primL(number q) "USAGE: primL(q); RETURN: list of the first primes p_1,...,p_r such that q>p_1*...*p_(r-1) and qt)&&(m>s))||((mm) { s=(x-1)^2; } else { s=(x+1)^2; } } if(t>m){return(x-1);} if(s==m){return(x+1);} return(x); } example { "EXAMPLE:"; echo = 2; ring R = 0,x,dp; intRoot(20); } proc squareRoot(number a, number p) "USAGE: squareRoot(a,p); RETURN: the square root of a in Z/p, p prime NOTE: assumes the Jacobi symbol is 1 or p=2. EXAMPLE:example squareRoot; shows an example " { if(p==2){return(a);} if((a mod p)==0){return(0);} if(powerN(a,p-1,p)!=1) { "p is not prime"; return(number(-5)); } number n=random(1,2147483647) mod p; if(n==0){n=n+1;} number j=Jacobi(n,p); if(j==0) { "p is not prime"; return(number(-5)); } if(j==1) { return(squareRoot(a,p)); } number q=p-1; number e=0; number two=2; number z,m,t; while((q mod 2)==0) { e=e+1; q=q/2; } number y=powerN(n,q,p); number r=e; number x=powerN(a,(q-1)/2,p); number b=a*x^2 mod p; x=a*x mod p; while(((b-1) mod p)!=0) { m=0;z=b; while(((z-1) mod p)!=0) { z=z^2 mod p; m=m+1; } t=powerN(y,powerN(two,r-m-1,p),p); y=t^2 mod p; r=m; x=x*t mod p; b=b*y mod p; } return(x); } example { "EXAMPLE:"; echo = 2; ring R = 0,x,dp; squareRoot(8315890421938608,32003); } proc solutionsMod2(matrix M) "USAGE: solutionsMod2(M); RETURN: an intmat containing a basis of the vector space of solutions of the linear system of equations defined by M over the prime field of characteristic 2 EXAMPLE:example solutionsMod2; shows an example " { def R=basering; ring Rhelp=2,z,(c,dp); matrix M=imap(R,M); matrix S=syz(M); setring(R); matrix S=imap(Rhelp,S); int i,j; intmat v[nrows(S)][ncols(S)]; for(i=1;i<=nrows(S);i++) { for(j=1;j<=ncols(S);j++) { if(S[i,j]==1){v[i,j]=1;} } } return(v); } example { "EXAMPLE:"; echo = 2; ring R = 0,x,dp; matrix M[3][3]=1,2,3,4,5,6,7,6,5; solutionsMod2(M); } proc powerX(int q, int i, ideal I) "USAGE: powerX(q,i,I); RETURN: the q-th power of the i-th variable modulo I ASSUME: I is a standard basis EXAMPLE:example powerX; shows an example " { if(q<=181){return(reduce(var(i)^int(q),I));} if((q mod 5)==0){return(reduce(powerX(q div 5,i,I)^5,I));} if((q mod 4)==0){return(reduce(powerX(q div 4,i,I)^4,I));} if((q mod 3)==0){return(reduce(powerX(q div 3,i,I)^3,I));} if((q mod 2)==0){return(reduce(powerX(q div 2,i,I)^2,I));} return(reduce(var(i)*powerX(q-1,i,I),I)); } example { "EXAMPLE:"; echo = 2; ring R = 0,(x,y),dp; powerX(100,2,std(ideal(x3-1,y2-x))); } //====================================================================== //=========================== Discrete Logarithm ======================= //====================================================================== //============== Shank's baby step - giant step ======================== proc babyGiant(number b, number y, number p) "USAGE: babyGiant(b,y,p); RETURN: the discrete logarithm x: b^x=y mod p NOTE: This procedure works based on Shank's baby step - giant step method. EXAMPLE:example babyGiant; shows an example " { int i,j,m; list l; number h=1; number x; //choose m minimal such that m^2>p for(i=1;i<=p;i++){if(i^2>p) break;} m=i; //baby-step: compute the table y*b^i for 1<=i<=m for(i=1;i<=m;i++){l[i]=y*b^i mod p;} //giant-step: compute b^(m+j), 1<=j<=m and search in the baby-step table //for an i with y*b^i=b^(m*j). If found then x=m*j-i number g=b^m mod p; while(j=1;j--) { if(x[i]==x[j]) { s=(e[j]-e[i]) mod (p-1); t=(f[i]-f[j]) mod (p-1); if(s!=0) { i=0; } else { e[1]=0; f[1]=random(0,2147483629) mod (p-1); x[1]=powerN(b,f[1],p); i=1; } break; } } } list w=exgcdN(s,p-1); number u=w[1]; number d=w[3]; number a=(t*u/d) mod (p-1); while(powerN(b,a,p)!=y) { a=(a+(p-1)/d) mod (p-1); } return(a); } example { "EXAMPLE:"; echo = 2; ring R = 0,x,dp; number b=2; number y=10; number p=101; rho(b,y,p); } //==================================================================== //====================== Primality Tests ============================= //==================================================================== //================================= Miller-Rabin ===================== proc MillerRabin(number n, int k) "USAGE: MillerRabin(n,k); RETURN: 1 if n is prime, 0 else NOTE: probabilistic test of Miller-Rabin with k loops to test if n is prime. Using the theorem: If n is prime, n-1=2^s*r, r odd, then powerN(a,r,n)=1 or powerN(a,r*2^i,n)=-1 for some i EXAMPLE:example MillerRabin; shows an example " { if(n<0){n=-n;} if((n==2)||(n==3)){return(1);} if((n mod 2)==0){return(0);} int i; number a,b,j,r,s; r=n-1; s=0; while((r mod 2)==0) { s=s+1; r=r/2; } while(iN . N is prime if and only if for each prime factor p of A we can find a_p such that a_p^(N-1)=1 mod N and gcd(a_p^((N-1)/p)-1,N)=1 EXAMPLE:example PocklingtonLehmer; shows an example " { number m=intRoot(N); if(size(#)>0) { list S=PollardRho(N-1,10000,1,#); } else { list S=PollardRho(N-1,10000,1); } int i,j; number A=1; number p,a,g; list PA; list re; while(im){break;} while(1) { p=p*S[i+1]; if(((N-1) mod p)==0) { A=A*p; } else { break; } } i++; } if(A<=m) { A=N/A; PA=list(S[size(S)]); } for(i=1;i<=size(PA);i++) { a=1; while(a0) { L=#; } for(i=1;i<=size(L);i++) { if((n mod L[i])==0) { re[size(re)+1]=L[i]; while((n mod L[i])==0) { n=n/L[i]; } } if(n==1){return(re);} } int e=size(re); //here the rho-algorithm starts number a,d,x,y; while(n>1) { a=random(2,2147483629); x=random(2,2147483629); y=x; d=1; i=0; while(i1) { re[size(re)+1]=d; while((n mod d)==0) { n=n/d; } break; } if(i==k) { re[size(re)+1]=n; n=1; } } } if(allFactors) //want to obtain all prime factors { i=e; while(iB) break; while(w*P[i]<=B) { w=w*P[i]; } k=k*w; } number a=random(2,2147483629); number d=gcdN(powerN(a,k,n)-1,n); if((d>1)&&(d2)&&(Jacobi(n,p)==-1)){i++;continue;}//n is no quadratic rest mod p j=1; while(j<=p) { if(j>2*c+1) break; f=S[j][2]; v=S[j][3]; s=0; while((f mod p)==0) { s++; f=f/p; } if(s) { S[j][2]=f; v[i]=s; S[j][3]=v; l=j; while(l+p<=2*c+1) { l=l+p; f=S[l][2]; v=S[l][3]; s=0; while((f mod p)==0) { s++; f=f/p; } S[l][2]=f; v[i]=s; S[l][3]=v; } } j++; } } list T; for(j=1;j<=2*c+1;j++) { if((S[j][2]==1)||(S[j][2]==-1)) { T[size(T)+1]=S[j]; } } //the system of equations for the exponents {l_s} for the f(s) such //product f(s)^l_s is a square (l_s are 1 or 0) matrix M[k+1][size(T)]; for(j=1;j<=size(T);j++) { if(T[j][2]==-1){M[1,j]=1;} for(i=1;i<=k;i++) { M[i+1,j]=T[j][3][i]; } } intmat G=solutionsMod2(M); //construction of x and y such that x^2=y^2 mod n and d=gcd(x-y,n) //y=square root of product f(s)^l_s //x=product s+m number x=1; number y=1; for(i=1;i<=ncols(G);i++) { kill v; intvec v; v[k]=0; for(j=1;j<=size(T);j++) { x=x*T[j][1]^G[j,i] mod n; if((T[j][2]==-1)&&(G[j,i]==1)){y=-y;} v=v+G[j,i]*T[j][3]; } for(l=1;l<=k;l++) { y=y*B[l]^(v[l]/2) mod n; } d=gcdN(x-y,n); if((d>1)&&(d0) { s=#[1]; } else { s=1; } while(k1) { if(K[size(K)][2]!=K[size(K)-1][2]) { d=gcdN(K[size(K)][2],K[size(K)-1][2]); if(ellipticMult(q,a,b,K[size(K)],d)[3]==0) { K[size(K)][2]=K[size(K)-1][2]; } } } i=int(m); break; } } i=i+1; Q=ellipticAdd(q,a,b,mP,Q); } k++; } if(size(K)>0) { int te=1; for(i=1;i<=size(K)-1;i++) { if(K[size(K)][2]!=K[i][2]) { if(ellipticMult(q,a,b,K[i],K[size(K)][2])[3]!=0) { te=0; break; } } } if(te) { return(K[size(K)][2]); } } return(ShanksMestre(q,a,b,s)); } example { "EXAMPLE:"; echo = 2; ring R = 0,z,dp; ShanksMestre(32003,71,602); } //==================== Schoof's algorithm ============================= proc Schoof(number N,number a, number b) "USAGE: Schoof(N,a,b); RETURN: the number of points of the elliptic curve defined by y^2=x^3+a*x+b over Z/N NOTE: algorithm of Schoof EXAMPLE:example Schoof; shows an example " { int pr=printlevel; //test for ellictic curve number D=4*a^3+27*b^2; number G=gcdN(D,N); if(G==N){ERROR("not an elliptic curve");} if(G!=1){ERROR("not a prime");} //=== small N // if((N<=500)&&(pr<5)){return(countPoints(int(N),a,b));} //=== the general case number q=intRoot(4*N); list L=primL(2*q); int r=size(L); list T; int i,j; for(j=1;j<=r;j++) { T[j]=(testElliptic(int(N),a,b,L[j])+int(q)) mod L[j]; } if(pr>=5) { "==================================================================="; "Chinese remainder :"; for(i=1;i<=size(T);i++) { " x =",T[i]," mod ",L[i]; } "gives t+ integral part of the square root of q (to be positive)"; chineseRem(T,L); "we obtain t = ",chineseRem(T,L)-q; "==================================================================="; } number t=chineseRem(T,L)-q; return(N+1-t); } example { "EXAMPLE:"; echo = 2; ring R = 0,z,dp; Schoof(32003,71,602); } /* needs 518 sec Schoof(2147483629,17,3567); 2147168895 */ proc generateG(number a,number b, int m) "USAGE: generateG(a,b,m); RETURN: m-th division polynomial NOTE: generate the so-called division polynomials, i.e., the recursively defined polynomials p_m=generateG(a,b,m) in Z[x, y] such that, for a point (x:y:1) on the elliptic curve defined by y^2=x^3+a*x+b over Z/N the point@* m*P=(x-(p_(m-1)*p_(m+1))/p_m^2 :(p_(m+2)*p_(m-1)^2-p_(m-2)*p_(m+1)^2)/4y*p_m^3 :1) and m*P=0 if and only if p_m(P)=0 EXAMPLE:example generateG; shows an example " { poly f; if(m==0){return(f);} if(m==1){return(1);} if(m==2){f=2*var(1);return(f);} if(m==3){f=3*var(2)^4+6*a*var(2)^2+12*b*var(2)-a^2;return(f);} if(m==4) { f=4*var(1)*(var(2)^6+5*a*var(2)^4+20*b*var(2)^3-5*a^2*var(2)^2 -4*a*b*var(2)-8*b^2-a^3); return(f); } if((m mod 2)==0) { f=(generateG(a,b,m/2+2)*generateG(a,b,m/2-1)^2 -generateG(a,b,m/2-2)*generateG(a,b,m/2+1)^2) *generateG(a,b,m/2)/(2*var(1)); return(f); } f=generateG(a,b,(m-1)/2+2)*generateG(a,b,(m-1)/2)^3 -generateG(a,b,(m-1)/2-1)*generateG(a,b,(m-1)/2+1)^3; return(f); } example { "EXAMPLE:"; echo = 2; ring R = 0,(x,y),dp; generateG(7,15,4); } proc testElliptic(int q,number a,number b,int l) "USAGE: testElliptic(q,a,b,l); RETURN: an integer t, the trace of the Frobenius NOTE: the kernel for the Schoof algorithm: looks for the t such that for all points (x:y:1) in C[l]={P in C | l*P=0},C the elliptic curve defined by y^2=x^3+a*x+b over Z/q with group structure induced by 0=(0:1:0), (x:y:1)^(q^2)-t*(x:y:1)^q -ql*(x:y:1)=(0:1:0), ql= q mod l, trace of Frobenius. EXAMPLE:example testElliptic; shows an example " { int pr=printlevel; def R=basering; ring S=q,(y,x),lp; number a=imap(R,a); number b=imap(R,b); poly F=y2-x3-a*x-b; // the curve C poly G=generateG(a,b,l); ideal I=std(ideal(F,G)); // the points C[l] poly xq=powerX(q,2,I); poly yq=powerX(q,1,I); poly xq2=reduce(subst(xq,x,xq,y,yq),I); poly yq2=reduce(subst(yq,x,xq,y,yq),I); ideal J; int ql=q mod l; if(ql==0){ERROR("q is not prime");} int t; poly F1,F2,G1,G2,P1,P2,Q1,Q2,H1,H2,L1,L2; if(pr>=5) { "==================================================================="; "q=",q; "l=",l; "q mod l=",ql; "the Groebner basis for C[l]:";I; "x^q mod I = ",xq; "x^(q^2) mod I = ",xq2; "y^q mod I = ",yq; "y^(q^2) mod I = ",yq2; pause(); } //==== l=2 ============================================================= if(l==2) { xq=powerX(q,2,std(x3+a*x+b)); J=std(ideal(xq-x,x3+a*x+b)); if(deg(J[1])==0){t=1;} if(pr>=5) { "==================================================================="; "the case l=2"; "the gcd(x^q-x,x^3+ax+b)=",J[1]; pause(); } setring R; return(t); } //=== (F1/G1,F2/G2)=[ql](x,y) ========================================== if(ql==1) { F1=x;G1=1;F2=y;G2=1; } else { G1=reduce(generateG(a,b,ql)^2,I); F1=reduce(x*G1-generateG(a,b,ql-1)*generateG(a,b,ql+1),I); G2=reduce(4*y*generateG(a,b,ql)^3,I); F2=reduce(generateG(a,b,ql+2)*generateG(a,b,ql-1)^2 -generateG(a,b,ql-2)*generateG(a,b,ql+1)^2,I); } if(pr>=5) { "==================================================================="; "the point ql*(x,y)=(F1/G1,F2/G2)"; "F1=",F1; "G1=",G1; "F2=",F2; "G2=",G2; pause(); } //==== the case t=0 : the equations for (x,y)^(q^2)=-[ql](x,y) === J[1]=xq2*G1-F1; J[2]=yq2*G2+F2; if(pr>=5) { "==================================================================="; "the case t=0 mod l"; "the equations for (x,y)^(q^2)=-[ql](x,y) :"; J; "the test, if they vanish for all points in C[l]:"; reduce(J,I); pause(); } //=== test if all points of C[l] satisfy (x,y)^(q^2)=-[ql](x,y) //=== if so: t mod l =0 is returned if(size(reduce(J,I))==0){setring R;return(0);} //==== test for (x,y)^(q^2)=[ql](x,y) for some point J=xq2*G1-F1,yq2*G2-F2; J=std(J+I); if(pr>=5) { "==================================================================="; "test if (x,y)^(q^2)=[ql](x,y) for one point"; "if so, the Frobenius has an eigenvalue 2ql/t: (x,y)^q=(2ql/t)*(x,y)"; "it follows that t^2=4q mod l"; "if w is one square root of q mod l"; "t =2w mod l or -2w mod l "; "-------------------------------------------------------------------"; "the equations for (x,y)^(q^2)=[ql](x,y) :"; xq2*G1-F1,yq2*G2-F2; "the test if one point satisfies them"; J; pause(); } if(deg(J[1])>0) { setring R; int w=int(squareRoot(q,l)); setring S; //=== +/-2w mod l zurueckgeben, wenn (x,y)^q=+/-[w](x,y) //==== the case t>0 : (Q1/P1,Q2/P2)=[w](x,y) ============== if(w==1) { Q1=x;P1=1;Q2=y;P2=1; } else { P1=reduce(generateG(a,b,w)^2,I); Q1=reduce(x*G1-generateG(a,b,w-1)*generateG(a,b,w+1),I); P2=reduce(4*y*generateG(a,b,w)^3,I); Q2=reduce(generateG(a,b,w+2)*generateG(a,b,w-1)^2 -generateG(a,b,w-2)*generateG(a,b,w+1)^2,I); } J=xq*P1-Q1,yq*P2-Q2; J=std(I+J); if(pr>=5) { "==================================================================="; "the Frobenius has an eigenvalue, one of the roots of w^2=q mod l:"; "one root is:";w; "test, if it is the eigenvalue (if not it must be -w):"; "the equations for (x,y)^q=w*(x,y)";I;xq*P1-Q1,yq*P2-Q2; "the Groebner basis"; J; pause(); } if(deg(J[1])>0){return(2*w mod l);} return(-2*w mod l); } //==== the case t>0 : (Q1/P1,Q2/P2)=(x,y)^(q^2)+[ql](x,y) ===== P1=reduce(G1*G2^2*(F1-xq2*G1)^2,I); Q1=reduce((F2-yq2*G2)^2*G1^3-F1*G2^2*(F1-xq2*G1)^2-xq2*P1,I); P2=reduce(P1*G2*(F1-xq2*G1),I); Q2=reduce((xq2*P1-Q1)*(F2-yq2*G2)*G1-yq2*P2,I); if(pr>=5) { "we are in the general case:"; "(x,y)^(q^2)!=ql*(x,y) and (x,y)^(q^2)!=-ql*(x,y) "; "the point (Q1/P1,Q2/P2)=(x,y)^(q^2)+[ql](x,y)"; "Q1=",Q1; "P1=",P1; "Q2=",Q2; "P2=",P2; pause(); } while(t<(l-1)/2) { t++; //==== (H1/L1,H2/L2)=[t](x,y)^q =============================== if(t==1) { H1=xq;L1=1; H2=yq;L2=1; } else { H1=x*generateG(a,b,t)^2-generateG(a,b,t-1)*generateG(a,b,t+1); H1=subst(H1,x,xq,y,yq); H1=reduce(H1,I); L1=generateG(a,b,t)^2; L1=subst(L1,x,xq,y,yq); L1=reduce(L1,I); H2=generateG(a,b,t+2)*generateG(a,b,t-1)^2 -generateG(a,b,t-2)*generateG(a,b,t+1)^2; H2=subst(H2,x,xq,y,yq); H2=reduce(H2,I); L2=4*y*generateG(a,b,t)^3; L2=subst(L2,x,xq,y,yq); L2=reduce(L2,I); } J=Q1*L1-P1*H1,Q2*L2-P2*H2; if(pr>=5) { "we test now the different t, 00) { d=#[1]; } while(i<=d) { L=ellipticRandomCurve(N); if(L[3]>1){return(L[3]);} //the discriminant was not invertible P=list(0,L[2],1); j=0; M=1; while(jB) break; while(w*S[j](4-th root(N) +1)^2 - m*P=0=(0:1:0) - (m/q)*P=(x:y:z) and z a unit in Z/N Then N is prime. EXAMPLE:example ECPP; shows an example " { list L,S,P; number m,q; int i; number n=intRoot(intRoot(N)); n=(n+1)^2; //lower bound for q while(1) { L=ellipticRandomCurve(N); //a random elliptic curve C m=ShanksMestre(N,L[1],L[2],3); //number of points of the curve C S=PollardRho(m,10000,1); //factorization of m for(i=1;i<=size(S);i++) //search for q between the primes { q=S[i]; if(n0;i--) { if(s[i]=="a"){v[i]=0;} if(s[i]=="b"){v[i]=1;} if(s[i]=="c"){v[i]=2;} if(s[i]=="d"){v[i]=3;} if(s[i]=="e"){v[i]=4;} if(s[i]=="f"){v[i]=5;} if(s[i]=="g"){v[i]=6;} if(s[i]=="h"){v[i]=7;} if(s[i]=="i"){v[i]=8;} if(s[i]=="j"){v[i]=9;} if(s[i]=="k"){v[i]=10;} if(s[i]=="l"){v[i]=11;} if(s[i]=="m"){v[i]=12;} if(s[i]=="n"){v[i]=13;} if(s[i]=="o"){v[i]=14;} if(s[i]=="p"){v[i]=15;} if(s[i]=="q"){v[i]=16;} if(s[i]=="r"){v[i]=17;} if(s[i]=="s"){v[i]=18;} if(s[i]=="t"){v[i]=19;} if(s[i]=="u"){v[i]=20;} if(s[i]=="v"){v[i]=21;} if(s[i]=="w"){v[i]=22;} if(s[i]=="x"){v[i]=23;} if(s[i]=="y"){v[i]=24;} if(s[i]=="z"){v[i]=25;} if(s[i]==" "){v[i]=26;} } for(i=1;i<=size(s);i++) { n=n+v[i]*t^(i-1); } return(n); } static proc numberToWord(number n) { int i,j; string v; list s; number t=27; number mm; number nn=n; while(nn>t) { j++; mm=nn mod t; s[j]=mm; nn=(nn-mm)/t; } j++; s[j]=nn; for(i=1;i<=j;i++) { if(s[i]==0){v=v+"a";} if(s[i]==1){v=v+"b";} if(s[i]==2){v=v+"c";} if(s[i]==3){v=v+"d";} if(s[i]==4){v=v+"e";} if(s[i]==5){v=v+"f";} if(s[i]==6){v=v+"g";} if(s[i]==7){v=v+"h";} if(s[i]==8){v=v+"i";} if(s[i]==9){v=v+"j";} if(s[i]==10){v=v+"k";} if(s[i]==11){v=v+"l";} if(s[i]==12){v=v+"m";} if(s[i]==13){v=v+"n";} if(s[i]==14){v=v+"o";} if(s[i]==15){v=v+"p";} if(s[i]==16){v=v+"q";} if(s[i]==17){v=v+"r";} if(s[i]==18){v=v+"s";} if(s[i]==19){v=v+"t";} if(s[i]==20){v=v+"u";} if(s[i]==21){v=v+"v";} if(s[i]==22){v=v+"w";} if(s[i]==23){v=v+"x";} if(s[i]==24){v=v+"y";} if(s[i]==25){v=v+"z";} if(s[i]==26){v=v+" ";} } return(v); } proc code(string s) "USAGE: code(s); s a string ASSUME: s contains only small letters and space COMPUTE: a number, RSA-coding of the string s RETURN: return RSA-coding of the string s as string EXAMPLE: code; shows an example " { ring r=0,x,dp; number p=398075086424064937397125500550386491199064362342526708406385189575946388957261768583317; number q=472772146107435302536223071973048224632914695302097116459852171130520711256363590397527; number n=p*q; number phi=(p-1)*(q-1); number e=1234567891; list L=exgcdN(e,phi); number d=L[1]; number m=wordToNumber(s); number c=powerN(m,e,n); string cc=string(c); return(cc); } example {"EXAMPLE:"; echo = 2; string s="i go to school"; code(s); } proc decodeString(string g) "USAGE: decodeString(s); s a string ASSUME: s is a string of a number, the output of code COMPUTE: a string, RSA-decoding of the string s RETURN: return RSA-decoding of the string s as string EXAMPLE: decodeString; shows an example " { ring r=0,x,dp; number p=398075086424064937397125500550386491199064362342526708406385189575946388957261768583317; number q=472772146107435302536223071973048224632914695302097116459852171130520711256363590397527; number n=p*q; number phi=(p-1)*(q-1); number e=1234567891; list L=exgcdN(e,phi); number d=L[1]; execute("number c="+g+";"); number f=powerN(c,d,n); string s=numberToWord(f); return(s); } example {"EXAMPLE:"; echo = 2; string s="78638618599886548153321853785991541374544958648147340831959482696082179852616053583234149080198937632782579537867262780982185252913122030800897193851413140758915381848932565"; string t=decodeString(s); t; } /* //=============================================================== //======= Example for DSA ===================================== //=============================================================== Suppose a file test is given.It contains "Oscar". //Hash-function MD5 under Linux md5sum test 8edfe37dae96cfd2466d77d3884d4196 //================================================================ ring R=0,x,dp; number q=2^19+21; //524309 number o=2*3*23*number(7883)*number(16170811); number p=o*q+1; //9223372036869000547 number b=2; number g=power(2,o,p); //8308467587808723131 number a=111111; number A=power(g,a,p); //8566038811843553785 number h =decimal("8edfe37dae96cfd2466d77d3884d4196"); //189912871665444375716340628395668619670 h= h mod q; //259847 number k=123456; number ki=exgcd(k,q)[1]; //50804 //inverse von k mod q number r= power(g,k,p) mod q; //76646 number s=ki*(h+a*r) mod q; //2065 //========== signatur is (r,s)=(76646,2065) ===================== //==================== verification ============================ number si=exgcd(s,q)[1]; //inverse von s mod q number e1=si*h mod q; number e2=si*r mod q; number rr=((power(g,e1,p)*power(A,e2,p)) mod p) mod q; //76646 //=============================================================== //======= Example for knapsack ================================ //=============================================================== ring R=(5^5,t),x,dp; R; // # ground field : 3125 // primitive element : t // minpoly : 1*t^5+4*t^1+2*t^0 // number of vars : 1 // block 1 : ordering dp // : names x // block 2 : ordering C proc findEx(number n, number g) { int i; for(i=0;i<=size(basering)-1;i++) { if(g^i==n){return(i);} } } number g=t^3; //choice of the primitive root findEx(t+1,g); //2091 findEx(t+2,g); //2291 findEx(t+3,g); //1043 intvec b=1,2091,2291,1043; // k=4 int z=199; intvec v=1043+z,1+z,2091+z,2291+z; //permutation pi=(0123) v; 1242,200,2290,2490 //(1101)=(e_3,e_2,e_1,e_0) //encoding 2490+2290+1242=6022 und 1+1+0+1=3 //(6022,3) decoding: c-z*c'=6022-199*3=5425 ring S=5,x,dp; poly F=x5+4x+2; poly G=reduce((x^3)^5425,std(F)); G; //x3+x2+x+1 factorize(G); //[1]: // _[1]=1 // _[2]=x+1 // _[3]=x-2 // _[4]=x+2 //[2]: // 1,1,1,1 //factors x+1,x+2,x+3, i.e. (1110)=(e_pi(3),e_pi(2),e_pi(1),e_pi(0)) //pi(0)=1,pi(1)=2,pi(2)=3,pi(3)=0 gives: (1101) */