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.

Directed graphs

A relation (V,V,E) is equivalent with a directed graph (V,E), where E⊆V×V is the set of edges and V is the set of vertices.

There is a path from the node x to the node y if there exist a sequence (x,x1),(x1,x2),...,(xn,y)∈E:

: path? \ x y -- flag | E --
  swap >zst zfence
  begin zover zover ztuck subimage      \ E s s s'
     union zdup dup smember             \ E s s" f
     if drop zdrop zdrop zdrop true exit
     then zswap zover zet=              \ E s" s=s"
     if drop zdrop zdrop false exit
     then
  again ;

To make a pair on zst-stack from two integers of the data stack:

: ipair \ m n -- | -- (m,n)
  2>r ( 2r> ) ;

Find the set s of all nodes in (V,E) without incoming arrows:

: sourceset \ (V,E) -- s
  unfence image diff ;

Find the set s of all nodes in (V,E) without outgoing arrows:

: sinkset \ (V,E) -- s
  unfence coimage diff ;

The next word merge the set on top of zst into the set on top of the xst leaving the result in xst.

: xzmerge \ s --
  xst zst setmove
  zswap zetmerge                   \ swap to build from right
  zst xst setmove ;

Given a directed graph (V,E), a topological sort is an ordered list including all the elements in V once, (x1,...xn), sorted so that if there is a path from xi to xj, then i<j. Such paths exist if and only if there are no directed cycles.

Kahn's algorithm, from Wikipedia.

: toposort \ (V,E) -- s
  0 >xst                           \ empty set in x
  zdup sourceset zst yst setmove   \ source nodes in y
  unfence znip                     \ drop V keep E
  begin yst@                       \ while source nodes left
  while ysplit yst> dup            \ remove node m
     zdup >zst zfence zdup xzmerge \ add m to the set in x
     subimage                      \ set of all n: m→n
     begin zst@                    \ while that set non empty
     while zsplit zst> zswap       \ remove node n, E tos
        2dup ipair zfence diff     \ E:=E\{(m,n)}
        dup zdup >zst zfence       \ build set of all nodes..
        subcoimage zst@ 0=         \ ..pointing at n
        if >yst yfence ymerge      \ add n to y-set if empty
        else drop                  \ else drop n
        then zdrop zswap           \ drop set, swap E back
     repeat zdrop drop             \ drop empty set and node m
  repeat yst> drop zst@ zdrop      \ drop empty set and E
  if xst setdrop 0 >zst            \ if |E|>0 flag with empty set
  else xst zst setmove             \ else move the x-set to zst
     zst> 1- >zst                  \ mark it as an ordered list
  then ;

A directed graph without directed cycles is called a dag (directed acyclic graph):

: dag? \ -- f | (V,E) --
  toposort zst@ 0= 0= zdrop ;

A loop is an edge from a node to itself and loopset gives the edge set to a graph in which all edges are loops. Can be used to filtrate the loops from digraphs.

: loopset \ V -- E
  { foreach ?do ( zst> dup ) loop } ;

Generate a random pair of nodes:

: randpair \ |V| -- | -- (m,n)
  dup random 1+ swap random 1+ ipair ;

Generate a random digraph with certain number of vertices and edges.

: rand-digraph \ |V| |E| -- | -- (V,E)
  { over 1+ 1 ?do i loop }
  0 >zst
  begin over rand-pair zfence union zdup cardinality over =
  until 2drop
  pair ;

: rand-noloop-digraph \ |V| |E| -- | -- (V,E)
  { over 1+ 1 ?do i loop }
  0 >zst
  begin over rand-pair zfence union
     zover loopset diff
     zdup cardinality over =
  until 2drop pair ;

: rand-acyclic-digraph \ m n -- | -- (V,E)
  begin 2dup rand-noloop-digraph zdup dag? 0=
  while zdrop
  repeat 2drop ;

10 20 rand-noloop-digraph  ok
zdup cr zet.
({1,2,3,4,5,6,7,8,9,10},{(5,2),(2,5),(10,3),(9,3),(9,10),(3,7),(4,10),(8,4),(3,8),(5,4),(3,6),(2,8),(6,10),(10,6),(7,6),(10,2),(7,4),(6,2),(3,9),(3,1)}) ok
toposort zet. 0  ok

Already the two first edges builds a directed loop.

10 20 rand-acyclic-digraph zdup cr zet.
({1,2,3,4,5,6,7,8,9,10},{(4,7),(1,10),(3,8),(4,8),(5,6),(8,2),(3,7),(2,10),(4,10),(2,6),(1,2),(9,2),(3,4),(7,2),(1,3),(1,7),(4,6),(9,4),(3,5),(5,4)}) ok
toposort zet. (9,1,3,5,4,7,8,2,6,10) ok

