From Newsgroup: comp.lang.forth

Using the standard program dc on linux to test for primality.
This hinges on | that does modular exponentiation.

ciforth idiosyncrasies:
Uses &c to denote characters
Uses ^c to denote control characters
Uses &, in large numbers
Uses RDROP
Uses real strings "ORANG UTAN"
WANT announces a facility that you need
REGRESS introduces a test, use it as an example.
The latter two can be regarded as comment

You cannot get around the prefix h that introduce huge numbers.
You could however replace
'' h10 '' by '' 10 >S ''
for single precision numbers.

A working /dev/random is used, but you could replace RAND
by a suitable random number generation.
Even `` 117 CONSTANT RAND '' would work here.

Communication to dc is established using the following interface
With this dummy file it is supposed to compile.
------------------ 8< --------------------------

\ At this point we have available
\ command (sc -- ) Get a string plus cr to dc.
\ cmd1 (char -- ) Get a single char command to dc.
\ h PREFIX : number preceded by h is put on the dc stack.
\ get : get one string from dc, possibly containing &\..
\ Define dummies and kill the regression test.
: command 2DROP ;
: cmd1 DROP ;
: h NAME 2DROP ; IMMEDIATE PREFIX
: get 0 0 ;
: REGRESS POSTPONE \ ; IMMEDIATE

------------------ 8< --------------------------

This is the actual program.
<op>h is the operator in dc
h<number> is a huge number to dc.

------------------ 8< --------------------------
\ Copyright (2026): Albert van der Horst {by GNU Public License}

\ An interface with the program dc that can handle h-numbers (huge).

WANT BOUNDS TRUE REGRESS

INCLUDE pipe.frt

: $>S command ;
: >S 0 <# #S #> command ;
: S>$ "ps_" command get ;
: S> S>$ EVALUATE ;
: S>d 0.0 S>$ >NUMBER 2DROP ;

REGRESS 1234 >S S> S: 1234
REGRESS "1234" $>S S>$ EVALUATE S: 1234

: +h &+ cmd1 ;
: -h &- cmd1 ;
: *h &* cmd1 ;
: /h &/ cmd1 ;
: MODh &% cmd1 ;
: /MODh &~ cmd1 &r cmd1 ;
: **h &^ cmd1 ;
: */h "s_ * l_ /" command ;
: |h &| cmd1 ;

REGRESS h13 h3 +h S> S: 16
REGRESS h13 h3 -h S> S: 10
REGRESS h13 h3 *h S> S: 39
REGRESS h13 h3 /h S> S: 4
REGRESS h13 h3 MODh S> S: 1
REGRESS h13 h3 /MODh S> S> S: 4 1
REGRESS h13 h2 **h S> S: 169
REGRESS h13 h2 h70 |h S> S: 29
REGRESS h12 h2 h6 */h S> S: 4

: DROPh "s_" command ;
: SWAPh &r cmd1 ;
: DUPh &d cmd1 ;
: OVERh "sAsBlBlAlB" command ;
REGRESS h1 h2 h3 h4 SWAPh S> S> S> S> S: 3 4 2 1
REGRESS h1 h2 h3 h4 DROPh S> S> S> S: 3 2 1
REGRESS h1 h2 h3 h4 DUPh S> S> S> S> S> S: 4 4 3 2 1
REGRESS h1 h2 h3 h4 OVERh S> S> S> S> S> S: 3 4 3 2 1

: .h S>$ TYPE ;
REGRESS h1 h2 +h .h S:

\ On dc the .S (f) command fails if the stack is empty.
\ So first DEPTH is added (z) and later removed (s_).
: .Sh "zfs_ %n" command get TYPE ;

\ Testing huge primes.
"/dev/random" 0 OPEN-FILE THROW CONSTANT devrand
: RAND _ DSP@ 7 devrand READ-FILE THROW 7 <> 1001 ?ERROR ;
REGRESS RAND . S:


: inithex "16o" command "16i" command ;

1,000,000,000 CONSTANT infinity

\ For X leave odd X' and POWER of 2. (H:X -- H:X' ) ( -- p)
\ X = X' * 2^p
: even-norm infinity 0 DO
DUPh h2 MODh S> 1 = IF I LEAVE THEN h2 /h
LOOP ;
REGRESS h26 even-norm S> S: 1 13

