## Fixpoints in nD

Help
Greg
2012-08-22
2013-11-04
• Greg - 2012-08-22

Hello again!

I now try to extend my analysis to 3dimensional systems,
but I have Problems to call the "fixedpoint_nD" method.
It took me ages to find out how to handle the already defined
Jacobian to that method (wonder why it not just check the 'haveJacobian' attribute and then just uses it)
, but still it won't take it :(

Heres is the code:

```from PyDSTool import *
from PyDSTool.Toolbox import phaseplane as pp
from PyDSTool.Toolbox.phaseplane import *
#==========SET UP SYSTEM==========================
#k2=4 !
pars = {'alpha' : 0.6, 's' : 0, 'beta' : 0.2, 'beta2' : 0.1}
icdict={'x' : 0.7 , 'y' : 0.1, 'x0' : 0.4}
# Set up model
Xeq = '1 - alpha*x*y - 0.1*x - s*x'
X0eq = '-alpha*0.1*x0*y + s*x'
Yeq = 'beta*x/(1+x) + beta2*x0/(1+x0) - 0.1*y'

DSargs = args(name='3dim')
DSargs.pars = pars
DSargs.varspecs = {'x': Xeq, 'y': Yeq, 'x0' : X0eq}
DSargs.ics = icdict         #initial conditions
DSargs.tdomain = [0,300]     # set the range of integration.
DSargs.xdomain = {'x': [0, 15], 'y': [0, 5.], 'x0' : [0,15]}
DSargs.fnspecs = {'Jacobian': (['t','x','x0','y'],"""[[-alpha*y-0.1-s,0,-alpha*x],[s,-alpha*0.1*y,-alpha*0.1*x0],[beta/(1+x)
-beta*x/pow((1+x),2),beta2/(1+x0)-beta2*x0/pow((1+x0),2),-0.1]]""")}
#==========================================0
ode = Generator.Dopri_ODEsystem(DSargs)
#==========================================0
print ode.funcspec.auxfns['Jacobian'][0]
#============PHASEPLANE=================================
#drawpp(ode)
fp_coord = pp.find_fixedpoints(ode, n=8, eps=1e-13, jac = 'Jacobian')
print fp_coord
fps=[]
#print dir(ode.auxfns)
print ode.haveJacobian
for i in range(len(fp_coord)):
fps.append(fixedpoint_nD(ode, Point(fp_coord[i]), coords=['x', 'x0','y'],jac=ode.auxfns['Jacobian'],
description='bottom', eps=1e-12))
```

That is the output:

"phaseplane.py", line 1857, in __init__
jac_test_arg[self.fp_coords]))
File "<string>", line 4, in Jacobian
File "<string>", line 2, in _auxfn_Jac
IndexError: invalid index to scalar variable."

I looked into phaseplane.py, and it seems like the Jacobian is treated as an 2-dim one?!
Any ideas would be great!

Greets,
Gregor

• Greg - 2012-08-22

uuh..it messed the code up, sorry for that.
Here is the plain repost:

from PyDSTool import *
from PyDSTool.Toolbox import phaseplane as pp
from PyDSTool.Toolbox.phaseplane import *

#==========SET UP SYSTEM==========================
#k2=4 !

pars = {'alpha' : 0.6, 's' : 0, 'beta' : 0.2, 'beta2' : 0.1}

icdict={'x' : 0.7 , 'y' : 0.1, 'x0' : 0.4}

# Set up model
Xeq = '1 - alpha*x*y - 0.1*x - s*x'
X0eq = '-alpha*0.1*x0*y + s*x'
Yeq = 'beta*x/(1+x) + beta2*x0/(1+x0) - 0.1*y'
#Beq = 'beta*A-gam*B'  #mdm protein

