Fixpoints in nD

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

    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

    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

    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

    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

    Greg - 2012-08-22

    :(

     
  • Rob Clewley

    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?

    Thanks for your patience.
    Rob

     
  • Greg

    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

    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

    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

    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

     

Get latest updates about Open Source Projects, Conferences and News.

Sign up for the SourceForge newsletter:





No, thanks