Menu

KeyError 'EP' occurring when plotting nullclines

Help
2021-03-10
2021-03-10
  • Brittany Howell

    Brittany Howell - 2021-03-10

    Hello. I have a 2D system of ODEs and I am trying to plot the nullclines to then vary my model parameters and plot the effects based on the Zero-crossing events tutorial.

    Based on the current parameter set, I have three steady states (nullclines are verified using XPP but I digress). My goal is to plot multiple nullclines as a parameter (or parameters) vary to see what my model does. I have currently gotten through the portion where the nullclines are plotted but based on the settings, my figure only gives the directional field and sort of white space where my nullclines should be and the KeyError 'EP' (meaning end point I would assume based on AUTO literature). I'm not sure if because I am setting [2] at the end of the fp_coord() function is working, and depending on the model parameter(s) that is(are) varied, I may get only one intersection by the nullclines as opposed to three (I want to plot this as well). I also want to extend my nullclines slightly into the negative by changing the domains of A1 and A2, but then I get a math domain error (although in some circumstances, my nullclines do slope into the negatives, so I would also be interested in fixing this.)

    Any assistance would be fantastic. Thanks so much. You will find my code attached.

    import PyDSTool as dst
    from PyDSTool import args
    import numpy as np
    from matplotlib import pyplot as plt
    from PyDSTool.Toolbox import phaseplane as pp
    
    pars = {'iota':3.09523809523809535, 'upsilon':27.5, 'vartheta':3.09523809523809535, 'varphi1':10, 'varphi2':10}
    icdict = {'A1':0.1, 'A2':0.1}
    A1str = '((1+((10*(A1**2.2))/((A1**2.2)+1)))*iota)-(upsilon*A1)+(varphi1*A2)'
    A2str = '((1+((10*(A2**2.2))/((A2**2.2)+1)))*vartheta)-(upsilon*A1)+(varphi2*A1)'
    
    event_A1_A2 = dst.makeZeroCrossEvent('A1-A2', 0,
                                       {'name': 'event_A1_A2',
                                       'eventtol':1e-6,
                                       'term': False,
                                       'active': True},
                                       varnames=['A1','A2'],
                                       targetlang='python')
    
    DSargs = args(name='twoDsignalmodel')
    DSargs.events = [event_A1_A2]
    DSargs.pars = pars
    DSargs.tdata = [0,3]
    DSargs.algparams = {'max_pts': 1000, 'init_step': 0.01, 'stiff': True}
    DSargs.varspecs = {'A1': A1str, 'A2': A2str}
    DSargs.xdomain = {'A1': [0.0,2.5], 'A2': [0.0,2.5]}
    DSargs.fnspecs = {'Jacobian': (['t','A1','A2'],
                                  """[[(10*2.2*(A1**(2.2-1)))/(((A1**2.2)+1)**2)*iota-upsilon, varphi1],
                                      [varphi2, (10*2.2*(A2**(2.2-1)))/(((A2**2.2)+1)**2)*vartheta-upsilon]]""")}
    DSargs.ics = icdict
    model = dst.Vode_ODEsystem(DSargs)
    
    traj = model.compute('test_traj')
    pts = traj.sample()
    evs = traj.getEvents('event_A1_A2')
    
    plt.figure(1)
    plt.plot(pts['t'], pts['A1'], 'r', linewidth=2)
    plt.plot(pts['t'], pts['A2'], 'b', linewidth=2)
    plt.show()
    
    plt.figure(2)
    pp.plot_PP_vf(model, 'A1', 'A2', scale_exp=-1)
    fp_coord= pp.find_fixedpoints(model, n=10, eps=1e-8)[2]
    fp = pp.fixedpoint_2D(model, dst.Point(fp_coord), eps=1e-8)
    nulls_A1, nulls_A2 = pp.find_nullclines(model, 'A1', 'A2', n=10, eps=1e-8, max_step=0.01, fps=[fp_coord])
    
    # plots fixed point(s)
    pp.plot_PP_fps(fp)
    
    # plots the nullclines
    plt.plot(nulls_A1[:,0], nulls_A1[:,1], 'g')
    plt.plot(nulls_A2[:,0], nulls_A2[:,1], 'y')
    
    # plots the trajectory
    plt.plot(pts['A1'], pts['A2'], 'k-o', linewidth=2)
    
    # plots the event points
    plt.plot(evs['A1'], evs['A2'], 'rs')
    
    plt.axis('tight')
    plt.title('Phase Plane')
    plt.xlabel('A1')
    plt.ylabel('A2')
    plt.show()
    
     
  • Brittany Howell

    Brittany Howell - 2021-03-10

    I'm also curious if this is needed to run the remainder of the tutorial, which is the part I need; to iterate a bunch of phase plane solutions while varying a parameter. Would this error have any effect on that or would it be sufficient to continue and vary a parameter and I should still get the multiple phase plane solutions as required? Thanks!

     
  • Brittany Howell

    Brittany Howell - 2021-03-10

    It is required for the part that I need, so I'm wondering how to fix this. Then when I continue the next part, I get an error of "invalid keys as arguments" for setting tdata=[0,20] so I am unsure how to get this working. Any help would be great. Thank you.

     

Log in to post a comment.