Relations

A relation is a triplet (A,B,R) where A and B are sets and R⊆A×B. A function from A to B is a relation (A,B,f) where

(x1,y),(x2,y)∈f ⇒ x1=x2 and x∈A ⇒ ∃y∈B: (x,y)∈f.

Some definitions:

\ remove brackets of object at top of stack

: unfence zst> drop ;

: domain \ (A,B,R) -- A
  unfence zdrop zdrop ;

: codomain \ (A,B,R) -- B
  unfence zdrop znip ;

: rel \ (A,B,R) -- R
  unfence znip znip ;

The set of all elements in the codomain that is related to some element in the domain:

\ y∈image(R) ⇔ ∃x:(x,y)∈R

: image \ R -- s
  { foreach ?do unfence zst> zst> drop loop } ;

The set of all elements in the domain that is related to some element in the codomain:

\ x∈coimage(R) ⇔ ∃y:(x,y)∈R

: coimage \ R -- s
  { foreach ?do unfence zst> drop zst> loop } ;

The image of a relation restricted to a subset s of the domain:

: subimage \ R s -- s'
  zst yst setmove
  { foreach
  ?do unfence zst> zst> yst zst setcopy smember 0=
     if drop then
  loop } yst setdrop ;

The coimage of a relation restricted to a subset s of the codomain:

: subcoimage \ R s -- s'
  zst yst setmove
  { foreach
  ?do unfence zst> zst> yst zst setcopy swap smember 0=
     if drop then
  loop } yst setdrop ;

Test if a relation (A,B,R) is a function:

: func? \ -- flag | (A,B,R) --
  unfence znip
  zst yst setmove true
  begin zst@
  while zsplit zst> yst zst setcopy >zst zfence
     subimage cardinality 1 = 0=
     if 0= zdrop yst setdrop exit then
  repeat zdrop yst setdrop ;

Evaluate f(x):

: eval \ x -- y | f --
  >zst zfence subimage unfence zst> ;

Making a ordered pair or triplet of the top bundles:

: pair \ s1 s2 -- (s1,s2)
  zswap zst@ 2 - zswap zst@ 2 - + 1- >zst ;

: triplet \ s1 s2 s3 -- (s1,s2,s3)
  zrot zst@ 2 - zrot zst@ 2 - zrot zst@ 2 - + + 1- >zst ;

The composition of two relations (A,B,R) and (B,C,S) is the relation (A,C,SR) defined by

(a,c)∈SR ⇔ ∃b∈B:(a,b)∈R & (b,c)∈S.

: composition \ (A,B,R) (B,C,S) -- (A,C,SR)
  0 >xst                                \ empty set on xst-stack
  unfence zrot zdrop zrot unfence       \ C S A B R
  zst yst setmove zdrop zswap           \ C A S
  zst yst setmove                       \ R S in yst
  zswap zover zover cartprod            \ A C A×C
  begin zst@                            \ while elements in top set
  while zsplit infence
     yzcopy1 zover zsplit znip subcoimage
     zst xst setmove
     yzcopy2 zover zsplit zdrop unfence subimage
     xst zst setmove intersection zst@ zdrop
     if unfence unfence zst> unfence >zst -5 >zst zfence
        xst zst setmove zetmerge zst xst setmove
     else zdrop
     then
  repeat zdrop yst setdrop yst setdrop
  xst zst setmove triplet ;

Set algebra

In bundles that represent sets all integers must be sorted, because the word smember stop searching for a single member when reaching an integer less than the integer to be tested. While zetmerge is faster than union it's safer to use union. Never the less, I have tried to replace union with zetmerge as much as possible, with considerable faster code. It seems to work.

: union \ -- | s s' -- sUs'
  zetmerge set-sort reduce ;

: intersection \ -- | s s' -- sΛs'
  0 >xst zst yst setmove
  begin zst@
  while zsplit zfence zdup zst> drop
     yst zst setcopy member
     if xst zst setmove zetmerge zst xst setmove
     else zdrop
     then
  repeat zdrop yst setdrop
  xst zst setmove reduce ;

: diff \ -- | s s' -- s\s'
  0 >xst zst yst setmove
  begin zst@
  while zsplit zfence zdup zst> drop
     yst zst setcopy member
     if zdrop
     else xst zst setmove zetmerge zst xst setmove
     then
  repeat zdrop yst setdrop
  xst zst setmove reduce ;