\ The number on the DC stack is one. (Z is how many digits).
: one? DUPh "Z" command S> 1 > IF DROPh FALSE ELSE S> 1 = THEN ;
REGRESS h1234 one? S: FALSE
REGRESS h0 one? S: FALSE
REGRESS h1 one? S: TRUE

\ The number on the DC stack is zero.
: zero? DUPh "Z" command S> 1 > IF DROPh FALSE ELSE S> 0 = THEN ;
REGRESS h1234 zero? S: FALSE
REGRESS h0 zero? S: TRUE
REGRESS h1 zero? S: FALSE

\ For dc containing a number: "it IS probably prime using one probe."
\ (H:n -- ) ( -- flag)
: ?Rabbin DUPh "sp" command \ local p
h1 -h even-norm >R
RAND >S SWAPh "lp" command |h
DUPh one? IF DROPh RDROP TRUE EXIT THEN

BEGIN DUPh h1 +h "lp" command -h zero? IF DROPh RDROP TRUE EXIT THEN
R@ WHILE R> 1- >R
DUPh *h "lp" command MODh \ Square
REPEAT DROPh RDROP FALSE ;
REGRESS h1001 ?Rabbin S: FALSE
REGRESS h1013 ?Rabbin S: TRUE

: test 1,000,000,000,000,000,100 1,000,000,000,000,000,001 DO
I 0 <# #S #> command ?Rabbin IF I . CR THEN 2 +LOOP ;
\ Expect 03 09 31 79

: post-google
h10 h100 **h
BEGIN DUPh ?Rabbin 0= WHILE h1 +h REPEAT
.h ;

------------------ 8< --------------------------

Examples:
Finding primes inside a range:

S[ ] OK test
1000000000000000003
1000000000000000009
1000000000000000031
1000000000000000079

Find the first prime past google (10^100)
S[ ] OK post-google 100000000000000000000000000000000000000000000000000000000000000000000\ 00000000000000000000000000000267

P.S. using the rabbin oracle, see Knuth TAOC.

Groetjes Albert
--
The Chinese government is satisfied with its military superiority over USA.
The next 5 year plan has as primary goal to advance life expectancy
over 80 years, like Western Europe.
--- Synchronet 3.21a-Linux NewsLink 1.2

From Newsgroup: comp.lang.forth

In article <nnd$58658c13$6e229668@5e1e42a0add12094>,
<albert@spenarnc.xs4all.nl> wrote:

Using the standard program dc on linux to test for primality.
This hinges on | that does modular exponentiation.

ciforth idiosyncrasies:
Uses &c to denote characters
Uses ^c to denote control characters
Uses &, in large numbers
Uses RDROP
Uses real strings "ORANG UTAN"
WANT announces a facility that you need
REGRESS introduces a test, use it as an example.
The latter two can be regarded as comment

You cannot get around the prefix h that introduce huge numbers.
You could however replace
'' h10 '' by '' 10 >S ''
for single precision numbers.

A working /dev/random is used, but you could replace RAND
by a suitable random number generation.
Even `` 117 CONSTANT RAND '' would work here.

Communication to dc is established using the following interface
With this dummy file it is supposed to compile.
------------------ 8< --------------------------

\ At this point we have available
\ command (sc -- ) Get a string plus cr to dc.
\ cmd1 (char -- ) Get a single char command to dc.
\ h PREFIX : number preceded by h is put on the dc stack.
\ get : get one string from dc, possibly containing &\..
\ Define dummies and kill the regression test.
: command 2DROP ;
: cmd1 DROP ;
: h NAME 2DROP ; IMMEDIATE PREFIX
: get 0 0 ;
: REGRESS POSTPONE \ ; IMMEDIATE

------------------ 8< --------------------------

This is the actual program.
<op>h is the operator in dc
h<number> is a huge number to dc.

------------------ 8< --------------------------
\ Copyright (2026): Albert van der Horst {by GNU Public License}

\ An interface with the program dc that can handle h-numbers (huge).

WANT BOUNDS TRUE REGRESS

INCLUDE pipe.frt

: $>S command ;
: >S 0 <# #S #> command ;
: S>$ "ps_" command get ;
: S> S>$ EVALUATE ;
: S>d 0.0 S>$ >NUMBER 2DROP ;

