Problem on solving periodic solution of ODE using orthogonal polynomial pseudo-arclength method

Specific Domains Numerics
1

a simple Duffing equation:
dx/dt=xt
dxt/dt=F0u-a1x1-a1tx1t-a3x1^3
du/dt=-omegav+u(1-u^2-v^2)
dv/dt=omegau+v(1-u^2-v^2)
parameters: [omega;F0;a1t;a3]=[2;0.01;1;0.008;0.75]
initial state point: [x0;xt0;u0;v0]=[0;0;1;0]

I am trying to solve iteratively referring to the implementation of Auto and Matcont.
before the first iteration, the discrete state points are given by ode45 (RK4).
after ode45, there are a series of state points [x0;xt0;u0;v0], …, [xN;xtN;uN;vN]
However, after some iterations, the phase diagram made of discrete state points becomes not smooth. And finally becomes the following diagram.
the newton corrector part is converged at first then failed when the phase diagram made of discrete points becomes not smooth at all.

image
image1300×770 24.6 KB

any sugguestions for following problem?

  1. proceeding iterations, the phase diagram made of discrete state points becomes less smooth.
  2. after iterations, the discrete state points don’t correspond to the limit cycle under that parmeter of the present iteration.

ps: my code is on matlab. I’m new on this difficult numerical method, I heard here are many experienced people :smile:

2

For numerical continuation in Julia, BifurcationKit.jl is the premier package. I recommend setting up your continuation using their tools.

For the particular problem you are facing in Matlab, a number of packages exist: Coco, Matcont, the list goes on. If the problem is that the solution becomes non-smooth, that’s typically an indication that your numerical method is unstable. You should consult those packages for help (hint: they do not use ode45, afaik).

3

thanks for your reply! I’ve seen BifurcationKit.jl package these days, and start to learn how to use it.
actually, ode45 is for generating initial periodic orbit. then start numerially conitinuation.

4

Hi

Can you put your equations in latex please to check that you can solve them using BifurcationKit?

In particular is the model time homogenous?

The pictures you are showing are from MatCont and it seems it has trouble near the period doubling point. You can try increase the number Ntst and put mesh adaptation

5

oh, sir, I sent to you an email haha, anyway, you’ve replied to me.
Acturally it’s a duffing equation as follows:

image
image804×149 4.21 KB

And I can solve it in Matcont. In Matcont the phase diagram, in every iteration, is always smooth. But the point is, results calculated by my own code are not smooth.
“In particular is the model time homogenous?”
I think you mean if the time interval is uniform? yes, I have N subintervals and m mesh points in each subinterval, to uniformly discrete time and the corresponding state variables of periodic orbit.

“Can you put your equations in latex please to check that you can solve them using BifurcationKit?
The pictures you are showing are from MatCont and it seems it has trouble near the period doubling point. You can try increase the number Ntst and put mesh adaptation”
The code is written by myself on matlab.
I’ve tried to increase the number Ntst and put mesh adaptation, but it seems to have no effect.
I’m trying to check my code about the part of generating Jacobian matrix of the system. and also I’ll try BifurcationKit these days.

2 Likes