: multincl \ -- |{s1,...,sn} x -- {s1U{x},...,snU{x}}
  0 >xst zfence zst yst setmove
  begin zst@
  while zsplit yst zst setcopy union zfence
     xst zst setmove zetmerge zst xst setmove
  repeat zdrop yst setdrop xst zst setmove ;

: powerset \ -- | s -- p(s)      Set of all subsets
  zst@ 0= if -2 >zst exit then
  zsplit zfence zst yst setmove recurse
  zdup yst zst setmove zst> drop multincl
  zetmerge ;

: cartprod \ -- | s s' -- s×s'   Cartesian product
  zst yst setmove
  zst xst setcopy xst> drop cardinality 0 0 >zst
  ?do xfence -1 xst+!
     yst setdup
     begin yst@
     while ysplit yfence -1 yst+!
        xst zst setcopy
        yst zst setmove vmerge
        zfence
        zetmerge
     repeat yst> drop xst setdrop
  loop yst setdrop ;

: infence \ -- |{x1,...,xn} -- {{x1},...,{xn}}
  0 >xst foreach
  ?do zfence zfence
     xst zst setmove zetmerge
     zst xst setmove
  loop xst zst setmove ;

For finite sets there is a recursive method to compute the set of all subsets of a certain order which is an analogy to Pascals rule for binomial coefficients:

Suppose A is a set (with n elements for the analogy), that f is a choice function on all finite sets and that p(A,k) is the set of all subsets of A with k elements. For a set of sets S, define S%x to be the set consisting of all sets X in S where x have been included. Then

p(A,k)=p(A\{f(A)},k)+(p(A\{f(A)},k-1)%f(A))

where + stands for union of two disjoint sets: non of the sets in the first set include f(A), which every set in the second set do.

The choice function here is the topmost element on the stack where the set is pushed. This will work because the order of the sets in the computation not will be changed.

\ Set of all subsets with k elements

\ p(A,k)=p(A\{f(A)},k)+(p(A\{f(A)},k-1)%f(A))
: power# \ k -- | s -- p(s,k)
  ?dup 0= if zdrop 0 >zst zfence exit then
  dup 1 = if drop infence exit then
  dup zdup cardinality =
  if drop zfence exit then
  dup 1 = if drop infence exit then
  zsplit zfence zst xst setmove
  dup zdup recurse
  zswap 1- recurse xst zst setmove
  zst> drop multincl
  zetmerge ;

Before trying to calculate p(A,k) it's a good idea to count the number of elements:

|p(A,k)|=choose(|A|,k).

\ http://rosettacode.org/wiki/Evaluate_binomial_coefficients#Forth
: choose \ n k -- nCk
  1 swap 0 ?do over i - i 1+ */ loop nip ;

: multiunion \ -- |{s1,...,sn} -- s1U...Usn
  foreach 0 >zst
  ?do zetmerge
  loop set-sort reduce ;

: zetcup \ -- |{s1,...,sn} s -- {s1Us,...,snUs}

  zst xst setmove 0 >yst foreach
  ?do xst zst setcopy union zfence
     yst zst setmove zetmerge zst yst setmove
  loop xst setdrop yst zst setmove ;

: zetcap \ -- |{s1,...,sn} s' -- {s1Λs',...,snΛs'}

  zst xst setmove 0 >yst foreach
  ?do xst zst setcopy intersection zfence
     yst zst setmove zetmerge zst yst setmove
  loop xst setdrop yst zst setmove ;

: zetunion \ -- |{ s1,...,sn} {t1,...,tm} -- {siUtj}ij
  0 >xst zst yst setmove foreach
  ?do yst zst setcopy
     zswap zetcup
     xst zst setmove union
     zst xst setmove
  loop yst setdrop xst zst setmove ;

The set of functions:

: functions \ -- | s s' -- fun(s,s')

  secobjad @ 0= if zdrop -2 >zst exit then
  secobjad @ -2 = if cartprod infence exit then
  zswap zsplit zfence zst xst setmove
  zover recurse zswap xst zst setmove
  zswap cartprod infence zetunion ;

utime { 1 2 3 } zdup functions 3 power# cardinality cr . utime d- d.
2925 -628052  ok

ZET creates and count this set with 2925 elements in 0.63 seconds.