REGRESS 1234 >S S> S: 1234
REGRESS "1234" $>S S>$ EVALUATE S: 1234

: +h &+ cmd1 ;
: -h &- cmd1 ;
: *h &* cmd1 ;
: /h &/ cmd1 ;
: MODh &% cmd1 ;
: /MODh &~ cmd1 &r cmd1 ;
: **h &^ cmd1 ;
: */h "s_ * l_ /" command ;
: |h &| cmd1 ;

REGRESS h13 h3 +h S> S: 16
REGRESS h13 h3 -h S> S: 10
REGRESS h13 h3 *h S> S: 39
REGRESS h13 h3 /h S> S: 4
REGRESS h13 h3 MODh S> S: 1
REGRESS h13 h3 /MODh S> S> S: 4 1
REGRESS h13 h2 **h S> S: 169
REGRESS h13 h2 h70 |h S> S: 29
REGRESS h12 h2 h6 */h S> S: 4

: DROPh "s_" command ;
: SWAPh &r cmd1 ;
: DUPh &d cmd1 ;
: OVERh "sAsBlBlAlB" command ;
REGRESS h1 h2 h3 h4 SWAPh S> S> S> S> S: 3 4 2 1
REGRESS h1 h2 h3 h4 DROPh S> S> S> S: 3 2 1
REGRESS h1 h2 h3 h4 DUPh S> S> S> S> S> S: 4 4 3 2 1
REGRESS h1 h2 h3 h4 OVERh S> S> S> S> S> S: 3 4 3 2 1

: .h S>$ TYPE ;
REGRESS h1 h2 +h .h S:

\ On dc the .S (f) command fails if the stack is empty.
\ So first DEPTH is added (z) and later removed (s_).
: .Sh "zfs_ %n" command get TYPE ;

\ Testing huge primes.
"/dev/random" 0 OPEN-FILE THROW CONSTANT devrand
: RAND _ DSP@ 7 devrand READ-FILE THROW 7 <> 1001 ?ERROR ;
REGRESS RAND . S:

: inithex "16o" command "16i" command ;

1,000,000,000 CONSTANT infinity

\ For X leave odd X' and POWER of 2. (H:X -- H:X' ) ( -- p)
\ X = X' * 2^p
: even-norm infinity 0 DO
DUPh h2 MODh S> 1 = IF I LEAVE THEN h2 /h
LOOP ;
REGRESS h26 even-norm S> S: 1 13

\ The number on the DC stack is one. (Z is how many digits).
: one? DUPh "Z" command S> 1 > IF DROPh FALSE ELSE S> 1 = THEN ;
REGRESS h1234 one? S: FALSE
REGRESS h0 one? S: FALSE
REGRESS h1 one? S: TRUE

\ The number on the DC stack is zero.
: zero? DUPh "Z" command S> 1 > IF DROPh FALSE ELSE S> 0 = THEN ;
REGRESS h1234 zero? S: FALSE
REGRESS h0 zero? S: TRUE
REGRESS h1 zero? S: FALSE

\ For dc containing a number: "it IS probably prime using one probe."
\ (H:n -- ) ( -- flag)
: ?Rabbin DUPh "sp" command \ local p
h1 -h even-norm >R
RAND >S SWAPh "lp" command |h
DUPh one? IF DROPh RDROP TRUE EXIT THEN

BEGIN DUPh h1 +h "lp" command -h zero? IF DROPh RDROP TRUE EXIT THEN
R@ WHILE R> 1- >R
DUPh *h "lp" command MODh \ Square
REPEAT DROPh RDROP FALSE ;
REGRESS h1001 ?Rabbin S: FALSE
REGRESS h1013 ?Rabbin S: TRUE

: test 1,000,000,000,000,000,100 1,000,000,000,000,000,001 DO
I 0 <# #S #> command ?Rabbin IF I . CR THEN 2 +LOOP ;
\ Expect 03 09 31 79

: post-google
h10 h100 **h
BEGIN DUPh ?Rabbin 0= WHILE h1 +h REPEAT
.h ;

------------------ 8< --------------------------

Examples:
Finding primes inside a range:

S[ ] OK test
1000000000000000003
1000000000000000009
1000000000000000031
1000000000000000079

Find the first prime past google (10^100)
S[ ] OK post-google >100000000000000000000000000000000000000000000000000000000000000000000\ >00000000000000000000000000000267

