Suppose there are numbers a b c d on the stack, then

locals| d c b a |

pop the values on the stack and store them in the locals. Normally there is no real need of locals in Forth, when factoring optimally, but when the stack is used for counted number series

n1 n2 ... nk k

locals is handy. And of course, locals could be used to uncomplicate algorithms.

The Pollard rho factoring routine for single cell numbers is fast, even in 64 bit systems, and can be used to define number theoretical functions. In BigZ the word

pollard# \ n -- p1 ... pk k

factorize the number and present it in a form that can be sorted by the word

sort \ n1 ... nk k -- m1 ... mk k

: maxprimefactor \ n -- p

pollard# sort

over >r drops r> ;

drops \ n1 ... nk k --

The radical of a number can be defined easily by factoring, dropping all copies of prime factors and multiply the rest of the factors:

: radical \ n -- r

1 swap \ just a value to be dropped at the end

pollard# sort \ p1 ... pk k sorted with largest on top of stack

1 swap 0 \ p1 ... pk 1 k 0

do undequ \ is two primes eual?

if nip \ drop the first of them (second on the stack)

else * \ multiply single prime

then

loop nip ; \ drop the number "1" used by undequ

The word

undequ \ a b c -- a b c flag

compares the second and the third values on the stack and flag is true if a=b else false.

{ 1000000 1001000 | all } ok

utime function radical transform-set utime d- d. -118123 ok

That is, transforming the set of the numbers {1000000,1000001,...,1000999} to the set of their radicals takes about a tenth of a second.

Now zdup cardinality . cr zet. gives

1000

{10,42,1034,1158,1954,3910,4119, ... ,1000995,1000997,1000999} ok

(Non of the numbers appears to have the same radical).

Also, it is easy to define a test for square free numbers

: sqrfree \ n -- flag

dup radical = ;

that's fairly fast

utime { 1 10000 | sqrfree } utime d- d. -2750161 ok

zdup cardinality . 6083 ok

cr zet.

{1,2,3,5,6,7,10,11,13,14,15, ... ,9993,9994,9995,9997,9998} ok

A nice word to analyse a sorted counted bundle of numbers of the stack is

: hist \ a1 ... ak k -- a1 ... ai i ak nk

2dup 0 locals| n k1 a k |

begin dup a = k1 and

while n 1+ to n

k1 1- to k1 drop

repeat k n - a n ;

that counts the uppermost copies of the same number, leaving the remaining counted bundle under the histogram value ak nk on the top of the stack, indicating nk copies of the number ak.

For example, define a function theta that gives the greatest factor of n that is a sum of two squares.

Facts:

any prime of the form 4n+1 can be written as a sum of two squares;

the product of two squaresums is a squaresum;

for primes p of the form 4n+3, p^2i is of the form 0²+b².

The word squsumfac gives the contribution from the prime factor pk.

: squsumfac \ pk nk -- fac fac=a²+b²

over 3mod4

if dup odd if 1- then

then ** ;

: theta \ n -- m

dup 1 = if exit then

1 locals| m |

oddpart dup 1 = \ r s flag, where n=s*2^r, s odd.

if swap lshift exit then

pollard# sort

begin hist squsumfac

m * to m dup 0=

until drop m swap lshift ;

The sets in BigZ can't have big number members (yet) but it might be interesting to create sets of single number factors of big numbers.

\ testing for small (fast) single number divisors

\ of big number w in the intervall n,m-1

: sfac \ w -- w ?v | m n -- f

beven if 2drop 2 bdup b2/ exit then

0 locals| flag |

0 locals| flag |

do bdup i pnr@ bs/mod 0=

if i pnr@ to flag leave then bdrop

loop flag ;

: sfacset \ b -- set

0 \ count of the number of elements

begin pi_plim 1 sfac ?dup

while >zst 2 - bnip

repeat bdrop >zst reduce ;

Testing a conjecture about divisibility of Fibonacci numbers:

: bsfib \ n -- F(n) single input and big output

dup 2 < if s>b exit then

bzero bone 1

do btuck b+ loop bnip ;

Conjecture:

dup 2 < if s>b exit then

bzero bone 1

do btuck b+ loop bnip ;

Conjecture:

Any prime number p<n divide some Fibonacci number F(m), 0<m≤n.

: fibtest \ m n -- flag

false locals| flag |

do i pnr@ 1+ 1

do j pnr@ i bsfib sfacset smember

if true to flag leave then

loop flag if leave then

loop flag ;

pi_plim . 1077871 ok

pi_plim . 1077871 ok

utime pi_plim 1 fibtest . utime d- d. -1 -12836944 ok (t<13 sec).

That is, the conjecture is true for all primes Pn where n<1077871.

That is, the conjecture is true for all primes Pn where n<1077871.

__________

Due to Wikipedia there is a formula:

p|F(p-i), where i=(5/p), the Legendre symbol.

: test \ p -- flag

dup

dup 5 over legendre -

bsfib

bs/mod bdrop 0= ;

100000 random nextprime test . -1 ok

Due to Wikipedia there is a formula:

p|F(p-i), where i=(5/p), the Legendre symbol.

: test \ p -- flag

dup

dup 5 over legendre -

bsfib

bs/mod bdrop 0= ;

100000 random nextprime test . -1 ok

## No comments:

## Post a Comment