@@ -0,0 +1,127 @@
#!/usr/bin/octave

printf("Hello, apple!\n")

# lsode_options("integration method", "stiff")
# lsode_options()

# s1_0=5.8;
# s2_0=0.9;
# s3_0=0.2;
# s4_0=0.2;
# n2_0=0.1;
# a3_0=2.4;
# s5_0=0.1;

J0 = 3; # mM/min
k1 = 100.0; # /mM*min
k2 = 6.0; # /mM*min
k3 = 16.0; # /mM*min
k4 = 100.0; # /mM*min
k5 = 1.28; # /min
k6 = 12.0; # /mM*min
k = 1.3; # /min
K = 13.0; # /min
q = 4.0;
K1 = 0.52; # mM
N = 1.0; # mM
A = 4.0; # mM
P = 0.1;

function xdot = g(x,t)
xdot = zeros(9,1); # a column vector
# vdot = zeros(NUMEQNS,1); # a column vector

J0 = 3; # mM/min
k1 = 100.0; # /mM*min
k2 = 6.0; # /mM*min
k3 = 16.0; # /mM*min
k4 = 100.0; # /mM*min
k5 = 1.28; # /min
k6 = 12.0; # /mM*min
k = 1.3; # /min
K = 13.0; # /min
q = 4.0;
K1 = 0.52; # mM
N = 1.0; # mM
A = 4.0; # mM
P = 0.1;

s1 = x(1);
s2 = x(2);
s3 = x(3);
s4 = x(4);
n2 = x(5);
a3 = x(6);
sX = x(7);

a2 = A-a3;
n1 = N-n2;
Ja = K*(s4-sX);

fA3=1+(a3/K1)^q;

v1 = k1*s1*a3/fA3;
v2 = k2*s2*n1;
v3 = k3*s3*a2;
v4 = k4*s4*n2;
v5 = k5*a3;
v6 = k6*s2*n2;
v7 = k*sX;

xdot(1)= J0 - v1;
xdot(2)= 2*v1 - v2 - v6;
xdot(3)= v2 - v3;
xdot(4)= v3 - v4 - Ja;
xdot(5)= v2 - v4 - v6;
xdot(6)= -2*v1 + 2*v3 - v5;
xdot(7)= P*Ja - v7;
xdot(8)= a2;
xdot(9)= n1;
endfunction

s1_0=1.2;
s2_0=0.2;
s3_0=0.05;
s4_0=0.11;
n2_0=0.08;
a3_0=2.5;
s5_0=0.078;
a2_0=A-a3_0;
n1_0=N-n2_0;
g0 = [s1_0, s2_0, s3_0, s4_0, n2_0, a3_0, s5_0, a2_0, n1_0]


# t = linspace (0,10.0,1000);
# y = lsode( "g", g0, t);

# # plot(t,y(:,1),'-',t,y(:,2),'-',t,y(:,3),'-',t,y(:,4),'-',t,y(:,5),'-',t,y(:,6),'-',t,y(:,7),'-')
# plot(t,y(:,1),'.',t,y(:,2),'.',t,y(:,3),'.',t,y(:,4),'.',t,y(:,5),'.',t,y(:,6),'.',t,y(:,7),'.')

# val = input("press any key to continue.");

function gen(numP, iCond, filename)
printf("Generating %d points f\n", numP)
t = linspace (0,10.0,numP);
y = lsode( "g", iCond, t);

printf("writing to file %s\n", filename)
fid = fopen(filename, 'wt');
fprintf(fid,"T\n");
fprintf(fid," S1 S2 S3 S4 N2 A3 S5 A2 N1\n");
for pnt = y
fprintf(fid,"%.6f %.6f %.6f %.6f %.6f %.6f %.6f %.6f %.6f %.6f\n", t, y);
end
fclose(fid);
# plot(t,y(:,1),'.',t,y(:,2),'.',t,y(:,3),'.',t,y(:,4),'.',t,y(:,5),'.',t,y(:,6),'.',t,y(:,7),'.')
printf("\n")
endfunction

