1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
|
#-----------------------------------------------------------------------------
# Name: Model.py
# Purpose: concentrate model related functions in a single module
#
# Author: Flavio Codeco Coelho
#
# Created: 2006/08/20
# RCS-ID: $Id: Model.py $
# Copyright: (c) 2004-2006
# Licence: GPL
# New field: Whatever
#-----------------------------------------------------------------------------
from scipy import integrate
from numpy import *
from string import *
class Model:
def __init__(self,equations,pars,inits, trange,):
"""
Equations: a function with the equations
inits: sequence of initial conditions
trange: time range for the simulation
"""
self.eqs = equations
self.pars = pars
self.Inits = inits
self.Trange = arange(0,trange,0.1)
self.compileParEqs()
def compileParEqs(self):
"""
compile equation and parameter expressions
"""
#Equations
eql = self.eqs.strip().split('\n')
pl = self.pars.strip().split('\n')
self.vnames = [v.split('=')[0] for v in eql if '=' in v]
self.pnames = [p.split('=')[0] for p in pl if '=' in p]
try:
self.ceqs = [compile(i.split('=')[-1],'<equation>','eval') for i in eql]
except SyntaxError:
dlg = wx.MessageDialog(self, 'There is a syntax error in the equation Box.\nPlease fix it and try again.',
'Syntax Error', wx.OK | wx.ICON_INFORMATION)
try:
dlg.ShowModal()
finally:
dlg.Destroy()
#Parameters
if self.pars.strip() =="":
self.cpars=[]
return
if self.pnames:
#in this case returns the compete expression, including the '='
try:
self.cpars = [compile(i,'<parameter>','exec') for i in pl]
except SyntaxError:
dlg = wx.MessageDialog(self, 'There is a syntax error in the parameter Box.\nPlease fix it and try again.',
'Syntax Error', wx.OK | wx.ICON_INFORMATION)
try:
dlg.ShowModal()
finally:
dlg.Destroy()
else:
try:
self.cpars = [compile(i,'<parameter>','eval') for i in pl]
except SyntaxError:
dlg = wx.MessageDialog(self, 'There is a syntax error in the parameter Box.\nPlease fix it and try again.',
'Syntax Error', wx.OK | wx.ICON_INFORMATION)
try:
dlg.ShowModal()
finally:
dlg.Destroy()
def Run(self):
"""
Do numeric integration
"""
t_courseList = []
t_courseList.append(integrate.odeint(self.Equations,self.Inits,self.Trange,
full_output=0, printmessg=0))
return (t_courseList,self.Trange)
def Equations(self,y,t):
"""
This function defines the system of differential equations, evaluating
each line of the equation text box as ydot[i]
returns ydot
"""
par = self.pars
#---Create Parameter Array----------------------------------------------------------------------------
pars = self.cpars#par.strip().split('\n')
Npar = len(pars)
if self.vnames:
exec('%s=%s'%(','.join(self.vnames),list(y)))
p = zeros((Npar),'d') #initialize p
if pars: #only if there is at least one parameter
for j in xrange(len(pars)):
if self.pnames:
exec(pars[j])
else:
p[j] = eval(pars[j]) #initialize parameter values
#---Create equation array----------------------------------------------------------------------------
eqs = self.ceqs#strip(self.eqs).split('\n')
Neq=len(eqs)
ydot = zeros((Neq),'d') #initialize ydot
for k in xrange(Neq):
ydot[k] = eval(eqs[k]) #dy(k)/dt
return ydot
|