DSargs = args(name='3dim')
DSargs.pars = pars
DSargs.varspecs = {'x': Xeq, 'y': Yeq, 'x0' : X0eq}
DSargs.ics = icdict         #initial conditions
DSargs.tdomain =      # set the range of integration.
DSargs.xdomain = {'x': , 'y': , 'x0' : }
DSargs.fnspecs = {'Jacobian': (,"""[,,[beta/(1+x)\
-beta*x/pow((1+x),2),beta2/(1+x0)-beta2*x0/pow((1+x0),2),-0.1]]""")}

#==========================================0
ode = Generator.Dopri_ODEsystem(DSargs)
#==========================================0

print ode.funcspec.auxfns
#============PHASEPLANE=================================

#drawpp(ode)

fp_coord = pp.find_fixedpoints(ode, n=8, eps=1e-13, jac = 'Jacobian')

print fp_coord
fps=
#print dir(ode.auxfns)
print ode.haveJacobian

for i in range(len(fp_coord)):
fps.append(fixedpoint_nD(ode, Point(fp_coord_), coords=,jac=ode.auxfns,
description='bottom', eps=1e-12))
_

• Greg - 2012-08-22

..? OK last try:

from PyDSTool import *
from PyDSTool.Toolbox import phaseplane as pp
from PyDSTool.Toolbox.phaseplane import *
from tools import *

#==========SET UP SYSTEM==========================
#k2=4 !

pars = {'alpha' : 0.6, 's' : 0, 'beta' : 0.2, 'beta2' : 0.1}

icdict={'x' : 0.7 , 'y' : 0.1, 'x0' : 0.4}

# Set up model
Xeq = '1 - alpha*x*y - 0.1*x - s*x'
X0eq = '-alpha*0.1*x0*y + s*x'
Yeq = 'beta*x/(1+x) + beta2*x0/(1+x0) - 0.1*y'
#Beq = 'beta*A-gam*B'  #mdm protein

DSargs = args(name='3dim')
DSargs.pars = pars
DSargs.varspecs = {'x': Xeq, 'y': Yeq, 'x0' : X0eq}
DSargs.ics = icdict         #initial conditions
DSargs.tdomain =      # set the range of integration.
DSargs.xdomain = {'x': , 'y': , 'x0' : }
DSargs.fnspecs = {'Jacobian': (,"""[,,
]""")}

#==========================================0
ode = Generator.Dopri_ODEsystem(DSargs)
#==========================================0

print ode.funcspec.auxfns
#============PHASEPLANE=================================

#drawpp(ode)

fp_coord = pp.find_fixedpoints(ode, n=8, eps=1e-13, jac = 'Jacobian')

print fp_coord
fps=
#print dir(ode.auxfns)
print ode.haveJacobian

for i in range(len(fp_coord)):
fps.append(fixedpoint_nD(ode, Point(fp_coord_), coords=,jac=ode.auxfns,
description='bottom', eps=1e-12))
_

• Greg - 2012-08-22

The line where the Jacobian gets defined gets screwed up here, I don't know why.
I try to give it line by line:

```DSargs.fnspecs = {'Jacobian': (['t','x','x0','y'],"""[[-alpha*y-0.1-s,0,-alpha*x],
[s,-alpha*0.1*y,-alpha*0.1*x0],
[beta/(1+x)-beta*x/pow((1+x),2),beta2/(1+x0)-beta2*x0/pow((1+x0),2),-0.1]]""")}
```

The rest should be ok!

• Greg - 2012-08-22

