I've reorganized and extended some top level words:
create-set \ m n xt -- set
filter-set \ set1 xt -- set2
build-set \ m n xt -- set
transform-set \ set1 xt -- set2
This words should be used with the word :| to define and form sets.
1 10 :| b a | a ^2 b ^2 + isprime ; create-set zdup cr zet.
{(1,1),(2,1),(4,1),(6,1),(1,2),(3,2),(5,2),(7,2),(2,3),(8,3),(1,4),(5,4),(9,4),(2,5),(4,5),(6,5),(8,5),(1,6),(5,6),(2,7),(8,7),(3,8),(5,8),(7,8),(4,9)} ok
:| b a | a ^2 b ^2 + ; transform-set zdup cr zet.
{2,5,13,17,29,37,41,53,61,73,89,97,113} ok
:| n | n 4 mod 3 = ; filter-set zet. 0 ok
The word :| define a nameless word and count the number of parameters.
: bl# \ ad n -- m count the number of spaces in the string
over + swap 0 -rot
do i c@ bl = -
loop ;
variable loc# \ the number of parameters
variable sta# \ the number of outputs on the stack
: locals# \ --
>in @ >r
[char] | parse bl# loc# !
r> >in !
1 sta# ! ; immediate
: :| \ --
:noname
postpone locals#
postpone locals| ;
The nameless words, represented by xt on stack, could be of two types:
1. Taking parameters and leaving a flag.
2. Taking parameters and leaving a non negative integer.
In the first case a set with the "dimension" stated by the parameters is the result. This works with CREATE-SET and FILTER-SET.
In the second case a set of values is the result. This works with BUILD-SET and TRANSFORM-SET.
1 10 :| x | x ^2 ; build-set cr zet.
{1,4,9,16,25,36,49,64,81} ok
Normally these nameless words leave one integer on the stack, but for TRANSFORM-SET there is an option (using 2; instead of ;) when the word leaves two non negative integers on the stack
{ 0 10 | all }
creates the set {0,1,2,3,4,5,6,7,8,9} which can be transformed to a two dimensional set
:| n | n n ^2 2; transform-set cr zet.
{(0,0),(1,1),(2,4),(3,9),(4,16),(5,25),(6,36),(7,49),(8,64),(9,81)} ok
The purpose with this words is to be able to quickly inspect conjectures and Diophantine equations.
a and b are coprime if there exists values x and y such that
ax+by=1
Since BigZ just deals with sets of non negative numbers we can search x and y in the equation
ax-by=1
1 10 :| b a y x | a x * b y * - 1 = ; create-set cr zet.
{(2,1,1,1),(3,2,1,1),(4,3,1,1),(5,4,1,1),(6,5,1,1),(7,6,1,1),(8,7,1,1),(9,8,1,1),(1,1,2,1),(2,3,2,1),(3,5,2,1),(4,7,2,1),(5,9,2,1),(1,2,3,1),(2,5,3,1),(3,8,3,1),(1,3,4,1),(2,7,4,1),(1,4,5,1),(2,9,5,1),(1,5,6,1),(1,6,7,1),(1,7,8,1),(1,8,9,1),(3,1,1,2),(5,2,1,2),(7,3,1,2),(9,4,1,2),(1,1,3,2),(3,4,3,2),(5,7,3,2),(1,2,5,2),(3,7,5,2),(1,3,7,2),(1,4,9,2),(4,1,1,3),(7,2,1,3),(2,1,2,3),(5,3,2,3),(8,5,2,3),(1,1,4,3),(4,5,4,3),(7,9,4,3),(2,3,5,3),(5,8,5,3),(1,2,7,3),(4,9,7,3),(2,5,8,3),(5,1,1,4),(9,2,1,4),(3,2,3,4),(7,5,3,4),(1,1,5,4),(5,6,5,4),(3,5,7,4),(1,2,9,4),(6,1,1,5),(3,1,2,5),(8,3,2,5),(2,1,3,5),(7,4,3,5),(4,3,4,5),(9,7,4,5),(1,1,6,5),(6,7,6,5),(3,4,7,5),(2,3,8,5),(4,7,9,5),(7,1,1,6),(5,4,5,6),(1,1,7,6),(7,8,7,6),(8,1,1,7),(4,1,2,7),(5,2,3,7),(2,1,4,7),(9,5,4,7),(3,2,5,7),(6,5,6,7),(1,1,8,7),(8,9,8,7),(4,5,9,7),(9,1,1,8),(3,1,3,8),(5,3,5,8),(7,6,7,8),(1,1,9,8),(5,1,2,9),(7,3,4,9),(2,1,5,9),(4,3,7,9),(8,7,8,9)} ok
The vectors in the set are of the form (x,y,a,b), that is the opposite of its appearance as parameters. (Which is logical in Forth with postfix notations and values on a stack, but still a bit awkward).
The idea is a one row code and for that purpose there is a need of short words:
: isp isprime ; \ n -- flag
: isq sqr ; \ is perfect square: n -- flag
: isqf sqrfree ; \ is square free: n -- flag
: isem bigomega 2 = ; \ is semi prime: n -- flag
: ispp smallomega 1 = ; \ is prime power: n -- flag
: 2sqs 2sqsum ; \ square sum: a b -- sum
: 3sqs 3sqsum ; \ square sum: a b c -- sum
: 4sqs 4sqsum ; \ square sum: a b c d -- sum
: cop coprime ; \ are coprime: m n -- flag
: div swap mod 0= ; \ divides: m n -- flag
: << \ i j k -- flag i<j<k
over > -rot < and ;
: <<= \ i j k -- flag i<=j<=k
over >= -rot <= and ;
: z. zet. ;
: fi postpone else postpone false postpone then ; immediate
\ a short version of ELSE FALSE THEN.
When inspecting a Diophantine equation there might be symmetries and trivial cases to weed out.
1 100 :| c b a | a b c << a b cop and if a b 2sqs c ^2 = else false then ; create-set
can be shortened to
1 100 :| c b a | a b c << a b cop and if a b 2sqs c ^2 = fi ; create-set
cr zet.
{(3,4,5),(5,12,13),(8,15,17),(7,24,25),(20,21,29),(12,35,37),(9,40,41),(28,45,53),(11,60,61),(33,56,65),(16,63,65),(48,55,73),(36,77,85),(13,84,85),(39,80,89),(65,72,97)} ok
Experimenting with sets gives a great opportunity to find and test conjectures:
1 50 :| b a | a b 2sqsum isprime ; create-set ok
:| b a | a b + ; transform-set cr z.
{2,3,5,7,9,11,13,15,17,19,21,23,25,27,29,31,33,35,37,39,41,43,45,47,49,51,53,55,57,59,61,63,65,67,69,71,73,75,77,79,83,85,87,89,91,93,95} ok
This suggests that all odd numbers larger than 1 can be written as a sum a+b where a²+b² is prime.
Define
: sqeq \ m n a b -- ma²+nb²
dup * rot * -rot
dup * * + ;
and the set {a+b<50|a,b>0 & ma²+nb² is prime} for some different values of m and n
5 value m
7 value n
1 50 :| b a | m n a b sqeq isprime ; create-set
:| b a | a b + ; transform-set cr zet.
{5,7,11,13,17,19,23,25,29,31,35,37,41,43,47,49,53,55,59,61,65,67,71,73,77,79,83,85,89,91} ok
3 to m
4 to n
1 50 :| b a | m n a b sqeq isprime ; create-set
:| b a | a b + ; transform-set cr zet.
{2,3,4,5,6,8,9,10,11,12,13,15,16,17,18,19,20,22,23,24,25,26,27,29,30,31,32,33,34,36,37,38,39,40,41,43,44,45,46,47,48,50,51,52,53,54,55,57,58,59,60,61,62,64,65,66,67,68,69,71,72,73,74,75,76,78,79,80,81,82,83,86,87,88,89,93,94,95,96} ok
This will show a pattern which suggest that
{a+b|a,b>0 & ma²+nb² is prime}={k>1|gcd(k,m+n)=1}
However, further tests will show that some changes are needed:
Routines to check the conjecture and print the table:
: prime_partition \ m n k -- flag
false locals| flag k |
k 1
?do 2dup k i - i sqeq isprime
if true to flag leave then
loop 2drop flag ;
\ 100000 value nx Gives the same table as
1000 value nx
: maximal_exception \ m n -- mx
2dup + 1 locals| mx m+n n m |
nx 3
do m+n i coprime
if m n i prime_partition 0=
if i to mx then
then
loop mx ;
: table \ m2 n2 --
locals| n2 | cr
2 spaces dup 1 do i 3 .r loop 1
do n2 1 cr i 2 .r
do j i coprime
if j i maximal_exception
else 0
then 3 .r
loop
loop ;
The source code for BigZ can be loaded here:
https://github.com/Lehs/BigZ/blob/master/bigzet.txt