version 1.1, 2001/10/02 11:16:59 |
version 1.2, 2002/09/11 07:26:45 |
|
|
squfof(n)= |
\\ return one non-trivial divisor of n > 1 using Shanks's SQUFOF |
|
squfof(n) = |
{ |
{ |
|
if (isprime(n) || issquare(n, &n), return(n)); |
|
|
if (isprime(n), return(n)); |
p = factor(n,0)[1,1]; |
if (issquare(n), return(sqrtint(n))); |
if (p != n, return(p)); |
|
|
p=factor(n,0)[1,1]; |
|
if (p!=n, return(p)); |
|
|
|
if (n%4==1, |
if (n%4==1, |
dd=n; d=sqrtint(dd); b=(((d-1)\2)<<1) + 1 |
D = n; d = sqrtint(D); b = (((d-1)\2) << 1) + 1 |
, |
, |
dd=n<<2; d=sqrtint(dd); b=(d\2)<<1 |
D = n << 2; d = sqrtint(D); b = (d\2) << 1 |
); |
); |
f=Qfb(1,b,(b^2-dd)>>2,0.); |
f = Qfb(1, b, (b^2-D)>>2); |
q=[]; lq=0; ii=0; l=sqrtint(d); |
l = sqrtint(d); |
|
|
|
q = []; lq = 0; i = 0; |
while (1, |
while (1, |
f=qfbred(f,3,dd,d); |
i++; |
ii++; a=component(f,1); |
f = qfbred(f, 3, D, d); |
if (!(ii%2), |
a = component(f, 1); |
if (issquare(a), |
if (!(i%2) && issquare(a, &as), |
as=sqrtint(a); j=1; |
j = 1; |
while (j<=lq, |
while (j<=lq, |
if (as==q[j], break); |
if (as == q[j], break); |
j++ |
j++ |
); |
|
if (j>lq, break) |
|
); |
); |
|
if (j > lq, break) |
); |
); |
|
|
if (abs(a)<=l, |
if (abs(a) <= l, |
q=concat(q,abs(a)); |
q = concat(q, abs(a)); |
print(q); lq++; |
print(q); lq++; |
) |
) |
); |
); |
|
|
gs=gcd(as,dd); |
print("i = ", i); print(f); |
print("i = ",ii); print(f); |
bb = component(f, 2); |
gs=gcd(gs,bb=component(f,2)); |
gs = gcd([as, bb, D]); |
if (gs>1, return (gs)); |
if (gs > 1, return (gs)); |
|
|
g=qfbred(Qfb(as,-bb,as*component(f,3),0.),3,dd,d); |
f = Qfb(as, -bb, as*component(f,3)); |
b=component(g,2); |
g = qfbred(f, 3, D, d); |
until (b1==b, |
b = component(g, 2); |
b1=b; g=qfbred(g,3,dd,d); |
until (b1 == b, |
b=component(g,2); |
b1 = b; g = qfbred(g, 3, D, d); |
|
b = component(g, 2); |
); |
); |
a=abs(component(g,1)); |
a = abs(component(g, 1)); |
if (!(a%2),a>>=1); |
if (a % 2, a, a>>1); |
a |
|
} |
} |