Ask Your Question
2

tree of vectors

asked 2017-11-07 12:01:32 -0600

Daniel L gravatar image

updated 2017-11-07 12:03:49 -0600

This is a more complicated version of apply functions iteratively

I have a starting "seed" vector, say e1 = matrix(4,1,(1,0,0,0)) Next, I have a collection of 4x4 matrices

thin_group = [T1,T2,T3,T4,T5]

where

$T1 = \left(\begin{array}{rrrr} -1 & 0 & 0 & 0 \\ 2 & 1 & 0 & 0 \\ 4 & 0 & 1 & 0 \\ 4 & 0 & 0 & 1 \end{array}\right), T2 = \left(\begin{array}{rrrr} 1 & 2 & 0 & 0 \\ 0 & -1 & 0 & 0 \\ 0 & 4 & 1 & 0 \\ 0 & 4 & 0 & 1 \end{array}\right), T3 = \left(\begin{array}{rrrr} 1 & 0 & 1 & 0 \\ 0 & 1 & 1 & 0 \\ 0 & 0 & -1 & 0 \\ 0 & 0 & 0 & 1 \end{array}\right), T4 = \left(\begin{array}{rrrr} 1 & 6 & 0 & 6 \\ 0 & 3 & 0 & 2 \\ 0 & 4 & 1 & 4 \\ 0 & -4 & 0 & -3 \end{array}\right), T5 = \left(\begin{array}{rrrr} 3 & 0 & 0 & 2 \\ 6 & 1 & 0 & 6 \\ 4 & 0 & 1 & 4 \\ -4 & 0 & 0 & -3 \end{array}\right)$

These matrices either fix e1 or send to a different vector. I have a "height" function, which just measures how large the vector entries are getting.

