The one line factoring from lehs bignum

[Albert van der Horst]

There is an interesting factoring algo in lehs bignum package.
http://wrap.warwick.ac.uk/54707/

I rewrote it to not use floating point square root.
The problem with the usual integer square root is the initial value.
Here we calculate square roots that are near to each other.
Then you can use the previous result as an initial value, meaning
Newton with mostly two iterations.
I use a modified SQRT to tel whether something is a square : SQ?

\ --------------------------------------------------
1,000,000,013 CONSTANT p

WANT SQ SQRT FACTOR GCD PRIME?
WANT REGRESS

\ For N return FLOOR of the square root of n.
VARIABLE T 1000 T !
: SQRT1 >R
T @ R@ OVER / + 2/ \ Minimize iterations.
BEGIN R@ OVER / OVER + 2/ .S 2DUP > WHILE
NIP REPEAT DROP RDROP DUP T ! ;
REGRESS 1001 SQRT S: 31
REGRESS 81 SQRT S: 9

\ For N return :" n IS a square"
VARIABLE T_s 1000 T_s !
: SQ? >R
T_s @ R@ OVER / + 2/ \ Minimize iterations.
BEGIN R@ OVER / OVER + 2/ .S 2DUP > WHILE
NIP REPEAT DUP T_s ! * R> = ;
REGRESS 1001 SQ? S: 0
REGRESS 81 SQ? S: -1

: hart DUP 1 DO
KEY? IF RDROP LEAVE THEN
DUP I M* THROW SQRT1 1+ DUP >R \ s
SQ OVER MOD \ m
DUP SQ? .S IF
SQRT \ t
R> SWAP - OVER GCD LEAVE
ELSE
R> 2DROP
THEN
.S LOOP ;
--------------------------------------
We can see the slow start and the speedy followup :

1,000,000,000,000,001 hart

....S[ 16 10 ]
S[ 10 9 ]
S[ 9 9 ]SQ? 81 SQ? S: -1\ PASSED
OK

S[ 1000000000000001 500000000500 250000001249 ]
S[ 1000000000000001 250000001249 125000002624 ]
S[ 1000000000000001 125000002624 62500005311 ]
S[ 1000000000000001 62500005311 31250010655 ]
S[ 1000000000000001 31250010655 15625021327 ]
S[ 1000000000000001 15625021327 7812542663 ]
S[ 1000000000000001 7812542663 3906335331 ]
S[ 1000000000000001 3906335331 1953295662 ]
S[ 1000000000000001 1953295662 976903808 ]
S[ 1000000000000001 976903808 488963725 ]
S[ 1000000000000001 488963725 245504433 ]
S[ 1000000000000001 245504433 124788839 ]
S[ 1000000000000001 124788839 66401188 ]
S[ 1000000000000001 66401188 40730579 ]
S[ 1000000000000001 40730579 32641078 ]
S[ 1000000000000001 32641078 31638660 ]
S[ 1000000000000001 31638660 31622780 ]
S[ 1000000000000001 31622780 31622776 ]
S[ 1000000000000001 31622776 31622776 ]
S[ 1000000000000001 25191728 1399544 699780 ]
S[ 1000000000000001 25191728 699780 349907 ]
S[ 1000000000000001 25191728 349907 174989 ]
S[ 1000000000000001 25191728 174989 87566 ]
S[ 1000000000000001 25191728 87566 43926 ]
S[ 1000000000000001 25191728 43926 22249 ]
S[ 1000000000000001 25191728 22249 11690 ]
S[ 1000000000000001 25191728 11690 6922 ]
S[ 1000000000000001 25191728 6922 5280 ]
S[ 1000000000000001 25191728 5280 5025 ]
S[ 1000000000000001 25191728 5025 5019 ]
S[ 1000000000000001 25191728 5019 5019 ]
S[ 1000000000000001 25191728 0 ]
S[ 1000000000000001 ]

.....
S[ 1000000000000001 1368575950 1368575902 ]
S[ 1000000000000001 1368575902 1368575902 ]
S[ 1000000000000001 2272263536 47689 47668 ]
S[ 1000000000000001 2272263536 47668 47668 ]
S[ 1000000000000001 2272263536 0 ]
S[ 1000000000000001 ]
S[ 1000000000000001 1368941245 1368941196 ]
S[ 1000000000000001 1368941196 1368941196 ]
S[ 1000000000000001 843790935 32684 29250 ]
S[ 1000000000000001 843790935 29250 29048 ]
S[ 1000000000000001 843790935 29048 29048 ]
S[ 1000000000000001 843790935 0 ]
S[ 1000000000000001 ]
S[ 1000000000000001 1369306442 1369306393 ]
S[ 1000000000000001 1369306393 1369306393 ]
S[ 1000000000000001 649281361 25700 25481 ]
S[ 1000000000000001 649281361 25481 25481 ]
S[ 1000000000000001 649281361 -1 ] OK
.S

S[ 1000000000000001 456426971 ] OK
/MOD
OK
.
2190931 OK
.
0 OK

The last part shows that we get at the square root at the second
iteration, typically.
The total number of divisions (in behalf of sqrt) for this
factorisation is about 8000, which is substantially less than
by trial division.

Groetjes Albert
--
Albert van der Horst, UTRECHT,THE NETHERLANDS
Economic growth -- being exponential -- ultimately falters.
albert@spe&ar&c.xs4all.nl &=n http://home.hccnet.nl/a.w.m.van.der.horst