Words for testing conjectures

Making code for testing conjectures can be cumbersome even in the case of easy programming. So far in Zet there is the possibility of make and calculate with sets of numbers interactively.

{ 1 1000 | prime } ok
{ 1 1000 | pairprime } { 1 1000 | notpairprime } union ok
zet= . -1 ok

Conditions so far are

: all dup = ;
: odd 1 and ;
: 1mod4 4 mod 1 = ;
: 3mod4 4 mod 3 = ;
: sqr dup sqrtf dup * = ;
: sqrfree dup radical = ;
: pairprime dup prime over 2 + prime rot 2 - prime or and ; 
: notpairprime dup prime swap pairprime 0= and ;
: semiprime bigomega 2 = ;  \ A product of two primes?
: uniprime smallomega 1 = ; \ Only divisional by one prime?
: biprime smallomega 2 = ;  \ Exact two different primes?

The construction { 1 10000 | pairprime } is fancy but slow and risk overflow in data stack. All the pairprimes in the intervall will first be created on the stack and then be moved to the zst-stack. It's better to check number for number and create the set directly on the zst-stack.

: intcond \ low hi xt -- | -- s   "intervall condition"
  loc{ xt } 
  swap 0 -rot
  do i xt execute 
     if i >zst 1+ then
  loop 2* negate >zst ;

utime 1 100000 ' pairprime intcond utime cr d- d. cardinality .

-35954 2447  ok