P.S. using the rabbin oracle, see Knuth TAOC.

Great disappointment. I thought that dc and python use the same
infinite precision libraries.
Apparently not. Using dc I am substantial slower than python.

Groetjes Albert
--
The Chinese government is satisfied with its military superiority over USA.
The next 5 year plan has as primary goal to advance life expectancy
over 80 years, like Western Europe.
--- Synchronet 3.21a-Linux NewsLink 1.2

From Newsgroup: comp.lang.forth

In article <nnd$27678a91$4a2a97f3@adcf8678945ca293>,
<albert@spenarnc.xs4all.nl> wrote:

In article <nnd$58658c13$6e229668@5e1e42a0add12094>,
<albert@spenarnc.xs4all.nl> wrote:

Using the standard program dc on linux to test for primality.
This hinges on | that does modular exponentiation.

<SNIP>

Great disappointment. I thought that dc and python use the same
infinite precision libraries.
Apparently not. Using dc I am substantial slower than python.

An investigating reveals that Python, in its kernel, implements infinite precision arithmetic using the Karatsuba algorithm. This uses Fourier transforms, and transforms a O(n^3) algorithm for exponentiation to
O(N. log(N)^m) where m small. Moreover, original thought for fp, there
is an integer version that Python uses. Division (needed for modulo computation) is approximately equal fast. (This plays a crucial role
in prime testing.)
I looked for libraries but it sits in the longobject.c file, sitting
in the Object directory where the ca 40 fundamental objects of Python
are defined.

(Note that the Great Internet Mersenne Prime Search uses floating point
taking advantage of fast fp engines.

Groetjes Albert

--
The Chinese government is satisfied with its military superiority over USA.
The next 5 year plan has as primary goal to advance life expectancy
over 80 years, like Western Europe.
--- Synchronet 3.21b-Linux NewsLink 1.2

From Newsgroup: comp.lang.forth

albert@spenarnc.xs4all.nl wrote:

In article <nnd$58658c13$6e229668@5e1e42a0add12094>,
<albert@spenarnc.xs4all.nl> wrote:

Using the standard program dc on linux to test for primality.
This hinges on | that does modular exponentiation.

ciforth idiosyncrasies:
Uses &c to denote characters
Uses ^c to denote control characters
Uses &, in large numbers
Uses RDROP
Uses real strings "ORANG UTAN"
WANT announces a facility that you need
REGRESS introduces a test, use it as an example.
The latter two can be regarded as comment

You cannot get around the prefix h that introduce huge numbers.
You could however replace
'' h10 '' by '' 10 >S ''
for single precision numbers.

A working /dev/random is used, but you could replace RAND
by a suitable random number generation.
Even `` 117 CONSTANT RAND '' would work here.

Communication to dc is established using the following interface
With this dummy file it is supposed to compile.
------------------ 8< --------------------------

\ At this point we have available
\ command (sc -- ) Get a string plus cr to dc.
\ cmd1 (char -- ) Get a single char command to dc.
\ h PREFIX : number preceded by h is put on the dc stack.
\ get : get one string from dc, possibly containing &\..
\ Define dummies and kill the regression test.
: command 2DROP ;
: cmd1 DROP ;
: h NAME 2DROP ; IMMEDIATE PREFIX
: get 0 0 ;
: REGRESS POSTPONE \ ; IMMEDIATE

------------------ 8< --------------------------

This is the actual program.
<op>h is the operator in dc
h<number> is a huge number to dc.

------------------ 8< --------------------------
\ Copyright (2026): Albert van der Horst {by GNU Public License}

\ An interface with the program dc that can handle h-numbers (huge).

WANT BOUNDS TRUE REGRESS

INCLUDE pipe.frt

: $>S command ;
: >S 0 <# #S #> command ;
: S>$ "ps_" command get ;
: S> S>$ EVALUATE ;
: S>d 0.0 S>$ >NUMBER 2DROP ;

REGRESS 1234 >S S> S: 1234
REGRESS "1234" $>S S>$ EVALUATE S: 1234

: +h &+ cmd1 ;
: -h &- cmd1 ;
: *h &* cmd1 ;
: /h &/ cmd1 ;
: MODh &% cmd1 ;
: /MODh &~ cmd1 &r cmd1 ;
: **h &^ cmd1 ;
: */h "s_ * l_ /" command ;
: |h &| cmd1 ;

REGRESS h13 h3 +h S> S: 16
REGRESS h13 h3 -h S> S: 10
REGRESS h13 h3 *h S> S: 39
REGRESS h13 h3 /h S> S: 4
REGRESS h13 h3 MODh S> S: 1
REGRESS h13 h3 /MODh S> S> S: 4 1
REGRESS h13 h2 **h S> S: 169
REGRESS h13 h2 h70 |h S> S: 29
REGRESS h12 h2 h6 */h S> S: 4

: DROPh "s_" command ;
: SWAPh &r cmd1 ;
: DUPh &d cmd1 ;
: OVERh "sAsBlBlAlB" command ;
REGRESS h1 h2 h3 h4 SWAPh S> S> S> S> S: 3 4 2 1
REGRESS h1 h2 h3 h4 DROPh S> S> S> S: 3 2 1
REGRESS h1 h2 h3 h4 DUPh S> S> S> S> S> S: 4 4 3 2 1
REGRESS h1 h2 h3 h4 OVERh S> S> S> S> S> S: 3 4 3 2 1

: .h S>$ TYPE ;
REGRESS h1 h2 +h .h S:

\ On dc the .S (f) command fails if the stack is empty.
\ So first DEPTH is added (z) and later removed (s_).
: .Sh "zfs_ %n" command get TYPE ;

\ Testing huge primes.
"/dev/random" 0 OPEN-FILE THROW CONSTANT devrand
: RAND _ DSP@ 7 devrand READ-FILE THROW 7 <> 1001 ?ERROR ;
REGRESS RAND . S:


: inithex "16o" command "16i" command ;

1,000,000,000 CONSTANT infinity

\ For X leave odd X' and POWER of 2. (H:X -- H:X' ) ( -- p)
\ X = X' * 2^p
: even-norm infinity 0 DO
DUPh h2 MODh S> 1 = IF I LEAVE THEN h2 /h
LOOP ;
REGRESS h26 even-norm S> S: 1 13

