home *** CD-ROM | disk | FTP | other *** search
/ Power-Programmierung / CD1.mdf / basic / compiler / ubasic / malm / bwppt2.ub < prev    next >
Encoding:
Text File  |  1990-08-22  |  1.3 KB  |  34 lines
  1.    10   *BWppt2(N,&F)
  2.    20   ' Baillie-Wagstaff strong Lucas pseudoprime test.  See their paper.
  3.    30   ' We assume n is odd and > 11.  D is chosen by method "A*".
  4.    40   ' F is 1 to indicate probable prime, 0 for composite, and -1 if
  5.    50   ' N is a square.
  6.    60   ' 5 June 1990.
  7.    70   local P=1,Pinv,Q,Qq,D=-3,Us=0,Um,Ub=1,Vs,Vb,Wn,S,Wnn,J=-1,I
  8.    80   dim Bit%(400)
  9.    90   Wn=isqrt(N):if res=0 then F=-1:return endif
  10.   100   repeat D=-sgn(D)*(abs(D)+2) until kro(D,N)<1
  11.   110   if kro(D,N)=0 then F=0:return endif
  12.   120   Q=(1-D)\4:Wn=N+1:if Q=-1 then Q=5:P=5 endif
  13.   130   Wnn=gcd(Q,N):if and{Wnn>1,Wnn<N} then F=0:return endif
  14.   140   S=0:repeat inc S:Wn=Wn\2:until odd(Wn)
  15.   150   Wnn=Wn
  16.   160   repeat inc J:Wnn=Wnn\2:Bit%(J)=res until Wnn=0
  17.   170   Wnn=(N+1)\2:Pinv=modinv(P,N)
  18.   180   for I=J to 0 step -1
  19.   190   Vs=(2*Ub-P*Us)@N:Vb=((P*Vs+D*Us)*Wnn)@N
  20.   200   Us=(Us*Vs)@N:Ub=(Ub*Vb)@N
  21.   210   Um=((Ub+Q*Us)*Pinv)@N
  22.   220   if Bit%(I) then Us=Um else Ub=Um endif
  23.   230   next I
  24.   240   if Us=0 then F=1:return endif
  25.   250   Vs=(2*Ub-P*Us)@N
  26.   260   if Vs=0 then F=1:return endif
  27.   270   Qq=modpow(Q,Wn,N)
  28.   280   for I=1 to S-1
  29.   290   Vs=(Vs*Vs-2*Qq)@N
  30.   300   if Vs=0 then F=1:cancel for:return endif
  31.   310   Qq=(Qq*Qq)@N
  32.   320   next I
  33.   330   F=0:return ' End of subroutine BWppt2.
  34.