## Wednesday, April 25, 2012

### Reaction Difusion Surface

Just had an amazing stroke of luck!!! found a working Reaction Diffusion code from an internet site
translating it from java to python was almost too easy

the site is:
http://cgjennings.ca/toybox/turingmorph/

and here are the results (on a plane for now)
and a closeup

and here is the python code, its pretty straightforward

enjoy

import math
import rhinoscriptsyntax as rs
import random as rnd

def InitMatrix (rndSeed):
global A0
global B0
rnd.seed(rndSeed)
for i in range(height):
thisA =[]
thisB =[]
for j in range(width):
thisA.append(rnd.random()*12+rnd.random()*2)
thisB.append(rnd.random()*12+rnd.random()*2)
A0.append(thisA)
B0.append(thisB)

def CopyMatrix(src):
res=[]
for i in range(len(src)):
list =[]
for j in range(len(src)):
list.append(src[i][j])
res.append(list)
return res

def NullMatrix():
res=[]
for i in range(height):
list =[]
for j in range(width):
list.append(0)
res.append(list)
return res

def Solve(CA,CB,iterations) :
An = NullMatrix()
Bn = NullMatrix()
global A0
global B0

###uses Euler's method to solve the diff eqns
for n in range(iterations):
for i in range(height):
### treat the surface as a torus by wrapping at the edges
iplus1 = i+1
iminus1 = i-1
if( i == 0 ) : iminus1 = height - 1
if( i == height - 1 ) : iplus1 = 0
for j in range (width):
jplus1 = j+1
jminus1 = j-1
if( j == 0 ) : jminus1 = width - 1
if( j == width - 1 ) : jplus1 = 0
###Component A
DiA = CA * ( A0[iplus1][j] - 2.0 * A0[i][j] + A0[iminus1][j] + A0[i][jplus1] - 2.0 * A0[i][j] + A0[i][jminus1] )
ReA = A0[i][j] * B0[i][j] - A0[i][j] - 12.0
An[i][j] = A0[i][j] + 0.01 * (ReA + DiA)
if( An[i][j] < 0.0 ) : An[i][j] = 0.0
###Component B
DiB = CB * ( B0[iplus1][j] - 2.0 * B0[i][j] + B0[iminus1][j]   + B0[i][jplus1] - 2.0 * B0[i][j] + B0[i][jminus1] )
ReB = 16.0 - A0[i][j] * B0[i][j]
Bn[i][j] = B0[i][j] + 0.01 * (ReB + DiB)
if( Bn[i][j] < 0.0 ) : Bn[i][j]=0.0

A0=CopyMatrix(An)
B0=CopyMatrix(Bn)

def ShowLines(which,spacing,Z):
mat = []
if which=='A' : mat=A0
if which=='B' : mat=B0

rs.EnableRedraw(False)
for i in range(height):
curve =[]
for j in range(width):
curve.append((i*spacing,j*spacing,mat[i][j]*Z))
rs.EnableRedraw()

###main call

width = 200
height = 200
seed = 1
A0 =[]
B0 = []

InitMatrix(seed)
Solve(1,16,500)
ShowLines('A',3,0.5)