\ The number on the DC stack is one. (Z is how many digits).
: one? DUPh "Z" command S> 1 > IF DROPh FALSE ELSE S> 1 = THEN ;
REGRESS h1234 one? S: FALSE
REGRESS h0 one? S: FALSE
REGRESS h1 one? S: TRUE

\ The number on the DC stack is zero.
: zero? DUPh "Z" command S> 1 > IF DROPh FALSE ELSE S> 0 = THEN ;
REGRESS h1234 zero? S: FALSE
REGRESS h0 zero? S: TRUE
REGRESS h1 zero? S: FALSE

\ For dc containing a number: "it IS probably prime using one probe."
\ (H:n -- ) ( -- flag)
: ?Rabbin DUPh "sp" command \ local p
h1 -h even-norm >R
RAND >S SWAPh "lp" command |h
DUPh one? IF DROPh RDROP TRUE EXIT THEN

BEGIN DUPh h1 +h "lp" command -h zero? IF DROPh RDROP TRUE EXIT THEN >> R@ WHILE R> 1- >R
DUPh *h "lp" command MODh \ Square
REPEAT DROPh RDROP FALSE ;
REGRESS h1001 ?Rabbin S: FALSE
REGRESS h1013 ?Rabbin S: TRUE

: test 1,000,000,000,000,000,100 1,000,000,000,000,000,001 DO
I 0 <# #S #> command ?Rabbin IF I . CR THEN 2 +LOOP ;
\ Expect 03 09 31 79

: post-google
h10 h100 **h
BEGIN DUPh ?Rabbin 0= WHILE h1 +h REPEAT
.h ;

------------------ 8< --------------------------

Examples:
Finding primes inside a range:

S[ ] OK test
1000000000000000003
1000000000000000009
1000000000000000031
1000000000000000079

Find the first prime past google (10^100)
S[ ] OK post-google >>100000000000000000000000000000000000000000000000000000000000000000000\ >>00000000000000000000000000000267

P.S. using the rabbin oracle, see Knuth TAOC.

Great disappointment. I thought that dc and python use the same
infinite precision libraries.
Apparently not. Using dc I am substantial slower than python.