utime { 1 100000 | prime } utime cr d- d. cardinality .
-6652182 9592  ok
utime { 1 10 | all } powerset utime cardinality cr . d- d.
512 -24502  ok
utime { 1 11 | all } powerset utime cardinality cr . d- d.
1024 -83816  ok
utime { 1 12 | all } powerset utime cardinality cr . d- d.
2048 -251767  ok
utime { 1 13 | all } powerset utime cardinality cr . d- d.
4096 -986219  ok
utime { 1 14 | all } powerset utime cardinality cr . d- d.
8192 -3944057  ok

{ 1 2 3 } zdup functions cr zet.

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

Nested sets

[Some errors are corrected and some other changes is done]

The previous posts about sets are obsolete from now and are removed from the code site, namely: Dynamic sets, Permutation groups, Subsets and subgroups, Moore about subgroups, A conjecture about groups, Topology.

The implementation of sets is similar to the previous, but now the sets are handled by the three implemented stacks xst, yst and zst, where zst is the parameterstack for sets. An other difference is that the bundles are coded to permit vectors as elements.

{ 0 1 2 3 ( 4 5 6 6 ) { 7 { 8 8 } } 9 0 } cr showz cr zet.
1 2 3 4 5 6 6 -9 7 8 -2 -6 9 0 -28
{1,2,3,(4,5,6,6),{7,{8}},9,0} ok

Negative integers indicate a bundle count, even for sets and odd for vectors. If the absolute value of these numbers are divided by 2 the number of integers in the bundle is obtained.

Stacks

Implementation of the three stacks:

: cs negate 2/ ;
: listflag 1 and ;

: objsize \ bc -- n
  dup 0< if cs 1+ else drop 1 then ;

cell negate constant -cell
: >stack ( n ad -- )  cell over +! @ ! ;
: stack> ( ad -- n )  dup @ @ -cell rot +! ;
: >stack> ( n ad -- m )  dup @ @ -rot @ ! ;
: stack@ ( ad -- n )  @ @ ;
: stack! ( n ad -- )  @ ! ;
: stack+! ( n ad -- )  @ +! ;

cell 1- log~ constant cellshift
: stack-depth ( ad -- n )  dup @ swap - cellshift rshift ;
: stack-cl ( ad -- )  dup ! ;
: stack-empty ( ad -- flag )  dup @ = ;

1 16 lshift cells allocate throw dup constant xst dup !
: >xst ( n -- )  xst >stack ;
: xst> ( -- n )  xst stack> ;
: >xst> ( n -- m )  xst >stack> ;
: xst@ ( -- n )  xst @ @ ;
: xst! ( n -- )  xst @ ! ;
: xst+! ( n -- )  xst @ +! ;
: >>xst ( xn ... x1 bc -- )  >r r@ cs 0 ?do >xst loop r> >xst ;
: xst>> ( -- x1 ... xn bc )  xst@ >r xst> cs 0 ?do xst> loop r> ;

1 20 lshift cells allocate throw dup constant yst dup !
: >yst ( n -- )  yst >stack ;
: yst> ( -- n )  yst stack> ;
: >yst> ( n -- m )  yst >stack> ;
: yst@ ( -- n )  yst @ @ ;
: yst! ( n -- )  yst @ ! ;
: yst+! ( n -- )  yst @ +! ;
: >>yst ( xn ... x1 bc -- )  >r r@ cs 0 ?do >yst loop r> >yst ;
: yst>> ( -- x1 ... xn bc )  yst@ >r yst> cs 0 ?do yst> loop r> ;

1 21 lshift cells allocate throw dup constant zst dup !
: >zst ( n -- )  zst >stack ;
: zst> ( -- n )  zst stack> ;
: >zst> ( n -- m )  zst >stack> ;
: zst@ ( -- n )  zst @ @ ;
: zst! ( n -- )  zst @ ! ;
: zst+! ( n -- )  zst @ +! ;
: >>zst ( xn ... x1 bc -- )  >r r@ cs 0 ?do >zst loop r> >zst ;
: zst>> ( -- x1 ... xn -n )  zst@ >r zst> cs 0 ?do zst> loop r> ;

: showx xst stack-depth if xst> >r recurse r> dup . >xst then ;
: showy yst stack-depth if yst> >r recurse r> dup . >yst then ;
: showz zst stack-depth if zst> >r recurse r> dup . >zst then ;

: >zet ( s -- | -- s)
  >>yst yst> dup >r cs 0
  ?do yst> >zst loop r> >zst ;

: zet> ( -- s | s -- )
  zst> dup >r cs 0
  ?do zst> >xst loop r> >xst xst>> ;

The word >>xst moves a bundle from the datastack to xst and xst>> do the reverse. But the order of the elements will be reversed. The words >zet and zet> from moves bundles between the datastack and the zst-stack without reversing the order.