:(

• Rob Clewley - 2012-08-22

Hi,

I think you have accidentally inserted a backspace character in your source code where it starts becoming "strike through", as that's what is upsetting the auto-formatter. If you put your code into a raw programming editor you should be able to see where to edit it out.

Anyway, a couple of things about the code. Your Jacobian should be 3 x 3 but it's only 2 x 3. It has to have entries for every dynamic variable.

That said, I've been working on the PP code recently and I noticed some problems with the n-D case. I tried fixing it already. So could you try updating your PyDSTool to the latest from the Head of the BZR repository (see the web page under Code Topics for the link, and http://doc.bazaar.canonical.com/latest/en/mini-tutorial/), making the Jac 3x3 then getting back to me if there are still problems?

Rob

• Greg - 2012-08-23

Hi Rob,

ok I opened it with an raw editor:

from PyDSTool import *
from PyDSTool.Toolbox import phaseplane as pp
from PyDSTool.Toolbox.phaseplane import *

#==========SET UP SYSTEM==========================
#k2=4 !

pars = {'alpha' : 0.6, 's' : 0.1, 'beta' : 0.2, 'beta2' : 0.1}

icdict={'x' : 0.7 , 'y' : 0.1, 'x0' : 0.4}

# Set up model
Xeq = '1 - alpha*x*y - 0.1*x - s*x'
X0eq = '-alpha*0.1*x0*y + s*x'
Yeq = 'beta*x/(1+x) + beta2*x0/(1+x0) - 0.1*y'
#Beq = 'beta*A-gam*B'  #mdm protein

DSargs = args(name='3dim')
DSargs.pars = pars
DSargs.varspecs = {'x': Xeq, 'y': Yeq, 'x0' : X0eq}
DSargs.ics = icdict         #initial conditions
DSargs.tdomain =      # set the range of integration.
DSargs.xdomain = {'x': , 'y': , 'x0' : }
DSargs.fnspecs = {'Jacobian': (,"""[,,
]""")}

#==========================================0
ode = Generator.Dopri_ODEsystem(DSargs)
#==========================================0

print ode.funcspec.auxfns
#============PHASEPLANE=================================

#drawpp(ode)

fp_coord = pp.find_fixedpoints(ode, n=8, eps=1e-13, jac = 'Jacobian')

print fp_coord
fps=
#print dir(ode.auxfns)
print ode.haveJacobian

for i in range(len(fp_coord)):
fps.append(fixedpoint_nD(ode, Point(fp_coord_), coords=,jac=ode.auxfns,
description='bottom', eps=1e-12))

So originally the Jacobian was 3x3 of course.
Unfortunately I have problems getting the bzr repos, I tried the comman:

_

``` bzr branch bzr://pydstool.bzr.sourceforge.net/bzrroot/pydstool
```

_
but recieved a time-out. I will try again later!
Thx so far!

Greg

_

• Rob Clewley - 2012-08-23

OK, I got it to run. I don't know why bzr didn't respond. It works for me. You can try pulling individual files from the web client, you should only need the updated phaseplane.py module anyway (although there are a few other recent fixes elsewhere that might be useful to pick up).

Here are some edits at the end of your file. You selected the jacobian in the wrong way. If you want to use the one made for the model then just use its method Jacobian directly. You can't pass the string 'Jacobian' either :)

```fp_coord = pp.find_fixedpoints(ode, n=8, eps=1e-13, jac=ode.Jacobian)
print fp_coord
fps=[]
#print dir(ode.auxfns)
print ode.haveJacobian()
for i in range(len(fp_coord)):
fps.append(fixedpoint_nD(ode, Point(fp_coord[i]), coords=['x', 'x0','y'],
jac=ode.Jacobian,
description='bottom', eps=1e-12))
```

Also, fp_coord is a tuple, so you had to select the 0 element of it to put into Point(). I'll be checking in a minor fix for phaseplane.py this afternoon so try grabbing it in a couple of hours and you should be good to go.
You also needed

• Greg - 2012-08-24

Ok, cool!
Maybe the firewall has sth. against me pulling bzr files?
I'll at home later today! As soon as I have the fixed pp file I'll let you know!

Thx,
Greg

• Greg - 2012-08-27

Hey,

yeah it works. Now I can analyze the evecs and evals in nD.
There is even a fixpoint classification for the nD case implemented as
I see now. So thx again for you help!

Greets,
Greg

PS: bzr worked fine today