If you want fast bignum arithmetic use gmp. The relevant functions
are: 'mpn_mul' (multiplication), 'mpn_tdiv_qr' (quotient with reminder), 'mpn_gcd' (GCD) and possibly 'mpn_add' and 'mpn_sub'. For the last two
it is not hard to write reasonably good implementation in assembly,
and for big numbers most execution time goes into multiplicative
operations, that is why they are not critical. Other are pretty complex
to get good speed. Of course, you can use more functions to get
somewhat better speed.

Note that converting to text or from text has similar cost to
multiplication, so if you use text interface for each operation
you will _not_ get full speed.
--
Waldek Hebisch
--- Synchronet 3.21b-Linux NewsLink 1.2

From Newsgroup: comp.lang.forth

In article <10mb9ek$342ld$1@paganini.bofh.team>,
Waldek Hebisch <antispam@fricas.org> wrote:

albert@spenarnc.xs4all.nl wrote:

In article <nnd$58658c13$6e229668@5e1e42a0add12094>,
<albert@spenarnc.xs4all.nl> wrote:

Using the standard program dc on linux to test for primality.
This hinges on | that does modular exponentiation.

ciforth idiosyncrasies:
Uses &c to denote characters
Uses ^c to denote control characters
Uses &, in large numbers
Uses RDROP
Uses real strings "ORANG UTAN"
WANT announces a facility that you need
REGRESS introduces a test, use it as an example.
The latter two can be regarded as comment

You cannot get around the prefix h that introduce huge numbers.
You could however replace
'' h10 '' by '' 10 >S ''
for single precision numbers.

A working /dev/random is used, but you could replace RAND
by a suitable random number generation.
Even `` 117 CONSTANT RAND '' would work here.

Communication to dc is established using the following interface
With this dummy file it is supposed to compile.
------------------ 8< --------------------------

\ At this point we have available
\ command (sc -- ) Get a string plus cr to dc.
\ cmd1 (char -- ) Get a single char command to dc.
\ h PREFIX : number preceded by h is put on the dc stack.
\ get : get one string from dc, possibly containing &\..
\ Define dummies and kill the regression test.
: command 2DROP ;
: cmd1 DROP ;
: h NAME 2DROP ; IMMEDIATE PREFIX
: get 0 0 ;
: REGRESS POSTPONE \ ; IMMEDIATE

------------------ 8< --------------------------

This is the actual program.
<op>h is the operator in dc
h<number> is a huge number to dc.

------------------ 8< --------------------------
\ Copyright (2026): Albert van der Horst {by GNU Public License}

\ An interface with the program dc that can handle h-numbers (huge).

WANT BOUNDS TRUE REGRESS

INCLUDE pipe.frt

: $>S command ;
: >S 0 <# #S #> command ;
: S>$ "ps_" command get ;
: S> S>$ EVALUATE ;
: S>d 0.0 S>$ >NUMBER 2DROP ;

REGRESS 1234 >S S> S: 1234
REGRESS "1234" $>S S>$ EVALUATE S: 1234

: +h &+ cmd1 ;
: -h &- cmd1 ;
: *h &* cmd1 ;
: /h &/ cmd1 ;
: MODh &% cmd1 ;
: /MODh &~ cmd1 &r cmd1 ;
: **h &^ cmd1 ;
: */h "s_ * l_ /" command ;
: |h &| cmd1 ;

REGRESS h13 h3 +h S> S: 16
REGRESS h13 h3 -h S> S: 10
REGRESS h13 h3 *h S> S: 39
REGRESS h13 h3 /h S> S: 4
REGRESS h13 h3 MODh S> S: 1
REGRESS h13 h3 /MODh S> S> S: 4 1
REGRESS h13 h2 **h S> S: 169
REGRESS h13 h2 h70 |h S> S: 29
REGRESS h12 h2 h6 */h S> S: 4

: DROPh "s_" command ;
: SWAPh &r cmd1 ;
: DUPh &d cmd1 ;
: OVERh "sAsBlBlAlB" command ;
REGRESS h1 h2 h3 h4 SWAPh S> S> S> S> S: 3 4 2 1
REGRESS h1 h2 h3 h4 DROPh S> S> S> S: 3 2 1
REGRESS h1 h2 h3 h4 DUPh S> S> S> S> S> S: 4 4 3 2 1
REGRESS h1 h2 h3 h4 OVERh S> S> S> S> S> S: 3 4 3 2 1