# gen(2000,g0,"yeast2000.txt");
# gen(20000,g0,"yeast20000.txt");
# gen(200000,g0,"yeast200000.txt");
gen(1000000,g0,"yeast1000000.txt");


printf("Goodbye!\n")

@@ -52,6 +52,9 @@ def Score(y_true, y_pred):
return (-1, "error: " + str(e))


def Integrate(model, init_vals, time_vals):

pass



@@ -3,6 +3,34 @@
import expand
import tests

### TODOS
#
# - benchmarks
# - functions
# - scikit learn
# - pandas DFs
# - more error metrics
# - diffeqs
# - models
# - evaluation RK4
# - other system types
# - invarients
# - hidden
# - pdes
# - algebra
# - growing / filtering policies
# - simplification / expansions
# - +C ???
# - networkx
# - logging
# - checkpointing
# - statistics
# - Ipython notebook examples
# - abstract expressions / memoization
# - when / where coefficients
# - domain alphabet
# - sub-expression frequencies in population
# - distributing to the cloud

def main():
print "hello pypge!\n"
@@ -48,7 +48,9 @@ def __init__(self,**kwargs):


def check_config(self):
if self.system_type == None or self.search_vars == None or self.usable_vars == None:
if self.system_type == None \
or self.search_vars == None \
or self.usable_vars == None:
return False
return True

@@ -67,7 +69,8 @@ def prepare(self):
self.nsga2_list = []
self.spea2_list = []

for i,e in enumerate(self.grower.first_exprs()):
first_exprs = self.grower.first_exprs()
for i,e in enumerate(first_exprs):
m = model.Model(e)
did_ins = self.memoizer.insert(m)
size = m.size()
@@ -1,8 +1,18 @@
pytest
pyyaml

scipy
numpy
pandas

scikit-learn
networkx

sympy
lmfit
deap

matplotlib
pygraphviz
pydot

@@ -0,0 +1,71 @@
#!/usr/bin/env python
"""
Program to plot the motion of a "springy pendulum".
(kindly taken from: http://julianoliver.com/share/free-science-books/comp-phys-python.pdf [page 102-103])
We actually have FOUR parameters to track, here:
L, L dot, theta, and theta dot.
So instead of the usual Nx2 array, make it Nx4.
Each 4-element row will be used for the state of
the system at one instant, and each instant is
separated by time dt. I'll use the order given above.
"""

import numpy as np
import scipy
from scipy.integrate import odeint

## Nx4
N = 1000 # number of steps to take
y = np.zeros([4])

Lo = 1.0 # unstretched spring length
L = 1.0 # Initial stretch of spring
vo = 0.0 # initial velocity
thetao = 0.3 # radians
omegao = 0.0 # initial angular velocity

y[0] = L # set initial state
y[1] = vo
y[2] = thetao
y[3] = omegao
time = np.linspace(0, 25, N)

k = 3.5 # spring constant, in N/m
m = 0.2 # mass, in kg
gravity = 9.8 # g, in m/s^2

def springpendulum(y, time):
"""
This defines the set of differential equations
we are solving. Note that there are more than
just the usual two derivatives!
"""
g0 = y[1]
g1 = (Lo+y[0])*y[3]*y[3] - k/m*y[0] + gravity*np.cos(y[2])
g2 = y[3]
g3 = -(gravity*np.sin(y[2]) + 2.0*y[1]*y[3]) / (Lo + y[0])
return np.array([g0,g1,g2,g3])


# Now we do the calculations.
answer = scipy.integrate.odeint(springpendulum, y, time)

# Now graph the results.
# rather than graph in terms of t, I'm going
# to graph the track the mass takes in 2D.
# This will require that I change L,theta data
# to x,y data.
xdata = (Lo + answer[:,0])*np.sin(answer[:,2])
ydata = -(Lo + answer[:,0])*np.cos(answer[:,2])


import matplotlib.pyplot as plt

plt.plot(xdata, ydata, 'r-')
plt.xlabel("Horizontal position")
plt.ylabel("Vertical position")

plt.show()