As with all projects there is a story with this one.

I was working on my RNN code, hoping to eventually get a GPU version of it running. Somewhere between formulas and implementations I couldn't help but notice convergence and initial conditions are very fragile things. So out of frustration I went and did some research on what could possibly make my RNN a little more stable. Amidst it all I came across some articles explaining the fractal structure of the convergence of recurrent networks. "Egads!" I thought to my self, "No wonder!" Both systems have a persisting state variable that is recursively sent through a (non)linear transform function. The weight matrix of the RNN works as the linear components of the IFS transformation matrix. The RNN activation function works as the final transformation to be applied to the result.

I embarked on designing a flame fractal renderer with the intentions of writing a renderer for the RNN basins of attraction. (Somewhat similar to all those Newton fractals you see out there: originally designed to solve for roots of quintic equations but in the end all they're truly used for is their pretty pictures.) Sadly my GPU greedy nature got the best of me and I went awry and ended up making a GPU-powered flame fractal renderer.

This thing is powered by FBO's and float textures to iterate all its points. It uses PBO's and VBO's to direct the point renderer in their direction. As a result it can do 16k points iterated 20 times a frame at 85fps. Looking at those million-particle simulators on the GPU, I'm sure there's room for improvement.

GPU Flame Fractal - [script]
Dependencies: GL Lua Shell

Thanks much the authors of the flame fractal paper is for all the help it has been in implementing this.

Derivation

The analogy of choosing which IFS branch at random to transform by is simulated by the varying input to the RNN as follows:
Lets say your RNN weight matrix looks like the following:

[new h1]       [w11 w12 w13 w14 w15] [h1]
[new h2]       [w21 w22 w23 w24 w25] [h2]
[  y1  ] = s ( [w31 w32 w33 w34 w35] [x1] )
[  y2  ]       [w41 w42 w43 w44 w45] [x2]
                                     [ 1]

x  = the i'th input signal
 i

h  = the i'th state signal
 i

y  = the i'th output signal
 i

w   = the weight matrix transforming the j'th input into the i'th output
 ij

s(x) = signal activation function 

Now the way the RNN works is you feed in input signals from an external source. The memory is the persistent state information. The two are combined into a vector, multiplied by the weight matrix, sent through some signal activation function, and end up in a vector of the system output signals and the new memory persistent state.

With IFS's, you have no input signals but instead your persistent data is sent through one of a various number of linear systems and transformation functions

[new h1]       [a b c]    [h1]
[new h2] = t ( [d e f]  * [h2] )
            i         i   [ 1]

t  (x) = IFS transform function
 i

For some i in [1,n]

Looks very similar doesn't it? The only difference is that, where the RNN had extra linear terms to represent its response to input, the IFS instead has different linear terms depending on what input is received.

Here's how they tie together:

Let t (x) = s(x) for all i in [1,n] 
     i

(i.e. all IFS transform functions are set to the RNN signal activation function) First we can eliminate the output of the RNN. All the IFS will show us is what happens to the state signals. We don't have to consider the RNN output for the IFS's sake. Our RNN now looks like this:

                                     [h1]
[new h1]       [w11 w12 w13 w14 w15] [h2]
[new h2] = s ( [w21 w22 w23 w24 w25] [x1] )
                                     [x2]
                                     [ 1]

Next we make a set of all possible input vectors which the RNN will accept.

      [x1]
X = { [x2] }
          i

And, for some given i'th input vector, we plug it into whats left of our RNN equation:

[new h1]       [w11 w12 (w13 * x1 + w14 * x2 + w15)] [h1]
[      ]       [                 i          i      ] [  ]
[      ]       [                                   ] [  ]
[new h2] = s ( [w21 w22 (w23 * x1 + w24 * x2 + w25)] [h2] )
               [                 i          i      ] [  ]
               [                                   ] [  ]
               [                                   ] [ 1]

Equating the RNN weight matrix with the IFS linear system gives us the following:

a = w11     b = w12     c = w13 * x1 + w14 * x2 + w15
                                    i          i

d = w21     e = w22     f = w23 * x1 + w24 * x2 + w25
                                    i          i

Take note that as the input signal vector changes only the translation changes while the coefficients of the state do not. All in all what does this tell us about the connection between RNNs and IFSs? It says that a RNN's state signals converge in the same fractal structure as an IFS does, so long as the IFS nonlinear transform function is the RNN signal activation function and so long as every branch in the IFS has the same state coefficients, i.e. the only differences in the linear coefficients of all the IFS branches is their translation.