Menu

PyCont: Convergence problem for a frequency response

Help
Anonymous
2014-05-03
2014-05-04
  • Anonymous

    Anonymous - 2014-05-03

    Hi,

    I'm trying to reproduce a frequency response plot for a duffing equation using Pycont - see here.

    However my code fails to converge, when the response becomes partly unstable (when damping is below beta<0.10).

    I hope someone can give me a hint of what I can do about this.

    Thanks,
    Henrik

    My code is found below:

    from PyDSTool import *
    import matplotlib.pyplot as plt
    
    # simulation parameters
    maxNumPoints = 1e5
    maxStepSize = 1e-2
    stepSize    = 1e-3
    minStepSize = 1e-6
    
    # physical parameters
    s = -1.0 # sigma = Omega-omega0 (detuning parameter)
    o = 1.0  # omega0
    g = 0.5  # gamma
    q = 0.2  # q
    b = 0.2  # beta
    
    # Equilibrium curve
    ustr = '-b*o*u+q/(2*o)*sin(v)' # amplitude
    vstr = 's*u-3*g/(8*o)*u*u*u+q/(2*o)*cos(v)' # phase
    
    DSargs = args(name='duffing')
    DSargs.pdomain = {'s': [-1.0, 1.0]}  # domain boundaries
    DSargs.pars = {'s': s, 'o': o, 'g': g, 'q': q, 'b': b}
    DSargs.varspecs = {'u': ustr, 'v': vstr}
    DSargs.ics = {'u': 0.01, 'v': 1.0}
    
    testDS = Generator.Vode_ODEsystem(DSargs)
    
    PyCont = ContClass(testDS)
    
    PCargs = args(name='EQ1', type='EP-C')
    PCargs.freepars = ['s']
    PCargs.StepSize = stepSize
    PCargs.MaxNumPoints = maxNumPoints
    PCargs.MaxStepSize = maxStepSize
    PCargs.MinStepSize = minStepSize
    PCargs.LocBifPoints = ['LP', 'BP', 'B']
    PCargs.StopAtPoints = ['B']
    PCargs.verbosity = 2
    PCargs.SaveEigen = True
    PCargs.Corrector = 'Natural'
    PyCont.newCurve(PCargs)
    
    print 'Computing first equilibrium curve...'
    start = clock()
    PyCont['EQ1'].forward()
    print 'done in %.3f seconds!' % (clock()-start)
    
    # Plots
    fig1 = plt.figure(num=1)
    PyCont.display(('s', 'u'), stability=True, linewidth=0.5)
    plt.grid(True)
    plt.title('amplitude')
    ymin, ymax = plt.ylim()
    plt.ylim([-ymax, ymax])
    fig1.show()
    
     

    Last edit: Anonymous 2014-05-03
  • Drew LaMar

    Drew LaMar - 2014-05-03

    Hi, Henrik. The MoorePenrose corrector is much more stable. I changed the corrector as follows:

     PCargs.Corrector = 'MoorePenrose'
    

    then modified beta to 0.05 to get the attached figure.

    • Drew
     
  • Anonymous

    Anonymous - 2014-05-04

    Thanks a lot!

    Can I find an overview of available correctors anywhere in the documentation?

     

Log in to post a comment.