: .h S>$ TYPE ;
REGRESS h1 h2 +h .h S:

\ On dc the .S (f) command fails if the stack is empty.
\ So first DEPTH is added (z) and later removed (s_).
: .Sh "zfs_ %n" command get TYPE ;

\ Testing huge primes.
"/dev/random" 0 OPEN-FILE THROW CONSTANT devrand
: RAND _ DSP@ 7 devrand READ-FILE THROW 7 <> 1001 ?ERROR ;
REGRESS RAND . S:


: inithex "16o" command "16i" command ;

1,000,000,000 CONSTANT infinity

\ For X leave odd X' and POWER of 2. (H:X -- H:X' ) ( -- p)
\ X = X' * 2^p
: even-norm infinity 0 DO
DUPh h2 MODh S> 1 = IF I LEAVE THEN h2 /h
LOOP ;
REGRESS h26 even-norm S> S: 1 13

\ The number on the DC stack is one. (Z is how many digits).
: one? DUPh "Z" command S> 1 > IF DROPh FALSE ELSE S> 1 = THEN ;
REGRESS h1234 one? S: FALSE
REGRESS h0 one? S: FALSE
REGRESS h1 one? S: TRUE

\ The number on the DC stack is zero.
: zero? DUPh "Z" command S> 1 > IF DROPh FALSE ELSE S> 0 = THEN ;
REGRESS h1234 zero? S: FALSE
REGRESS h0 zero? S: TRUE
REGRESS h1 zero? S: FALSE

\ For dc containing a number: "it IS probably prime using one probe."
\ (H:n -- ) ( -- flag)
: ?Rabbin DUPh "sp" command \ local p
h1 -h even-norm >R
RAND >S SWAPh "lp" command |h
DUPh one? IF DROPh RDROP TRUE EXIT THEN

BEGIN DUPh h1 +h "lp" command -h zero? IF DROPh RDROP TRUE EXIT THEN >>> R@ WHILE R> 1- >R
DUPh *h "lp" command MODh \ Square
REPEAT DROPh RDROP FALSE ;
REGRESS h1001 ?Rabbin S: FALSE
REGRESS h1013 ?Rabbin S: TRUE

: test 1,000,000,000,000,000,100 1,000,000,000,000,000,001 DO
I 0 <# #S #> command ?Rabbin IF I . CR THEN 2 +LOOP ;
\ Expect 03 09 31 79

: post-google
h10 h100 **h
BEGIN DUPh ?Rabbin 0= WHILE h1 +h REPEAT
.h ;

------------------ 8< --------------------------

Examples:
Finding primes inside a range:

S[ ] OK test
1000000000000000003
1000000000000000009
1000000000000000031
1000000000000000079

Find the first prime past google (10^100)
S[ ] OK post-google >>>100000000000000000000000000000000000000000000000000000000000000000000\ >>>00000000000000000000000000000267

P.S. using the rabbin oracle, see Knuth TAOC.

Great disappointment. I thought that dc and python use the same
infinite precision libraries.
Apparently not. Using dc I am substantial slower than python.

If you want fast bignum arithmetic use gmp. The relevant functions
are: 'mpn_mul' (multiplication), 'mpn_tdiv_qr' (quotient with reminder), >'mpn_gcd' (GCD) and possibly 'mpn_add' and 'mpn_sub'. For the last two
it is not hard to write reasonably good implementation in assembly,
and for big numbers most execution time goes into multiplicative
operations, that is why they are not critical. Other are pretty complex
to get good speed. Of course, you can use more functions to get
somewhat better speed.

Note that converting to text or from text has similar cost to
multiplication, so if you use text interface for each operation
you will _not_ get full speed.

I have extracted the preprocessed c-source from python.
The problem is not so much the entry address, but the data structure
that should be delivered to the procedures. Python insists that
single precision and multiple precision integers are the same objects.

I look into gpm. Maybe that is better.

Of course I do not converse numbers to readable output, certainly
not decimal.
dc and python has a single operator for modular exponentiation.

--
Waldek Hebisch

Groetjes Albert
--
The Chinese government is satisfied with its military superiority over USA.
The next 5 year plan has as primary goal to advance life expectancy
over 80 years, like Western Europe.
--- Synchronet 3.21b-Linux NewsLink 1.2