All words beginning with set is called with one or two stack addresses.

: setdrop \ ad --
  dup @ @ cs cells cell+ negate swap +! ;

: setdup \ ad --
  >r
  r@ @ @ cs cells                 \ n'
  r@ @ over -                     \ n' ad1
  r@ @ cell+                      \ n' ad1 ad2
  rot cell+ dup r> +! cmove ;

: setover \ ad --
  dup >r @ @ cs cells cell+       \ nr of bytes 1'st set
  r@ @ swap -                     \ ad to 2'nd set
  dup @ cs cells cell+ dup >r -   \ ad to 3'rd set
  cell+ r> r@ @ cell+             \ ad to move to
  swap dup r> +! cmove ;

: setcopy loc{ ad1 ad2 -- }
  ad1 @ @ cs cells             \ n'
  ad1 @ over - swap cell+      \ ad1-n' n
  ad2 @ cell+ over ad2 +! swap cmove ;

: setmove \ ad1 ad2 --
  swap dup rot setcopy setdrop ;

The three words below gets the addresses and the counts for the first, second and third sets on the zst-stack.

: adn1 zst@ cs cells zst @ over - swap cell+ ;
: adn2 adn1 drop cell- dup @ cs cells tuck - swap cell+ ;
: adn3 adn2 drop cell- dup @ cs cells tuck - swap cell+ ;

All words beginning with z acts on the zst-stack, and the words below manipulates sets.

: zdup  zst setdup ;
: zdrop  zst setdrop ;
: zover  adn2 tuck zst @ cell+ swap cmove zst +! ;
: zswap  zover adn2 adn3 rot + move zdrop ;
: znip  zswap zdrop ;
: ztuck  zswap zover ;
: zrot  zst>> zswap >>zst zswap ;

Output of sets

The output is built up backwards in a buffer which is printed out.

0 value addr1

: addr1- \ --
  addr1 1- to addr1 ;

: fillad$ \ addr n --
  dup 1- negate addr1 + dup to addr1 swap move addr1- ;

: n>addr1 \ n --
  0 <# #s #> fillad$ ;

: a>addr1 \ c --
  addr1 c! addr1- ;

Defining the cardinality here for the definition of foreach that removes the fences around the top set on the zst-stack and prepare for a do-loop.

: cardinality \ -- n | s --
  zst> cs dup >xst 0
  ?do zst@ 0<
     if zst@ dup cs negate xst+! >r zdrop r> cs 1+
     else zst> drop 1
     then
  +loop xst> ;

: foreach \ -- n 0 | s -- z1...zn
  zdup cardinality zst> drop 0 ;

: closep \ -- bc asc
  zst@ dup listflag if [char] ) else [char] } then ;

