2007-07-27, 10:32 | #1 |
Jul 2005
2^{3}·5 Posts |
Modification of Fermat's method
Greetings,
I believe I've discovered a way to replace one of the loops in Fermat's factorization method with a binary search. It does what was intended and seems to work, but is nothing more than a curiousity because the multiplication needed in the binary search is so much more computationally expensive than addition that the method is slower than Fermat's traditional method. The motivation for the technique described below was a reply I received to a question on the Yahoo PrimeNumbers forum some time ago at: http://tech.groups.yahoo.com/group/p.../message/18687 The specific implementation of Fermat's method I used is described in Knuth's "Art of Computer Programming" series, in Volume 2, "Seminumerical Algorithms", Algorithm C - Fermat's method, in section 4.5. Since I'm a programmer rather than a mathematician I would rather simply post the program itself rather than try to express it in terms of equations and botch things up. It's written in the Maple programming language. The derivation of the predictor equation "f(n)= r-n*(y+n-1)" used in function "f2" below is another subject altogether, but is outlined in the Yahoo PrimeNumbers message mentioned above. I will be happy to provide the Maple code for the derivation if anyone is interested. The "inner loop" that was replaced with a binary search is this code from the function "f1": while (r > 0) do r:= r - y: y:= y + 2: od: #================================================ # Maple Code Follows #================================================ Code:
# Ron S. Dotson (c) July 27, 2007 restart: # NOTE: "N" is the number to be factored kernelopts(ASSERT=true): #=============== Fermat Classic ================ f1:= proc(r0,x0,y0) local r, x, y, u, v, i, j; r:= r0: x:= x0: y:= y0: i:= 0: j:= 0: while (r <> 0) do while (r < 0) do r:= r + x - y: y:= y + 2: x:= x + 2: j:= j+1: od: while (r > 0) do r:= r - y: y:= y + 2: i:= i+1: od: od: u:= (x-y)/2; v:= (x+y-2)/2; return [u,v,i,j]: end proc: #=============== Fermat Modification ================ f2:= proc(r0,x0,y0) local f,r,x,y,y1,rNext,rUpper,rLower,rHalfRange, rEst, yEst, u,v,i,j; f:= proc(n, r, y) local fn; #p:= r-n*(y+n-1); fn:= r-n*(y+n-1); #f:= unapply(p,n): return(fn); end proc: r:= r0: x:= x0: y:= y0: i:= 0: j:= 0: rLower:= 0: rNext:= 0: rUpper:= 0: rEst:= 0: while (r <> 0) do y1:= y: if (r < 0) then r:= r + x: x:= x + 2: j:= j+1: fi: #printf("% 4d: rLower= % 4d, rNext= % 6d, ", i,rLower,rNext); #printf("r= % 6d, y= % 6d, rUpper= % 6d, rEst= % 6d, x= % 6d\n" ,r,y,rUpper,rEst,x); while (r > 0) do rUpper:= trunc(r/y)+1; # first estimate is an upper bound on "r" rLower:= 0; rHalfRange:= trunc((rUpper-rLower)/2): rNext:= rHalfRange: #printf("% 4d: rLower= % 4d, rNext= % 6d, ", i,rLower,rNext); #printf("r= % 6d, y= % 6d, rUpper= % 6d, rEst= % 6d\n" ,r,y,rUpper,rEst); while (rHalfRange > 0) do i:= i+1: rEst:= f(rNext, r, y): if (rEst < 0) then ASSERT(rNext<rUpper, "rNext NOT < rUpper"); rUpper:= rNext; # decrease upper bound else # rEst > 0 ASSERT(rNext>rLower, "rNext NOT > rLower"); rLower:= rNext; # increase lower bound fi: rHalfRange:= trunc((rUpper-rLower)/2): rNext:= rLower + rHalfRange: #printf("% 4d: rLower= % 4d, rNext= % 6d, ", i,rLower,rNext); #printf("r= % 6d, y= % 6d, rUpper= % 6d, rEst= % 6d\n" ,r,y,rUpper,rEst); od: #printf("\n"); if f(rLower, r, y)=0 then rNext:= rLower; else rNext:= rUpper; fi; r:= f(rNext, r, y); y:= y1 + 2*rNext; #ASSERT(r=0, "r NOT = 0"); #printf("% 4d: rLower= % 4d, rNext= % 6d, ", i,rLower,rNext); #printf("r= % 6d, y= % 6d, rUpper= % 6d, rEst= % 6d\n" ,r,y,rUpper,rEst); od: od: u:= (x-y)/2; v:= (x+y-2)/2; return [u,v,i,j]: end proc: #============ Begin Program =============== Digits:= 7: f1Total:= 0.0: f2Total:= 0.0: f1Totali:=0: f2Totali:= 0: for i from 1 to 10 do a:= nextprime(trunc(rand()/10^(12-Digits))); # rand() = 12 digits b:= nextprime(trunc(rand()/10^(12-Digits))); N:= a*b; printf("\n%2d Factors: a= %d, b= %d, N= %d\n", i, a, b, N); sqrtN:= trunc(evalf[800](sqrt(N))): # floor(sqrt(N)) y:= 1: x:= 2*sqrtN+1: # 2 times the integer square root of N, plus 1 r:= sqrtN^2 - N: # error in using sqrtN as a factor (negative) ASSERT(r<0,"Initialization1: r NOT < 0"); L1:= []; L2:= []; st:= time(): L1:= f1(r,x,y): et:= time(): f1secs:= et-st; f1Total:= f1Total + f1secs: f1Totali:= f1Totali + L1[3]: printf(" f1secs = % 5.1f secs, X-iterations= % 6d, Y-iterations= % 6d\n", f1secs,L1[4],L1[3]); st:= time(): L2:= f2(r,x,y): et:= time(): f2secs:= et-st; f2Total:= f2Total + f2secs: f2Totali:= f2Totali + L2[3]: printf(" f2secs = % 5.1f secs, X-iterations= % 6d, Y-iterations= % 6d\n", f2secs,L2[4],L2[3]); u:= L1[1]: v:= L1[2]: if (u*v <> N) then printf("%2d f1 FAILURE, a= %d, b= %d, N= %d\n", i,a,b,N); fi: x:= L2[1]: y:= L2[2]: if (x*y <> N) then printf("%2d f2 FAILURE, a= %d, b= %d, N= %d\n", i,a,b,N); fi: od: printf("TOTALS: f1= % 5.0f secs (%d Y-iterations), f2=% 5.0f secs (%d Y-iterations)\n", f1Total,f1Totali, f2Total,f2Totali); |
Thread Tools | |
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
Help with series convergence in Fermat method of factoring | EdH | Other Mathematical Topics | 52 | 2021-01-29 21:17 |
B1 and B2 in P-1 method | Miszka | Math | 13 | 2013-12-27 20:23 |
New Method | Unregistered | Miscellaneous Math | 14 | 2013-05-24 10:55 |
Elliptic Curve Method factoring - Fermat numbers | philmoore | Math | 131 | 2006-12-18 06:27 |
Fermat's factoring method with Gauss's inprovement | Carlos | Programming | 0 | 2005-09-11 12:50 |