A set of 2447 primes is created in about 0.04 seconds. This construction is also possible to use in definitions, then using ['] instead of '.

To filtrate a set on the zst-stack:

: setcond \ xt -- | s -- s'       "set condition"
  loc{ xt } 0
  foreach
  do zst> dup xt execute
     if >xst 1+ else drop then
  loop dup 0
  do xst> >zst 
  loop 2* negate >zst ;

{ 1 100 | prime } ' 1mod4 setcond cr zet.

{5,13,17,29,37,41,53,61,73,89,97} ok

It's also nice to be able to create the image of a function:

: intimage \ low hi xt -- | -- s  "intervall image"
  loc{ xt } 
  swap 2dup
  do i xt execute >zst
  loop - 2* negate >zst
  set-sort reduce ;

: setimage \ xt -- | s -- s'      "set image"
  loc{ xt } 0
  foreach 
  do zst> xt execute >xst 1+
  loop dup 0
  do xst> >zst
  loop 2* negate >zst
  set-sort reduce ;

Functions so far are: 

log~ ( n -- nr ) where nr=1+²log n

random ( u1 -- u2 ) where 0≤u2<u1
nextprime ( numb -- prime )

prevprime ( numb -- prime )

sqrtf ( m -- n ) "floor"

sqrtc ( m -- n ) "ceiling"

radical ( n -- r )
totients ( n -- t )

bigomega ( n -- b )

smallomega ( n -- s )
ufaculty ( u -- u! )
pnr@ ( n -- p ) prime number n

pi ( x -- n ) number of primes ≤ x

Functions and conditions both must have the stackdiagram ( m -- n ), but the concept will be generalized.

1 20 ' radical intimage zet. {1,2,3,5,6,7,10,11,13,14,15,17,19} ok

Some test functions:

: square dup * ;                \ x → x²
: sqr>prime square nextprime ;  \ x → nextprime(x²)
: sqr<prime square prevprime ;  \ x → prevprime(x²)
: foo dup totients mod ;        \ x → x(mod ϕ(x)) Euler's totient.

{ 1 100 | all } ' foo setimage cr zet.
{0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,21,22,23,25,27,31,33,35,39} ok

1 100 ' square intimage ' foo setimage cr zet.

{0,3,5,7,11,13,17,19,20,23,27,28,29,31,37,41,43,44,47,52,53,59,61,67,68,71,73,76,79,80,83,89,92,97,105,112,116,124,125,148,164,172,176,180,188,189,208,243,252,272,304,320,343,368,385,396,429,448,468,500,585,704,720,825,945,969,1008,1105,1197,1280,1309,1372,1540,1620,1701,1725,1729,1785,2185,2187,2625,2697,3069,3861} ok

Hmm, it seems like all odd primes less than 100 belongs to the image...

1 10000 ' square intimage ' foo setimage ok
1 10000 ' prime intcond ok
zswap diff zet. {2} ok

So I asked Mathematics stack exchange about it. (: 

Well, it might be sound to expect non dramatic explanations to conjectures, especially conjectures concerning primes.

To check relations n R m there is a need for testing subsets of Cartesian products, sets of pairs of integers.

: paircond \ xt -- | s -- s'
  loc{ xt } 0
  foreach
  do zdup zet> drop xt execute
     if zst xst setmove 1+ else zdrop then
  loop 6 * negate >xst
  xst zst setmove ;

{ 1 10 | all } zdup cartprod ' = paircond cr zet.

{(1,1),(2,2),(3,3),(4,4),(5,5),(6,6),(7,7),(8,8),(9,9)} ok

: pairimage \ xt -- | s -- s'
  loc{ xt } 0
  foreach
  do 1+ zet> drop xt execute >xst
  loop dup 0 
  do xst> >zst
  loop 2* negate >zst
  set-sort reduce ;

{ 2 10 | all } zdup cartprod ' * pairimage cr zet.

{4,6,8,9,10,12,14,15,16,18,20,21,24,25,27,28,30,32,35,36,40,42,45,48,49,54,56,63,64,72,81} ok

Some conditions and functions N²→N are =, <>, <, >, >=, <=, +, *, /, mod, **, ugcd, -, invmod, legendre, jacobi, kronecker, gnorm, choose, where m ≥ n for m n - and m n must be coprime for m n invmod. the gnorm of two integers m n is the norm of the gaussian integer m+in, that is the number m²+n².

: coprime ugcd 1 = ;

: divide swap mod 0= ; 

{ 1 10 | all } zdup cartprod ' coprime paircond cr zet.

{(1,1),(1,2),(1,3),(1,4),(1,5),(1,6),(1,7),(1,8),(1,9),(2,1),(2,3),(2,5),(2,7),(2,9),(3,1),(3,2),(3,4),(3,5),(3,7),(3,8),(4,1),(4,3),(4,5),(4,7),(4,9),(5,1),(5,2),(5,3),(5,4),(5,6),(5,7),(5,8),(5,9),(6,1),(6,5),(6,7),(7,1),(7,2),(7,3),(7,4),(7,5),(7,6),(7,8),(7,9),(8,1),(8,3),(8,5),(8,7),(8,9),(9,1),(9,2),(9,4),(9,5),(9,7),(9,8)} ok

{ 1 10 | all } zdup cartprod ' divide paircond cr zet.

{(1,1),(1,2),(1,3),(1,4),(1,5),(1,6),(1,7),(1,8),(1,9),(2,2),(2,4),(2,6),(2,8),(3,3),(3,6),(3,9),(4,4),(4,8),(5,5),(6,6),(7,7),(8,8),(9,9)} ok

{ 1 10 | all } zdup cartprod ' coprime paircond ' gnorm pairimage cr zet.

{2,5,10,13,17,25,26,29,34,37,41,50,53,58,61,65,73,74,82,85,89,97,106,113,130,145} ok

Simple graphs

The idea of using stacks for sets is not a bad idea - in my opinion. Besides from the simplified garbage collecton there is the benefit that the data can be accessed in arrays (in the stacks), which at least can be used to define some fast primitive routines. The set routines are not slow. When enumerating all the cards in a deck from 0 to 51, Zet can calculate the set of all 1326 possible hold cards in Texas hold´em poker in a few hundredths of a second. On my Android:

{ 0 52 | all } utime 2 power# utime 2swap d- cr d. cardinality cr .
29099 
1326

Here utime counts in μs.

Formally a simple graph is a set of vertices V and a subset of all unordered pairs of vertices E. Visually, a simple graph is a collection of verticies joined by zero or one edges and which consist no loops.

A subgraph (V,E) of a simple graph (V',E') is a graph such that V is a subset of V' and E of E'.

: subgraph \ -- flag | (V,E) (V',E') -- 
  unfence zrot unfence 
  zrot subset
  zswap subset and ;

There is a maximal subgraph for each subset V generated by E'.

\ E = intersection of E' and power#(2,V)
: edges~ \ E' V -- E
  2 power# intersection ;

But this straightforward implementation is inefficient. About 20 times faster is:

: edges \ E' V -- E
  0 >xst 
  zst yst setmove 
  foreach \ {u,v}∈E'
  ?do zdup unfence yzcopy1 member 
     if yzcopy1 member
        if zfence xzmerge
        else zdrop
        then
     else zst> drop zdrop
     then
  loop yst setdrop xst zst setmove ;

To make a random simple graph with v vertices and with an edge between two vertices in m cases of n:

: randgraph \ m n v -- | -- (V,E)
  loc{ m n v } 
  0 >xst
  { v 0 do i 1+ loop } 
  zdup 2 power# foreach \ {u,v}
  do n random m < 
     if zfence xzmerge
     else zdrop
     then
  loop xst zst setmove pair ;

The word extend creates a superset to the graph created by edges where all edges connected to V is submitted plus all edges connected to the submitted points.

\ V={x∈V'|y∈V" & {x,y}∈E'}
\ E={{x,y}∈E'|x∈V & y∈V}
: extend \ E' V" -- (V,E)
  zswap zst yst setmove 
  zst xst setcopy 
  foreach \ v∈V"
  do zst> yzcopy1
     begin zst@ 
     while zsplit zdup dup smember
        if xzmerge
        else zdrop
        then
     repeat zet> 2drop
  loop xst zst setmove 
  set-sort reduce
  yst zst setmove 
  zover edges pair ;

Counts all isolated points in a graph:
  
: isolated-vertices# \ -- n | (V,E) --
  unfence 0 dup loc{ flag }
  zst yst setmove
  foreach
  do zst> yzcopy1 true to flag
     begin zst@ flag and
     while zsplit dup smember 0= to flag
     repeat zdrop drop flag -
  loop yst zst setmove zdrop ;

4 5 9 randgraph zdup cr zet.

({1,2,3,4,5,6,7,8,9},{{8,9},{7,9},{6,9},{5,9},{4,9},{3,9},{2,9},{1,9},{6,8},{4,8},{1,8},{6,7},{5,7},{3,7},{2,7},{1,7},{5,6},{4,6},{3,6},{2,6},{1,6},{4,5},{3,5},{2,5},{3,4},{2,4},{1,4},{1,3}}) ok

isolated-vertices# . 0  ok

Counts all isolated components in a graph:

  
: components# \ -- n | (V,E) --
  zdup 0 >xst
  unfence 
  znip zst yst setcopy
  foreach
  do begin yzcopy1 zover
        extend unfence zdrop ztuck zet=
     until zfence xzmerge
  loop yst setdrop 
  xst zst setmove reduce cardinality
  isolated-vertices# + ;

Due to the formula for vertices, edges and components for a forest, that is, a graph without circuits, v=e+c:

: forest? \ -- flag | (V,E) -- 
  zdup unfence 
  cardinality \ e
  cardinality \ v
  components# \ c
  rot + = ;

4 5 9 randgraph zdup cr zet.
({1,2,3,4,5,6,7,8,9},{{7,9},{6,9},{5,9},{4,9},{3,9},{2,9},{6,8},{5,8},{4,8},{3,8},{2,8},{1,8},{5,7},{4,7},{3,7},{1,7},{5,6},{4,6},{3,6},{2,6},{1,6},{4,5},{3,5},{2,5},{1,5},{3,4},{2,4},{1,4},{2,3},{1,2}}) ok

forest? . 0  ok

\ Using set-sort to sort a vector
: vector-sort \ s -- s'
  set-sort zst> 1- >zst ;

\ check if E is a cycle 
: cycle \ -- flag | E --
  zdup multiunion
  zdup cardinality true loc{ v flag }
  zover zdup cardinality v = 0=
  if triplet zdrop false exit
  then pair components# 1 >
  if zdrop false exit
  then 0 >xst foreach
  do xzmerge
  loop xst zst setmove
  zet> cs sort 2 - 0
  do over = flag and to flag
     over > flag and to flag
  +loop = flag and ;

: clear-table \ s --
  pad 0 foreach
  do zst> max
  loop cells erase ;
 
: cyc!check \ n -- flag
  cells pad + 1 over +! @ 2 > ;

\ Test if (V,E) is 2-regular 
: 2-regular \ -- flag | (V,E) --
  unfence zswap clear-table
  begin zst@
  while zsplit unfence
     zst> cyc!check if zst> drop zdrop false exit then
     zst> cyc!check if zdrop false exit then
  repeat zdrop true ;

4 5 9 randgraph zdup cr zet.
({1,2,3,4,5,6,7,8,9},{{8,9},{7,9},{6,9},{4,9},{3,9},{2,9},{7,8},{6,8},{5,8},{4,8},{2,8},{1,8},{6,7},{4,7},{2,7},{1,7},{4,6},{2,6},{1,6},{4,5},{3,5},{1,5},{3,4},{1,4},{1,3},{1,2}}) ok

2-regular . 0  ok

Fast generation of the symmetric and alternating groups

From Wikipedia I got this algorithm, how to generate all permutation in alphabetical order:

1. Find the largest index k such that a[k]<a[k+1]. If no such index
   exists, the permutation is the last.
2. Find the largest index l greater than k such that a[k]<a[l].
3. Swap the value of a[k] with that of a[l].
4. Reverse the sequence from a[k+1] up to and including the final
   element a[n].

First a word that reverse the order of all n characters starting at address ad.

: reverse-string \ ad n --
  2dup + 1- loc{ ad1 n ad2 } n 2/ 0
  ?do ad1 i + c@ ad2 i - c@ 
     ad1 i + c! ad2 i - c!
  loop ; 

Then the 1'st part of the algorithm, returning the address corresponding to the index k if it exists or else return 0.

: lex-perm1 \ ad n -- a1
  0 loc{ a1 } 2 - over + 
  do i c@ i 1+ c@ <
     if i to a1 leave then -1
  +loop a1 ;

Find the largest address a2 greater than a1 such that [a1]<[a2].

: lex-perm2 \ ad n a1 -- a2
  0 loc{ a1 a2 } 1- over +
  do a1 c@ i c@ <
     if i to a2 leave then -1
  +loop a2 ;

Swap the values at addresses a1 and a2.

: lex-perm3 \ a1 a2 --
  over c@ over c@
  swap rot c!
  swap c! ;

Reverse the order of the last characters, from address a1 to the end.

: lex-perm4 \ ad n a1 -- 
  reverse from a1+1 to ad+n-1 
  1+ -rot            \ a1+1 ad n
  + over -           \ a1+1 ad+n-(a1+1) 
  reverse-string ; 

Calculate the next permutation:
  
: nextp \ ad n -- 
  2dup 2dup          \ ad n ad n ad n
  lex-perm1 dup 0=
  if 2drop 2drop drop exit 
  then dup >r        \ ad n ad n a1
  lex-perm2 r>       \ ad n a2 a1
  tuck swap          \ ad n a1 a1 a2
  lex-perm3          \ ad n a1
  lex-perm4 ;

Create the string 123...n:

: n>str \ n -- ad n
  dup 0 do i 49 + pad i + c! loop pad swap ;

Create a vector on the z-stack from the string.

: str>vect \ ad n -- | -- s
  loc{ ad n } n dup 0
  do ad i + c@ 15 and >zst loop 2* 1+ negate >zst ;

Fast calculation of the symmetry group of n! permutations.

: sym \ n -- | -- s
  n>str loc{ ad n }
  n dup ufaculty dup 0
  do ad n str>vect 
     ad n nextp
  loop swap 1+ * 2* negate >zst ;

utime 7 sym cardinality . utime d- d. 5040 -3931  ok

What would take hours with straight forward generation is now done in 4 milliseconds.

Next word calculates how many components in the vector s that is greater than the number m:

: perm> \ m -- n | s --
  loc{ m } 0
  foreach do zst> m > + loop negate ;

This is used to calculate the number of pairs of components in the vector s that is unsorted:
  
: #perm \ -- n | s -- 
  0
  begin zst@ -3 <
  while zsplit zst> zdup perm> +
  repeat zdrop ;

Which determine if the vector correspond to an odd permutation:

: oddperm \ -- flag | s --
  #perm 1 and ; 

: alt \ n -- | -- s
  n>str loc{ ad n }
  n dup ufaculty dup 0
  do ad n str>vect zdup oddperm
     if zdrop then ad n nextp
  loop swap 1+ * negate >zst ;

utime 7 alt cardinality . utime d- d. 2520 -35424  ok

To filter out the odd permutations takes some time, so the alternating group of n!/2 even permutations runs in 35 ms.

What is left is to figure out how to generate general groups fast. And to write a manual!

Tutorial 2: Euler project 1 - 5

https://projecteuler.net/problem=1

If we list all the natural numbers below 10 that are multiples of 3 or 5, we get 3, 5, 6 and 9. The sum of these multiples is 23.

Find the sum of all the multiples of 3 or 5 below 1000.

: multiple-of-3  3 mod 0= ;
: multiple-of-5  5 mod 0= ;

: setsum \ -- n | s --
  0 foreach do zst> + loop ;

{ 1 1000 | multiple-of-3 } { 1 1000 | multiple-of-5 } union setsum .

or

: multiple-of-3-or-5  dup multiple-of-3 swap multiple-of-5 or ;

{ 1 1000 | multiple-of-3-or-5 } setsum . 


https://projecteuler.net/problem=2

Each new term in the Fibonacci sequence is generated by adding the previous two terms. By starting with 1 and 2, the first 10 terms will be:

1, 2, 3, 5, 8, 13, 21, 34, 55, 89, ...

By considering the terms in the Fibonacci sequence whose values do not exceed four million, find the sum of the even-valued terms.

: next-fib-pair \ m n -- n m+n
  tuck + ;

: euler2 \ -- sum
  0 loc{ sum } 1 2
  begin dup 1 and 0=
     if dup sum + to sum 
     then next-fib-pair dup 4000000 >
  until 2drop sum ;

euler2 . 


https://projecteuler.net/problem=3

The prime factors of 13195 are 5, 7, 13 and 29.

What is the largest prime factor of the number 600851475143 ?

600851475143 pollard# sort over . drops 


https://projecteuler.net/problem=4

A palindromic number reads the same both ways. The largest palindrome made from the product of two 2-digit numbers is 9009 = 91 × 99.

Find the largest palindrome made from the product of two 3-digit numbers.

: palindrome \ n -- flag
  s>d <# #s #> 2dup + 1- true 
  loc{ add1 nr add2 flag } 
  nr 2/ 0
  do add1 i + c@ 
     add2 i - c@ = 0=
     if false to flag leave then
  loop flag ;

Just for curiosity:

{ 10000 1000000 | palindrome } cardinality . 1800  ok

: split-fact \ n -- i j where ij=n and i+j is minimal
  dup sqrtf 
  begin 2dup mod
  while 1-
  repeat tuck / ;

: euler4 \ -- n
  10000 999999 
  do i palindrome
     if i split-fact
        100 1000 within swap
        100 1000 within and
        if i leave then
     then -1
  +loop ;

euler4 . 


https://projecteuler.net/problem=5

2520 is the smallest number that can be divided by each of the numbers from 1 to 10 without any remainder.

What is the smallest positive number that is evenly divisible by all of the numbers from 1 to 20?

: divisible20 \ m -- flag
  true loc{ m flag } 21 1
  do m i mod
     if false to flag leave then
  loop flag ;

: euler5 \  -- n
  -1 1
  do i divisible20
     if i leave then
  loop ;

utime euler5 . utime d- d. xxxxxxxxx -36728698  ok

That takes more than 36.7 seconds. However, the wanted number is the product of all numbers in the intervall 1...20 that are of the form p^n, where p is a prime and p^n<=20<p^(n+1).

20 value numb
: pn \ m -- flag 
  dup uniprime 0= if drop false exit then 
  dup pollard# over >r drops r> * numb > ; 

: setmul \ -- n | s -- 
  1 foreach do zst> * loop ; 

utime { 2 21 | pn } setmul . utime d- d.  xxxxxxxxx -190  ok

which takes 190 micro seconds.

Some groups

First, I have improved the word psubgroups. Stupidly enough it calculated the whole (known) group. Now it works 5-10 times faster and calculate the set of subgroups of Sym(4) in a second. To calculate the set of subgroups to Sym(5) takes about 30 minutes, though, so enhancements must still be made.

Now handling groups with order around 100 and subgroups of groups with order around 30 is okay. A little more with a fast system...

Some finite groups

\ cyclic group of permutations of 1...n
: cyc \ n -- | -- s
  pcirc pgen ;

6 cyc zet. {(6,1,2,3,4,5),(5,6,1,2,3,4),(4,5,6,1,2,3),(3,4,5,6,1,2),(2,3,4,5,6,1),(1,2,3,4,5,6)} ok

\ symetric group of permutations of 1...n, n<6
: sym \ n -- | -- s
  dup 2 >
  if dup pcirc zfence proll zfence zmerge generate
  else 2 = if ( 2 1 ) pgen else ( 1 ) pgen then
  then ;

The n-th symmetry group is the set of all bijections of {1,...,n}.

4 sym cr zet.
{(4,2,3,1),(1,2,3,4),(1,3,4,2),(2,3,4,1),(2,4,3,1),(1,4,3,2),(2,1,3,4),(4,1,3,2),(3,4,1,2),(2,4,1,3),(3,4,2,1),(1,4,2,3),(3,1,4,2),(2,1,4,3),(3,2,4,1),(1,2,4,3),(4,3,1,2),(2,3,1,4),(4,2,1,3),(3,2,1,4),(4,3,2,1),(1,3,2,4),(4,1,2,3),(3,1,2,4)} ok

\ dihedral group of permutations of 1...n
: dih \ n -- | -- s
  dup >r pcirc zfence
  ( 1 r> ?do i -1 +loop ) zfence
  zetmerge generate ;

6 dih cr zet.
{(6,5,4,3,2,1),(6,1,2,3,4,5),(2,1,6,5,4,3),(2,3,4,5,6,1),(4,5,6,1,2,3),(4,3,2,1,6,5),(3,2,1,6,5,4),(3,4,5,6,1,2),(5,4,3,2,1,6),(5,6,1,2,3,4),(1,2,3,4,5,6),(1,6,5,4,3,2)} ok

There is also an other tradition where the dihedral group is denoted after the order of the group, but I think it's more consequent to denote it after the number of permutation elements.

Any permutation can be factorized in simple so called 2-cycles: (n m) where n maps to m and m maps to n. Example: (2,3,1)=(1 2)(1 3). Certain permutations can be factorized in an even number of 2-cycles and some can not. The product of even permutations is of course an even permutation and those permutation forms a subgroup Alt(S) of Sym(S). Both Sym(S) and Alt(S) can be generated by two elements, while their subgroups might not.

\ alternating group of permutations of 1...n, n<6
: alt \ n -- | -- s   n>2
  dup 3 = if drop ( 2 3 1 ) pgen exit then
  dup 1 and
  if >r
     { r@ pcirc
     ( r@ 2 - 1 do i loop r@ 1- r@ r> 2 - )
     } generate
  else >r
     { ( r@ 2 do i loop 1 r@ )
     ( r@ 2 - 1 do i loop r@ 1- r@ r> 2 - )
     } generate
  then ;

4 alt cr zet.
{(3,4,1,2),(3,1,2,4),(3,2,4,1),(4,1,3,2),(2,1,4,3),(1,3,4,2),(2,4,3,1),(1,4,2,3),(2,3,1,4),(4,3,2,1),(1,2,3,4),(4,2,1,3)} ok

\ quaternion group Q8={±1,±i,±j,±k} as group of permutations of 1..8
\ q8 \ -- s
: { ( 2 4 6 7 3 8 1 5 ) ( 3 5 4 8 7 2 6 1 ) } generate ;

The product of two permutation groups given as a permutation group

\ extend, to the right, bijection v to permute n elements 
: rext \ n -- | v -- v'
  >r ( r> zst@ cs do i 1+ loop )
  1 zst+! zswap 1 zst+! zswap zmerge -1 zst+! ;

\ extend to the left
: lext \ n -- | v -- v'
  dup >r ( r> zst@ cs - 1+ 1 do i loop ) 1 zst+!
  zswap zst@ tuck cs - loc{ x y }
  1 zst+! foreach
  do zst> y + loop x 1+ >>zst
  zmerge -1 zst+! ;

\ extend all functions in a set to the right 
: multirext \ n -- | s -- s'
  0 >xst foreach
  do dup rext zfence xzmerge
  loop drop xst zst setmove ;

3 sym 4 multirext cr zet.
{(3,2,1,4),(1,2,3,4),(2,3,1,4),(1,3,2,4),(3,1,2,4),(2,1,3,4)} ok

\ extend all to the left 
: multilext \ n -- | s -- s'
  0 >xst foreach
  do dup lext zfence xzmerge
  loop drop xst zst setmove ;

5 cyc 6 multilext cr zet.
{(1,2,3,4,5,6),(1,3,4,5,6,2),(1,4,5,6,2,3),(1,5,6,2,3,4),(1,6,2,3,4,5)} ok

\ the product of two groups s and s'
: gprod \ s s' -- sxs'
  ord zswap ord +
  dup multirext zswap
  multilext union generate ;

2 cyc zdup gprod zet. {(2,1,3,4),(1,2,3,4),(1,2,4,3),(2,1,4,3)} ok

3 cyc 4 alt gprod cr zet.
{(2,3,1,4,5,6,7,8,9,10,11,15,12,14,13),(2,3,1,4,5,6,7,8,9,10,11,15,14,13,12),(2,3,1,4,5,6,7,8,9,10,11,12,14,15,13),(2,3,1,4,5,6,7,8,9,10,11,13,14,12,15),(2,3,1,4,5,6,7,8,9,10,11,14,15,12,13),(2,3,1,4,5,6,7,8,9,10,11,13,12,15,14),(2,3,1,4,5,6,7,8,9,10,11,14,12,13,15),(2,3,1,4,5,6,7,8,9,10,11,13,15,14,12),(2,3,1,4,5,6,7,8,9,10,11,12,15,13,14),(2,3,1,4,5,6,7,8,9,10,11,14,13,15,12),(2,3,1,4,5,6,7,8,9,10,11,12,13,14,15),(1,2,3,4,5,6,7,8,9,10,11,12,14,15,13),(1,2,3,4,5,6,7,8,9,10,11,12,15,13,14),(1,2,3,4,5,6,7,8,9,10,11,12,13,14,15),(1,2,3,4,5,6,7,8,9,10,11,14,15,12,13),(1,2,3,4,5,6,7,8,9,10,11,13,15,14,12),(1,2,3,4,5,6,7,8,9,10,11,15,12,14,13),(1,2,3,4,5,6,7,8,9,10,11,13,14,12,15),(1,2,3,4,5,6,7,8,9,10,11,15,14,13,12),(1,2,3,4,5,6,7,8,9,10,11,13,12,15,14),(1,2,3,4,5,6,7,8,9,10,11,14,12,13,15),(1,2,3,4,5,6,7,8,9,10,11,14,13,15,12),(3,1,2,4,5,6,7,8,9,10,11,12,14,15,13),(3,1,2,4,5,6,7,8,9,10,11,12,15,13,14),(3,1,2,4,5,6,7,8,9,10,11,12,13,14,15),(3,1,2,4,5,6,7,8,9,10,11,14,15,12,13),(3,1,2,4,5,6,7,8,9,10,11,13,15,14,12),(3,1,2,4,5,6,7,8,9,10,11,15,12,14,13),(3,1,2,4,5,6,7,8,9,10,11,13,14,12,15),(3,1,2,4,5,6,7,8,9,10,11,15,14,13,12),(3,1,2,4,5,6,7,8,9,10,11,13,12,15,14),(3,1,2,4,5,6,7,8,9,10,11,14,12,13,15),(3,1,2,4,5,6,7,8,9,10,11,15,13,12,14),(3,1,2,4,5,6,7,8,9,10,11,14,13,15,12),(2,3,1,4,5,6,7,8,9,10,11,15,13,12,14),(1,2,3,4,5,6,7,8,9,10,11,15,13,12,14)} ok

Pseudo isomorphism test

\ the set of cyclic subgroups of the group s
: pcsubs \ s -- s'
  0 >xst
  foreach
  do pgen zfence xzmerge
  loop xst zst setmove reduce ;

\ flag true if not equal cardinality
: card<> \ -- flag | s s' -- s s'
  zover cardinality zdup cardinality = 0= ;

\ sort a list of non-negative integers
: vect-sort \ v -- v'
  set-sort zst> 1- > zst ;

\ compute vector of orders of all cyclic subgroups in s
: pscan \ s -- v
  0 foreach do cardinality swap 1+ loop
  sort 2* 1+ negate >zet ;

: pseudoiso \ -- flag | s s' --
  card<>
  if zdrop zdrop false exit then
  pcsubs zswap pcsubs card<>
  if zdrop zdrop false exit then
  pscan zswap pscan vector= ;

4 dih pcsubs pscan zet. (1,2,2,2,2,2,4) ok

4 dih 8 cyc pseudoiso . 0  ok

3 sym 3 dih pseudoiso . -1  ok

Due to Mathematics Stack Exchange the psudoiso test holds for all subgroups of Sym(7) but already in Sym(8) there are counterexamples. A counterexample in Sym(16) is also:

4 cyc zdup gprod 2 cyc q8 gprod pseudoiso . -1  ok

Tutorial: play around 1

In addition to the two Forth stacks s and r, there are three stacks x, y and z in Zet. The main stack (set parameter stack) for bundles is z. Most of the algebraic is done in z, for example:

reduce ( -- ) that eliminates copies of members in sorted sets in z;
zdup zdrop zover zswap znip ztuck and zrot manipulates bundles in z;
cardinality ( -- n | s -- ) that counts the number of elements in sets or components vectors;
foreach ( -- n 0 | s -- z1...zn ) that "appends" a set and prepare for a do loop for each element;
zet. ( s -- ) that prints the set/vector on z;
subset zet= member that examines sets and vectors on z;
union intersection diff powerset cartprod etc that works on z;
set-sort ( -- | m1...mk -- n1---nk ) 
zmerge ( s s' -- s" )

The other two stacks, x and y, are help stacks that compensate the lack of variables. Some operations working on x and y are:

setdup ( ad -- | obj -- obj obj )
setdrop ( ad -- | obj -- )
setover ( ad -- | obj1 obj2 -- obj1 obj2 obj1
setmove ( ad1 ad2 -- )
setcopy ( ad1 ad2 -- )
_fence ( ad -- | obj -- {obj} )
_split ( ad -- | s -- s' obj )  ad=yst,zst 

And there are also some special words for stack interaction:

xzmerge ( s -- ) takes set from z and merge so set in x
yzcopy1 ( -- s ) copy set from top y to z
yzcopy2 ( -- s ) copy set from next after top y to z

There are two main methods to penetrate a set, to use foreach or to use zsplit (or ysplit). When all elements are to be penetrated foreach is handy, but no objects under the top set on z can be reached under the penetration. When using zsplit, that splits a set in the top element and the rest of the set, the sets under can be reached, and the penetration can be abrupt without stack problems.

The facility with | in

{ 1 100 | prime } cr zet.{2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97} ok

does only work in interpretation mode and can't be compiled. 

1 value num

: coprime \ n -- flag 

  num ugcd 1 = ;

{ 1 30 dup to num | coprime } zet. {1,7,11,13,17,19,23,29} ok

17 num invmod . 23  ok

17 23 num u*mod . 1  ok

About adding all numbers in two sets:

: setint+ \ n -- | s -- s'

  0 >xst                       \ empty set on x

  foreach                      \ for each number in s

  ?do zst> over + >zst         \ element to datastack and back

     zfence xzmerge            \ merge {x+y} to the set on x

  loop drop xst zst setmove ;  \ drop n and move the set to z

: set+ \ s s' -- s"

  0 >xst 

  zst yst setmove

  foreach

  ?do zst>                     \ element to data stack

     yzcopy1 setint+           \ add element to all elements in s'

     xzmergered                \ merge this set to the set on x

  loop yst setdrop 

  xst zst setmove ;

{ 3 100 | prime } zdup set+ cr zet.

{6,8,10,12,14,16,18,20,22,24,26,28,30,32,34,36,38,40,42,44,46,48,50,52,54,56,58,60,62,64,66,68,70,72,74,76,78,80,82,84,86,88,90,92,94,96,98,100,102,104,106,108,110,112,114,116,118,120,122,124,126,128,130,132,134,136,138,140,142,144,146,148,150,152,154,156,158,160,162,164,166,168,170,172,176,178,180,186,194} ok

Hmm..! Goldbach seems to be right...

The first failure is 174.

{ 3 200 | prime } zdup set+ cr zet.

{6,8,10,12,14,16,18,20,22,24,26,28,30,32,34,36,38,40,42,44,46,48,50,52,54,56,58,60,62,64,66,68,70,72,74,76,78,80,82,84,86,88,90,92,94,96,98,100,102,104,106,108,110,112,114,116,118,120,122,124,126,128,130,132,134,136,138,140,142,144,146,148,150,152,154,156,158,160,162,164,166,168,170,172,174,176,178,180,182,184,186,188,190,192,194,196,198,200,202,204,206,208,210,212,214,216,218,220,222,224,226,228,230,232,234,236,238,240,242,244,246,248,250,252,254,256,258,260,262,264,266,268,270,272,274,276,278,280,282,284,286,288,290,292,294,296,298,300,302,304,306,308,310,312,314,316,318,320,322,324,326,328,330,332,334,336,338,340,342,344,346,348,350,352,354,356,358,360,362,364,366,370,372,374,376,378,380,382,384,386,388,390,392,394,396,398} ok

First failure at 368.

{ 3 300 | prime } zdup set+ cr zet.

{6,8,10,12,14,16,18,20,22,24,26,28,30,32,34,36,38,40,42,44,46,48,50,52,54,56,58,60,62,64,66,68,70,72,74,76,78,80,82,84,86,88,90,92,94,96,98,100,102,104,106,108,110,112,114,116,118,120,122,124,126,128,130,132,134,136,138,140,142,144,146,148,150,152,154,156,158,160,162,164,166,168,170,172,174,176,178,180,182,184,186,188,190,192,194,196,198,200,202,204,206,208,210,212,214,216,218,220,222,224,226,228,230,232,234,236,238,240,242,244,246,248,250,252,254,256,258,260,262,264,266,268,270,272,274,276,278,280,282,284,286,288,290,292,294,296,298,300,302,304,306,308,310,312,314,316,318,320,322,324,326,328,330,332,334,336,338,340,342,344,346,348,350,352,354,356,358,360,362,364,366,368,370,372,374,376,378,380,382,384,386,388,390,392,394,396,398,400,402,404,406,408,410,412,414,416,418,420,422,424,426,428,430,432,434,436,438,440,442,444,446,448,450,452,454,456,458,460,462,464,466,468,470,472,474,476,478,480,482,484,486,488,490,492,494,496,498,500,502,504,506,508,510,512,514,516,518,520,522,524,526,528,532,534,538,540,542,544,546,548,550,552,554,556,558,560,562,564,566,570,574,576,586} ok

Fail(300)=530. The first failure is a kind of measure of the probability of Goldbachs conjecture to be true.

Additional words for permutation groups:

\ Is s a permutation subgroup of s'?

: psub? \ -- flag | s s' --

  zover zswap subset 0= 

  if zdrop false exit 

  then permgroup? ;

\ Is s a normal permutation subgroup of s'?

: pnsub? \ -- flag | s s' --

  zover zover psub? 0= 

  if zdrop zdrop false exit then

  zswap zst yst setmove

  begin zst@

  while zsplit yzcopy1 zover prcoset 

     zswap yzcopy1 plcoset zet= 0=    \ false

     if zdrop yst setdrop false exit then

  repeat zdrop          \ dropping the empty set left in z

  yst setdrop true ;

\ s' is the set of normal subgroups of s

: pnsubgroups \ s -- s'

  zst yst setcopy

  psubgroups

  0 >xst

  begin zst@

  while zsplit zdup yzcopy1 pnsub?

     if zfence xzmerge else zdrop then

  repeat zdrop yst setdrop xst zst setmove ;

{ ( 4 1 2 3 ) ( 2 1 3 4 ) } generate zdup cardinality . 24  ok

pnsubgroups cr zet.

{{(2,1,3,4),(3,4,2,1),(3,1,4,2),(2,3,4,1),(1,4,2,3),(1,3,4,2),(1,2,3,4),(4,2,3,1),(4,1,2,3),(3,2,1,4),(4,1,3,2),(4,2,1,3),(4,3,1,2),(2,1,4,3),(2,3,1,4),(2,4,1,3),(3,1,2,4),(3,4,1,2),(3,2,4,1),(4,3,2,1),(2,4,3,1),(1,2,4,3),(1,3,2,4),(1,4,3,2)},{(4,3,2,1),(1,3,4,2),(1,2,3,4),(4,2,1,3),(3,1,2,4),(2,3,1,4),(2,4,3,1),(3,2,4,1),(1,4,2,3),(4,1,3,2),(2,1,4,3),(3,4,1,2)},{(4,3,2,1),(1,2,3,4),(3,4,1,2),(2,1,4,3)},{(1,2,3,4)}} ok

A simple implementation of permutation groups

I have big problems with groups. If you in GAP write 

gap> s7 := Group ( (1, 2, 3, 4, 5, 6, 7), (1, 2) );

GAP immediately respond 

Group([ (1,2,3,4,5,6,7), (1,2) ])

And if you then write

gap> Elements ( s7 );

GAP prints the list of all 7! elements of Sym(7) in a second or so.

In my simple implementation of permutation groups it takes about two seconds to generate Sym(5), and Sym(6) seems to be out of reach. The big discrepancy comes from that GAP has been developed at the universities for twenty years or so and that I so far has almost no knowledge of computational algebra. 

In GAP a group is determined by its generators while in ZET it's determined by the set of permutations. To me it's interesting to examine subgroups of Sym(5) in this form, but it would also be interesting to try to make a more general and perhaps more efficient implementation later on.

In GAP i.e. (1,2,3,4) denotes a cycle but here the same vector will denote the identity permutation and the cycle in GAP would correspond to (2,3,4,1) here.

\ The number of permutations in a set of permutations
: ord \ -- n | s -- s
  zst> zst> 2dup >zst >zst
  cs 1+ swap cs swap / ;

\ The number of elements to be permuted in v
: numb \ -- n | v --
  zst@ cs zdrop ;

\ j=v(i)
: pmaps \ i -- j | v --
  zdrop cells zst @ + @ ;

\ composition of permutations as functions
: permcomp \ v1 v2 -- v1v2
  ( zst@ cs 1+ 1
  do zover zover i pmaps pmaps
  loop ) znip znip ;

\ generation of cyclic permutation group
: pgen \ v -- s
  zst yst setcopy -1 1
  do zdup yzcopy1 permcomp zdup yzcopy1 vector=
     if numb 1+ i * 2* negate >zst leave then
  loop yst setdrop ;

\ right coset
: prcoset \ s v -- s'
  0 >xst
  zst yst setmove
  foreach
  ?do yzcopy1 permcomp zfence xzmerge
  loop yst setdrop xst zst setmove ;

\ left coset
: plcoset \ v s -- s'
  0 >xst
  zswap zst yst setmove
  foreach
  ?do yzcopy1 zswap permcomp zfence xzmerge
  loop yst setdrop xst zst setmove ;

\ componentwise composition of permutation sets
: permset* \ -- | s1 s2 -- s3
  0 >xst
  zst yst setmove
  foreach
  ?do yzcopy1 plcoset
  xzmergered
  loop yst setdrop
  xst zst setmove ;

: permgroup? \ -- flag | s --
  zdup zdup permset* zet= ;

\ Generation of standard permutations
: pidentity \ n -- | -- v
  >r ( r> 1+ 1 ?do i loop ) ;

: pcirc \ n -- | -- v
  >r ( r> dup 1 ?do i loop )  ;

: proll \ n -- | -- v
  >r ( r@ 1- dup 1 do i loop r> ) ;

\ The number of element to be permuted in permutations in s
: perm# \ -- n | s -- s
  zst> zst> tuck >zst >zst cs ;

\ Calculate the inverse permutation
: pinv \ v -- v'
  zdup adn2 drop adn1 -rot loc{ a2 a1 } cell/ 1
  do i dup 1- cells a2 + @ 1- cells a1 + ! loop znip ;

\ add the inverses to all permutations in s
: adinv \ s -- s'
  0 >xst zdup xzmerge foreach
  do pinv zfence xzmerge
  loop xst zst setmove reduce ;

\ generates the group s' from the generators in s
: generate \ s -- s'
  zst yst setcopy 0 >xst foreach
  ?do pgen xzmerge
  loop xst zst setmove reduce 1
  begin yzcopy1 zswap permset*
     yzcopy1 permset* ord tuck =
  until yst setdrop drop ;

\ generate set of groups s' from set of generators s
: multigen \ s -- s'
  0 >xst foreach
  ?do generate zfence xzmerge
  loop xst zst setmove reduce ;

\ Set of all subgroups to s
: psubgroups \ s -- s'
  perm# pidentity zfence zfence
  zst yst setmove foreach
  do yst zst setmove zdup zrot multincl
     multigen union zst yst setmove
  loop yst zst setmove ;

{ ( 4 1 2 3 ) ( 2 1 3 4 ) } generate  ok

ord . 24  ok

psubgroups  ok

zdup  ok

cardinality . 30  ok

zet. {{(4,3,2,1),(1,2,3,4)},{(1,2,3,4)},{(1,3,4,2),(1,4,2,3),(1,2,3,4)},{(2,1,4,3),(1,2,3,4)},{(4,1,3,2),(2,4,3,1),(1,2,3,4)},{(3,4,1,2),(1,2,3,4)},{(4,3,2,1),(1,2,3,4),(3,4,1,2),(2,1,4,3)},{(3,1,2,4),(2,3,1,4),(1,2,3,4)},{(1,3,2,4),(1,2,3,4)},{(3,2,1,4),(1,2,3,4)},{(1,2,4,3),(1,2,3,4)},{(4,2,3,1),(1,2,3,4)},{(4,3,2,1),(1,2,3,4),(4,2,3,1),(1,3,2,4)},{(3,4,1,2),(1,2,3,4),(1,4,3,2),(3,2,1,4)},{(1,3,4,2),(1,4,2,3),(1,2,3,4),(1,4,3,2),(1,2,4,3),(1,3,2,4)},{(1,4,3,2),(1,2,3,4)},{(4,3,1,2),(2,1,4,3),(3,4,2,1),(1,2,3,4)},{(1,2,4,3),(1,2,3,4),(3,2,1,4),(4,2,1,3),(3,2,4,1),(4,2,3,1)},{(4,3,2,1),(1,3,4,2),(1,2,3,4),(4,2,1,3),(3,1,2,4),(2,3,1,4),(2,4,3,1),(3,2,4,1),(1,4,2,3),(4,1,3,2),(2,1,4,3),(3,4,1,2)},{(3,2,4,1),(4,2,1,3),(1,2,3,4)},{(3,1,4,2),(3,4,1,2),(1,3,2,4),(1,2,3,4),(2,1,4,3),(4,3,2,1),(2,4,1,3),(4,2,3,1)},{(2,4,1,3),(3,1,4,2),(4,3,2,1),(1,2,3,4)},{(3,1,2,4),(2,3,1,4),(1,2,3,4),(2,1,3,4),(1,3,2,4),(3,2,1,4)},{(4,1,3,2),(2,4,3,1),(1,2,3,4),(2,1,3,4),(1,4,3,2),(4,2,3,1)},{(2,1,4,3),(1,2,3,4),(2,1,3,4),(1,2,4,3)},{(2,1,3,4),(1,2,3,4)},{(4,3,2,1),(1,2,3,4),(2,1,3,4),(3,4,2,1),(1,2,4,3),(4,3,1,2),(3,4,1,2),(2,1,4,3)},{(2,1,3,4),(3,4,2,1),(3,1,4,2),(2,3,4,1),(1,4,2,3),(1,3,4,2),(1,2,3,4),(4,2,3,1),(4,1,2,3),(3,2,1,4),(4,1,3,2),(4,2,1,3),(4,3,1,2),(2,1,4,3),(2,3,1,4),(2,4,1,3),(3,1,2,4),(3,4,1,2),(3,2,4,1),(4,3,2,1),(2,4,3,1),(1,2,4,3),(1,3,2,4),(1,4,3,2)},{(4,1,2,3),(3,4,1,2),(2,3,4,1),(1,2,3,4)},{(4,1,2,3),(2,1,4,3),(3,2,1,4),(1,2,3,4),(4,3,2,1),(3,4,1,2),(2,3,4,1),(1,4,3,2)}} ok

My purpose with the nested sets was to investigate the possibility of simple dynamic data structures without the need of garbage collection. The stacks are administrated as usual and rest data on the stacks comes from faulty programming. So in a way the word q, that resets the stacks in case of error (manually now, but should be automatic), and the word drop, replace the garbage collection systems used in traditional programming with dynamic data. Except from with the present primitive implementation of groups it turned out to be surprisingly efficient.