: openp \ bc -- asc
  listflag if [char] ( else [char] { then ;

: list$ \ n1...nk -k ad -- ad n
  dup to addr1 false loc{ addr2 flag }
  closep a>addr1
  foreach
  do flag if [char] , a>addr1 then zst@ 0<
     if addr1 recurse 2drop
     else zst> n>addr1
     then flag 0= if true to flag then
  loop openp a>addr1
  addr1 1+ addr2 over - 1+ ;

\ Corrected code:
1 20 lshift dup allocate throw swap cell - + constant printbuf

: zet. \ -- | s --   prints top set on zst stack
  zst@ 0=
  if zst> .
  else printbuf list$ type
  then ;

: set. \ ad --  prints top set on xst or yst stack
  zst setcopy zet. ;

Analyzing sets

The next word analyse a bundle cell: 0 integer, 1 vector, 2 set.

: ?obj \ x -- 0,1,2
  dup 0<
  if listflag
     if 1 else 2 then
  else drop 0
  then ;

Splitting a set is to put the top element at top of stack just before the rest of the set, which eventually is the empty set 0.

: _split \ ad --   ad=yst,zst
  dup >r @ cell- @ 0< 0=
  if r@ stack> 2 + r@ stack> swap r@ >stack r> >stack exit then
  r@ stack>
  r@ xst setmove
  xst@ cs 1+ 2* + r@ >stack
  xst r> setmove ;

: ysplit \ -- | s -- s' x  in yst stack
  yst _split ;

: zsplit \ -- | s -- s' x
  zst _split ;

Set equal, subset and membership

The word zet= is defined by subset, member is defined by zet= and subset is defined by member, by recursion. Next word examines if the integer n is a member of the set s.

A change is made in this section. All sets of integers are sorted and smember use that for faster exit when 'not member'.

: zetmerge \ -- | s s' -- s"
  zst yst setmove
  yst@ zst> +
  yst zst setmove
  zst! ;

: vmerge \ -- | v v'-- v"
  zst yst setmove
  yst@ zst> + 1+
  yst zst setmove
  zst! ;

: _fence \ ad -- | x -- {x}
  dup >r stack@ ?obj
  case 0 of -2 r@ >stack endof
       1 of r@ stack@ 1- r@ >stack endof
       2 of r@ stack@ 2 - r@ >stack endof
  endcase rdrop ;

: xfence xst _fence ;
: yfence yst _fence ;
: zfence zst _fence ;
 
: set-sort \ -- | s -- n1...nk -2k
  0 loc{ counter } 0 >xst 0 >yst
  foreach
  ?do zst@ ?obj
     case 0 of counter 1+ to counter zst> endof
          1 of zfence xst zst setmove zetmerge zst xst setmove endof
          2 of zfence yst zst setmove zetmerge zst yst setmove endof
     endcase
  loop counter sort 2* negate >zet
  xst zst setmove zetmerge
  yst zst setmove zetmerge ;
 
: smember \ n -- flag | s --
  zst@ cs false loc{ m flag }
  foreach
  ?do zst@ 0<
     if m zst@ cs 1+ - to m zdrop
     else m 1- to m dup zst> 2dup >
        if false to flag 2drop
           m cells negate zst +! leave
        then =
        if true to flag
           m cells negate zst +! leave
        then
     then
  loop drop flag ;

Equality for vectors:

: vect= \ s -- flag | s' --
\ non empty list not including non empty sets
  dup zst@ = 0=
  if zdrop cs 0 ?do drop loop false exit
  then true loc{ flag } zst> drop cs 0
  ?do flag
     if zst> = 0= if false to flag then
     else zst> 2drop
     then
  loop flag ;

: vector= \ -- flag | s s' --
  zet> vect= ;

Examines if the vector s is a member in the set s':

: vmember \ -- flag | s s' --
  zswap zst yst setmove
  zst@ cs false loc{ m flag }
  foreach
  ?do zst@ ?obj
    case 0 of m 1 - to m zst> drop endof
         1 of m zst@ cs 1+ - to m
              yst zst setcopy vector=
              if true to flag
                 m cells negate zst +! leave
              then endof
         2 of m zst@ cs 1+ - to m
              zst@ cs 1+ cells negate zst +! endof
    endcase
  loop yst setdrop flag ;

Get the count/integer of the second object of the zst-stack:

: secobjad \ -- ad | x y -- x y
  zst @ zst@ 0< if zst@ cs 1+ cells - else cell - then ;

Move the second object of zst-stack to datastack:

: routout \ -- x | x s -- s
  secobjad dup @ swap dup cell+ swap zst@ cs 1+ cells move
  zst> drop ;

0 value 'subset 
: subset \ -- flag | s s' --
  'subset execute ;

: zet= \ -- flag | s s' --
  zover zover subset
  if zswap subset
  else zdrop zdrop false
  then ;

Examines if s is a set-member:

: zet-member \ -- flag | s s' --
  zswap zst yst setmove
  begin zst@                         \ set not empty?
  while zsplit zst@ ?obj 2 =         \ element is a set?
     if yst zst setcopy zet= 
        if yst setdrop zdrop true exit then
     else zst@ ?obj if zdrop else zst> drop then
     then
  repeat yst setdrop zdrop false ;

: member \ -- flag | x s --
  secobjad @ ?obj
  case 0 of routout smember endof
       1 of vmember endof
       2 of zet-member endof
  endcase ;

:noname \ -- flag | s s' --          \ the subset code
  zst @ cell - 2@ or 0=
  if zdrop zdrop true exit then      \ true if both sets are empty
  zswap zst yst setmove
  begin yst@                         \ set is not empty?
  while ysplit yst@ ?obj
     if yst zst setmove zover member
     else yst> zdup smember
     then 0= if yst setdrop zdrop false exit then
  repeat yst> drop zdrop true ; to 'subset

Merge two sets on zst-stack:

: zetmerge \ -- | s s' -- s"
  zst yst setmove
  yst@ zst> +
  yst zst setmove
  zst! ;

Merge two vectors on zst-stack:

: vmerge \ -- | v v'-- v"
  zst yst setmove
  yst@ zst> + 1+
  yst zst setmove
  zst! ;

: _fence \ ad -- | x -- {x}
  dup >r stack@ ?obj
  case 0 of -2 r@ >stack endof
       1 of r@ stack@ 1- r@ >stack endof
       2 of r@ stack@ 2 - r@ >stack endof
  endcase rdrop ;

: xfence xst _fence ;
: yfence yst _fence ;
: zfence zst _fence ;

The important word that reduce multiple members in a set at top of zst-stack:

: reduce \ -- | s -- s'
  0 >yst foreach
  ?do zfence zdup zst> drop
     yst zst setcopy member
     if zdrop
     else yst zst setmove
        zetmerge zst yst setmove
     then
  loop yst zst setmove ;

Input of sets

0 create match ,
true value sort?

: { \ --
  1 match +! depth >xst true to sort? ;

\ Integer sorting is included
: } \ x1...xk --
  depth xst> - 2* negate
  -1 match +! >zet sort?
  if set-sort then reduce match @
  if zet> then true to sort? ;

Next word resets everything and should be automatic on error somehow.

: q  xst stack-cl yst stack-cl zst stack-cl 0 match ! abort ;

: (  { ;

: ) \ x1...xk --
  depth xst> - 2* 1+ negate
  -1 match +! >zet match @ if zet> then ;

Integer conditions

\ n -- flag
: pairprime dup prime over 2 + prime rot 2 - prime or and ; 
: odd 1 and ;  \ n -- flag
: 1mod4 4 mod 1 = ;  \ n -- flag
: 3mod4 4 mod 3 = ;  \ n -- flag
: sqr dup sqrtf dup * = ;
: all dup = ;
: sqrfree dup radical = ;     \ square free test
: semiprime bigomega 2 = ;    \ number is product of two primes?
: uniprime smallomega 1 = ;   \ number is power of single prime?
: biprime smallomega 2 = ;    \ number has two different pfactors?

: 2sqrsum dup 0               \ number sum of two squares?
  ?do dup i dup * - dup
     0< if drop false leave then
     sqr if true leave then
  loop nip ;

: | \ m n -- x1...xk
  swap ' loc{ xt }
  ?do i xt execute if i then
  loop false to sort? ;

{ 10000 20000 | pairprime } cardinality . 274  ok
53 >zst { 1 100 | prime } member . -1 ok

{ 100 200 | pairprime } { 100 200 | prime } subset . -1  ok

{ 1000 2000 | sqrfree } cardinality . 607  ok
{ 2000 3000 | sqrfree } cardinality . 609  ok
{ 3000 4000 | sqrfree } cardinality . 609  ok
{ 8000 9000 | sqrfree } cardinality . 608  ok

( 1 1 ) { { 0 } ( 1 1 ) } member . -1  ok

{ 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

{ 1 100 | uniprime } { 1 100 | semiprime } union cr zet.
{2,3,4,5,6,7,8,9,10,11,13,14,15,16,17,19,21,22,23,25,26,27,29,31,32,33,34,35,37,38,39,41,43,46,47,49,51,53,55,57,58,59,61,62,64,65,67,69,71,73,74,77,79,81,82,83,85,86,87,89,91,93,94,95,97} ok

Topology

[Attention, this code is not included in the code site!]

The blog has become a little messy and I will try to clean it up. These simple implementations of sets and groups are limited but interesting. I believe in the idea of bundles on a stack, but the data stack is not perfect for the purpose: for the sake of portability the data stack can only be manipulated with a restricted number of words. With a special stack for bundles, stack manipulations can be made much more efficient.

While going on with the current implementation for a while, I will try to find out a faster implementation. With better algorithms: until now I have used sloppy straight forward brutal force algorithms or worse! Also the implementation of permutations and groups must be enhanced.

Topological spaces

Topology is very easy to understand in spite of it's often formalistic presentation for the students. In set theory the question is: does an element belong to the set or not. In topology the question is if the element is close to the set or not.

That's what a topological structure λ on a set X submit, the possibility to decide whether a point x in X belongs to the closure of a subset A of X or not. The set λ of open subsets of X brings the method: the point x is in the closure of A if and only if every open set O in λ that contain x, also contain a point in A.

Two subsets of X are obligatory in λ and that's ø and X. Beyond that there are two rules for λ to be a topology on X. (1) Given a subset B of λ, then the union of all sets in B should belong to λ. (2) given two sets O,Ô in λ, then the intersection of O and Ô should belong to λ.

Examples:

• The topology of the real numbers (a topology which is a necessary part of the definition of the set of real numbers) consists of all unions of limited open intervalls, that is, open intervalls of the type {x|a<x<b} for real numbers a,b. Any union of such intervalls belongs to this topology. The point 1 belongs to the closure of {x|0<x<1}, without actually being a member of that set, because every open set that contains 1 also contains a point in {x|0<x<1}.
• The set {ø,X} is a topology on X called the trivial topology and the set of all subsets of X is a topology called the discrete topology.
• Any infinite set X has an obvious topology consisting of all complements of finite subsets of X plus the empty set.
• Any set of subsets of a set X generates a minimal topology on X which contain those subsets.

A subset C of X is said to be closed if X\C is open, that is belong to λ. The closure of a subset A of X is the smallest closed set C such that A is a subset of C. The most interesting infinite topological spaces are the Hausdorff spaces, where all singleton sets {x} are closed. But for finite topological spaces all Hausdorff spaces are discrete. Finite spaces are determined by the closure function of all the singleton sets, because with that function it's possible to decide which points that are close to which sets.

In stack diagram for setcup s is a set of sets and s' is a set of points. And s" is the set of the unions of all elements of s ans s':

{ { 1 2 } { 2 3 4 } } { 3 4 5 } setcup cr set.
{{1,2,3,4,5},{2,3,4,5}} ok

: setcup ( s s' -- s")
  >>xst 0 >yst foreach
  ?do nxst@ union fence yst>> union >>yst
  loop nxstdrop yst>> ;


Similar for the word:

: setcap ( s s' -- s")
  >>xst 0 >yst foreach
  ?do nxst@ intersection fence yst>> union >>yst
  loop nxstdrop yst>> ;

The word capgen generates all intersections:

: capgen ( s -- s')
  ndup >>yst ndup >>xst foreach
  ?do nyst@ nswap setcap xst>> union >>xst
  loop nystdrop xst>> ;

The word cupgen generates the set of all unions os s. When first apply capgen and then cupgen all the non trivial open sets in the smallest topology that includes s are generated. And the word topology completes with ø and X, which might or might not be generated by capgen and cupgen.

: cupgen ( s -- s' )
  ndup >>yst ndup >>xst foreach
  ?do nyst@ nswap setcup xst>> union >>xst
  loop nystdrop xst>> ;

{ { 1 2 } { 2 3 4 } { 3 4 5 } } capgen ndup cr set.
{{2},0,{3,4},{1,2},{2,3,4},{3,4,5}} ok
cupgen ndup cr set.
{{1,2,3,4},{1,2,3,4,5},{2,3,4,5},{2},0,{3,4},{1,2},{2,3,4},{3,4,5}} ok
{ 1 2 3 4 5 } { { 1 2 } { 2 3 4 } { 3 4 5 } } topology cr set.
{{1,2,3,4},{2,3,4,5},{2},{3,4},{1,2},{2,3,4},{3,4,5},{1,2,3,4,5},0} ok

: topology \ X s -- λ
  nover nswap capgen
  cupgen nswap incl 0 incl ;

In the figure the rectangle correspond to X and the colored shapes to the set of sets generating a topology. The colored surfaces are the primitive intersections, and the topology is the set of all combinations of unions of those intersections.

The next word computes the closure to the singleton {x}.

: singlecl \ λ x -- s  
  0 false loc{ x y flag } >>yst 0 >xst
  nyst@ nunion foreach
  ?do to y true to flag nyst@ foreach
     ?do ndup y member
        if x member 0=
           if false to flag then
        else ndrop
        then
     loop flag if xst>> y incl >>xst then
  loop nystdrop xst>> ;

And the λ-closure for any set:

: closure ( λ s -- s')
  nswap >>xst 0 >yst foreach
  ?do >r nxst@ r> singlecl yst>> union >>yst
  loop nxstdrop yst>> ;

{ 1 2 3 4 5 } { { 1 2 } { 2 3 4 } { 3 4 5 } } topology   ok
{ 1 3 5 } closure cr set.
{1,4,3,5} ok

So except for the members 1,3,5 also 4 is close to {1,3,5}.
The word opsubs computes all open sets that are subsets of s.

: opsubs ( λ s -- s' )
  >>yst 0 >xst foreach
  ?do ndup nyst@ subset
     if fence xst>> union >>xst
     else ndrop
     then
  loop nystdrop xst>> ;

The interior of a set s is the set of all points which are members of an open set that is a subset of s. That is, s' is the union of all open subsets of s.

: interior ( λ s --- s')
  opsubs nunion ;

{ 1 2 3 4 5 } { { 1 2 } { 2 3 4 } { 3 4 5 } } topology ndup  ok
{ 1 3 5 } interior cr set.
0  ok
{ 2 3 } interior cr set.
{2} ok