def h(X): return max(abs(X[0,0],abs(X[1,0]), abs(X[2,0]),abs(X[3,0]))

I want to develop the following procedure:

hmax = 50 Begin applying T1,...,T5 to e1. As soon as I get a different vector (T1*e1 = (-1,2,4,4)), measure the height of that vector. If the height is below hmax (h(-1,2,4,4) = 4 < hmax = 50), store that vector somewhere. Now apply another matrix from thin_group to (-1,2,4,4), (but not T1 again). If the new vector has a height below hmax, store it, and apply another matrix from thin_group, just not the previous one. Continue this process until we reach a vector with a height larger than hmax. Once this happens, go back to the previous vector and apply a different matrix from thin_group i.e. go up a different branch of the tree. Go up this branch until hmax is exceeded, then go up a different branch, etc.. Eventually, this process will generate all thin_group images of e1 with hmax < 50. (this is probably about 50,000 vectors)

This procedure should "remember" the vectors that were counted/stored, so I'm thinking that I may be able to index each vector in the tree with a sequence (a1,a2,a3,...,an) where aj=1,2,3,4,5 to denote which matrix T1,...,T5 we applied previously.

Any/all ideas or partial solutions greatly appreciated. Thanks!

edit retag flag offensive close merge delete

2 answers

Sort by » oldest newest most voted
4

answered 2017-11-07 14:22:51 -0600

Sébastien gravatar image

updated 2017-11-07 16:37:37 -0600

I think you want to use RecursivelyEnumeratedSet in sage which consist in applying iteratively a successor function to a set of initial seeds:

sage: v = matrix(4,1,(1,0,0,0))
sage: v.set_immutable()
sage: seed = [v]
sage: T1 = matrix( ZZ, 4, [ -1,0,0,0,  2, 1,0,0,  4,0, 1,0,  4, 0,0, 1 ] )
sage: T2 = matrix( ZZ, 4, [  1,2,0,0,  0,-1,0,0,  0,4, 1,0,  0, 4,0, 1 ] )
sage: T3 = matrix( ZZ, 4, [  1,0,1,0,  0, 1,1,0,  0,0,-1,0,  0, 0,0, 1 ] )
sage: T4 = matrix( ZZ, 4, [  1,6,0,6,  0, 3,0,2,  0,4, 1,4,  0,-4,0,-3 ] )
sage: T5 = matrix( ZZ, 4, [  3,0,0,2,  6, 1,0,6,  4,0, 1,4, -4, 0,0,-3 ] )
sage: thin_group = [T1, T2, T3, T4, T5]
sage: def succ(u):
....:     ans = []
....:     for m in thin_group:
....:         m_u = m*u
....:         m_u.set_immutable()
....:         ans.append(m_u)
....:     return ans
....: 
sage: R = RecursivelyEnumeratedSet(seed, succ, max_depth=4)
sage: R
A recursively enumerated set (breadth first search) with max_depth=4

To get a finite set, you may want to use max_depth as above or hmax as you suggest:

sage: def the_rec_enum_set(hmax):
....:     def succ(u):
....:         ans = []
....:         for m in thin_group:
....:             m_u = m*u
....:             m_u.set_immutable()
....:             if max(map(abs, m_u)) < hmax:
....:                 ans.append(m_u)
....:         return ans
....:     R = RecursivelyEnumeratedSet(seed, succ)
....:     return R
....: 
sage: the_rec_enum_set(20)
A recursively enumerated set (breadth first search)

The recursively enumerated set can be enumerated in different maner (look at the documentation). What I like the most is to create the generated graph:

sage: G = the_rec_enum_set(20).to_digraph(); G
Looped multi-digraph on 10 vertices
sage: G = the_rec_enum_set(50).to_digraph(); G
Looped multi-digraph on 37 vertices
sage: G = the_rec_enum_set(500).to_digraph(); G
Looped multi-digraph on 775 vertices

The above numbers confirm what dan_fulea obtain.

Then, you may want to use my optional package to visualize it (with dot2tex and graphviz installed):

sage -pip install slabbe
sage -i dot2tex

and then:

sage: from slabbe import TikzPicture
sage: G = the_rec_enum_set(20).to_digraph()
sage: tikz = TikzPicture.from_graph(G)
sage: tikz.pdf()        # or tikz.png() or tikz.svg()

image description

Looking at the picture, your structure looks to be symmetric (y is a successor of x if and only if x is a successor of y). If so, use

sage: R = RecursivelyEnumeratedSet(seed, succ, structure='symmetric')

which allows to enumerate all elements without keeping all of them in memory. Even better, if your structure is really a tree or a forest (that is if you are able to avoid backward arrows and any loop in the successor function), then use

sage: R = RecursivelyEnumeratedSet(seed, succ, structure='forest')

which will allow you to enumerate the vectors in parallel using fast map reduce algorithms written by Florent Hivert without keeping the elements in memory. In any case, first consult the documentation I have written:

sage: RecursivelyEnumeratedSet?

One information that is missing in the graph is the label of the edges. One way to solve this is to momentarily hide the label information inside the arrival vertex. Once, the graph is made, we put back the information into the edge. I usually do like this:

sage: seed = [(v,None)]
sage: def the_rec_enum_set(hmax):
....:     def succ(node):
....:         u,label = node
....:         ans = []
....:         for i,m in zip([1,2,3,4,5], thin_group):
....:             m_u = m*u
....:             m_u.set_immutable()
....:             if max(map(abs, m_u)) < hmax:
....:                 ans.append((m_u,i))
....:         return ans
....:     R = RecursivelyEnumeratedSet(seed, succ)
....:     return R
....: 
sage: R = the_rec_enum_set(20)
sage: G = R.to_digraph()
sage: G
Looped multi-digraph on 22 vertices
sage: from slabbe.graph import digraph_move_label_to_edge
sage: GG = digraph_move_label_to_edge(G)
sage: GG
Looped digraph on 10 vertices
sage: TikzPicture.from_graph(GG).pdf()

graph with labels indicating which matrix was used

I intend to move this feature directly in Sage once I have time. Tell me if you have any other comment, I intend to improve this module in a near future.

edit flag offensive delete link more

Comments

Very nice code, thank you. Can I ask, in line 4, why did you use random matrices?

Daniel L gravatar imageDaniel L ( 2017-11-07 15:49:01 -0600 )edit

That was my first attempt. I edited my answer later on using your matrices instead of random ones, but I forgot at some point to edit that line. It should be ok now.

Sébastien gravatar imageSébastien ( 2017-11-07 15:52:26 -0600 )edit

Ok got it. In some cases, my thin_group may have loops and in other cases not. The visualization is very nice!

Daniel L gravatar imageDaniel L ( 2017-11-07 16:18:12 -0600 )edit

On the last line tikz.pdf(), I don't get any output. I have dot2tex, graphviz, and slabbe installed. Any idea why this may be so?

Daniel L gravatar imageDaniel L ( 2017-11-09 14:31:28 -0600 )edit

Maybe it depends on the way you use Sage? On the terminal, it works for me. Also in the ipython notebook, you may try:

tikz = TikzPicture.from_graph(G)
path_to_file = tikz.png()
from ipython.display import Image
Image(path_to_file)
Sébastien gravatar imageSébastien ( 2017-11-11 08:25:27 -0600 )edit
2

answered 2017-11-07 14:40:36 -0600

dan_fulea gravatar image

My generating statistics differs, so i hope i understood the question well. Here is the quick code:

T1 = matrix( ZZ, 4, [ -1,0,0,0,  2, 1,0,0,  4,0, 1,0,  4, 0,0, 1 ] )
T2 = matrix( ZZ, 4, [  1,2,0,0,  0,-1,0,0,  0,4, 1,0,  0, 4,0, 1 ] )
T3 = matrix( ZZ, 4, [  1,0,1,0,  0, 1,1,0,  0,0,-1,0,  0, 0,0, 1 ] )
T4 = matrix( ZZ, 4, [  1,6,0,6,  0, 3,0,2,  0,4, 1,4,  0,-4,0,-3 ] )
T5 = matrix( ZZ, 4, [  3,0,0,2,  6, 1,0,6,  4,0, 1,4, -4, 0,0,-3 ] )

T_List = ( T1, T2, T3, T4, T5 )
v = vector( ZZ, 4, [ 1,0,0,0 ] )
v = v.column()

def ht( w ):
    return max( [ abs(entry) for entry in w ] )

def fillList( T_List, vector_List, height ):
    v_List = vector_List[:]    # and we will append
    v_New  = vector_List[:]    # use them to generate new vectors, this is the "last wave"

    count = 0    # give me some small control...
    while v_New and count < 100 :
        count += 1    # this will protect me against the usual programming errors
        w_New = []    # and we will append, this is the "new wave" 
        for v in v_New:
            for T in T_List:
                w = T*v
                if w not in v_List and ht(w) < height:
                    w_New.append( w )
        v_List += w_New
        L = len( w_New )
        print "%s new entr%s added..." % ( L, 'y was' if L == 1 else 'ies were' )
        v_New = w_New[:]

    print "A list with %s vectors was generated in %s steps" % ( len( v_List ), count )
    return v_List


for height in ( 20, 50, 500 ):
    print 60*'-'
    print "Generation for height %s" % height
    v_List = fillList( T_List, [ v, ], height )

Results:

------------------------------------------------------------
Generation for height 20
2 new entries were added...
5 new entries were added...
2 new entries were added...
0 new entries were added...
A list with 10 vectors was generated in 4 steps
------------------------------------------------------------
Generation for height 50
2 new entries were added...
8 new entries were added...
16 new entries were added...
7 new entries were added...
2 new entries were added...
1 new entry was added...
0 new entries were added...
A list with 37 vectors was generated in 7 steps
------------------------------------------------------------
Generation for height 500
2 new entries were added...
8 new entries were added...
32 new entries were added...
116 new entries were added...
232 new entries were added...
183 new entries were added...
85 new entries were added...
42 new entries were added...
21 new entries were added...
15 new entries were added...
10 new entries were added...
6 new entries were added...
5 new entries were added...
5 new entries were added...
5 new entries were added...
2 new entries were added...
1 new entry was added...
1 new entry was added...
1 new entry was added...
1 new entry was added...
1 new entry was added...
0 new entries were added...
A list with 775 vectors was generated in 22 steps
sage:

As i said, possibly i am missing something... (Have to post and catch that train.)

edit flag offensive delete link more

Comments

Dan: very nice code. I'm still working on checking it with what I have here, but it appears to agree with my hand calculation for vectors of ht < 20. Thank you

Daniel L gravatar imageDaniel L ( 2017-11-07 15:47:21 -0600 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2017-11-07 12:01:32 -0600

Seen: 65 times

Last updated: Nov 07