Hopf with positive Lyapunov coef. - but followed by stable limit cycles
Numerical Bifurcation Analysis Toolbox in Matlab
Brought to you by:
hilmeijer,
willy_govaerts
Dear all,
during equilibrium continuation MatCont informs me of a Hopf bifurcation very near to the parameter value that literature on this particular system suggests. Now the limit cycle at a parameter value 1% higher is stable if perturbed (see figure) also in accordance with bibliography, which seems to be a supercritical bifurcation.
However, the Lyapunov coefficient at the Hopf point is positive:
label = H , x = ( 0.000127 -0.000000 -0.000106 0.000000 466.900186 )
First Lyapunov coefficient = 8.584801e+10
How can this be possibly explained? Shouldn't the Lyapunov coefficient at the bifurcation be negative when stable limit cycles follow?
Thanks
Kind regards
lanast
Last edit: lanast 2020-03-03
Theoretically, yes, BUT
Your Lyapunov coefficient is also extremely large. Do you actually trust it? Are the computation settings reliable? Besides, how can we judge that if the system is not provided?
Cheers, Hil
From: lanast lanast@users.sourceforge.net
Sent: Tuesday, March 3, 2020 7:32 PM
To: [matcont:discussion]
Subject: [matcont:discussion] Hopf with positive Lyapunov coef. - but followed by stable limit cycles
Hello together,
During equilibrium continuation MatCont informs me of a Hopf bifurcation very near to the parameter value that literature on this particular system suggests. Now the limit cycles for a parameter value 1% higher than at is stable (see figure, also in accordance with bibliography, which seems to be a supercritical bifurcation
However, the Lyapunov coefficient at the Hopf point is positive:
label = H , x = ( 0.000127 -0.000000 -0.000106 0.000000 466.900186 )
First Lyapunov coefficient = 8.584801e+10
How can this be explained? Shouldn't the Lyapunov coefficient at the bifurcation be negative when stable limit cycles follow?
Kind regards
lanast
Hopf with positive Lyapunov coef. - but followed by stable limit cycles
Sent from sourceforge.net because you indicated interest in https://sourceforge.net/p/matcont/discussion/762214/
To unsubscribe from further messages, please visit https://sourceforge.net/auth/subscriptions/
Thanks you for this adivce Hill,
You are right, the magnitude of the coefficient is very high.
The system is actually a mass under the influece of gravity and nonlinear hydrodynamic forces. I am using the command line version, here is the dynamic system file. Inside the .m-file I am calling the "hydrodyn" at each iteration.
function out = RotorFile
out{1} = @init;
out{2} = @funeval;
out{3} = [];
out{4} = [];
out{5} = [];
out{6} = [];
out{7} = [];
out{8} = [];
out{9} = [];
% --------------------------------------------------------------------------
function dydt = funeval(t,kmrgd,parOm)
xj = kmrgd(1);
yj = kmrgd(3);
xjt= kmrgd(2);
yjt= kmrgd(4);
[FY1,FZ1] = hydrodyn(xj,yj,xjt,yjt);
dydt(1) = kmrgd(2);
dydt(2) = 2FZ1/(M);
dydt(3) = kmrgd(4);
dydt(4) = - gi + 2FY1/(M);
% --------------------------------------------------------------------------
function [tspan,y0,options] = init
handles = feval(RotorFile);
y0=[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0];
options = odeset('Jacobian',[],'JacobianP',[],'Hessians',[],'HessiansP',[]);
tspan = [0 10];
% --------------------------------------------------------------------------
function jac = jacobian(t,kmrgd,parOm)
% --------------------------------------------------------------------------
function jacp = jacobianp(t,kmrgd,parOm)
% --------------------------------------------------------------------------
function hess = hessians(t,kmrgd,parOm)
% --------------------------------------------------------------------------
function hessp = hessiansp(t,kmrgd,parOm)
%---------------------------------------------------------------------------
function tens3 = der3(t,kmrgd,parOm)
%---------------------------------------------------------------------------
function tens4 = der4(t,kmrgd,parOm)
%---------------------------------------------------------------------------
function tens5 = der5(t,kmrgd,parOm)**
The settings are following:
opt=contset;
opt=contset(opt,'Singularities' , 1);
opt=contset(opt,'InitStepsize', 0.01);
opt=contset(opt,'MinStepsize' , 0.01);
opt=contset(opt,'MaxStepsize' , 0.1);
opt=contset(opt,'MaxNewtonIters', 3);
opt=contset(opt,'MaxCorrIters', 10);
opt=contset(opt,'MaxTestIters', 10);
opt=contset(opt,'MaxNumPoints', 10000);
opt=contset(opt,'Increment', 1e-08);
opt=contset(opt,'Eigenvalues', 1);
Since the dynamic system is kind of cascaded, it would be tricky to run it if I send you the .m-files.
Kind regards
Lysandros
Last edit: lanast 2020-03-03
Difficult to judge,
The function hydrodyn hides the nonlinearities, and I have no idea about values.
But you're not using symbolic derivatives, so it may really be that the computation is unreliable for your system.
Hil
From: lanast lanast@users.sourceforge.net
Sent: Tuesday, March 3, 2020 8:26 PM
To: [matcont:discussion]
Subject: [matcont:discussion] Hopf with positive Lyapunov coef. - but followed by stable limit cycles
Thanks you for this adivce Hill,
You are right, the magnitude of the coefficient is very high.
The system is actually a mass under the influece of gravity and nonlinear hydrodynamic forces. I am using the command line version, here is the dynamic system file. Inside the .m-file I am calling the "hydrodyn" at each iteration.
function out = RotorFile
out{1} = @init;
out{2} = @funeval;
out{3} = [];
out{4} = [];
out{5} = [];
out{6} = [];
out{7} = [];
out{8} = [];
out{9} = [];
% --------------------------------------------------------------------------
function dydt = funeval(t,kmrgd,parOm)
xj = kmrgd(1);
yj = kmrgd(3);
xjt= kmrgd(2);
yjt= kmrgd(4);
[FY1,FZ1] = hydrodyn(xj,yj,xjt,yjt);
dydt(1) = kmrgd(2);
dydt(2) = 2FZ1/((2MJ+MD));
dydt(3) = kmrgd(4);
dydt(4) = - gi + 2FY1/((2MJ+MD));
% --------------------------------------------------------------------------
function [tspan,y0,options] = init
handles = feval(RotorFile);
y0=[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0];
options = odeset('Jacobian',[],'JacobianP',[],'Hessians',[],'HessiansP',[]);
tspan = [0 10];
% --------------------------------------------------------------------------
function jac = jacobian(t,kmrgd,parOm)
% --------------------------------------------------------------------------
function jacp = jacobianp(t,kmrgd,parOm)
% --------------------------------------------------------------------------
function hess = hessians(t,kmrgd,parOm)
% --------------------------------------------------------------------------
function hessp = hessiansp(t,kmrgd,parOm)
%---------------------------------------------------------------------------
function tens3 = der3(t,kmrgd,parOm)
%---------------------------------------------------------------------------
function tens4 = der4(t,kmrgd,parOm)
%---------------------------------------------------------------------------
function tens5 = der5(t,kmrgd,parOm)**
The settings are following:
opt=contset;
opt=contset(opt,'Singularities' , 1);
opt=contset(opt,'InitStepsize', 0.01);
opt=contset(opt,'MinStepsize' , 0.01);
opt=contset(opt,'MaxStepsize' , 0.1);
opt=contset(opt,'MaxNewtonIters', 3);
opt=contset(opt,'MaxCorrIters', 10);
opt=contset(opt,'MaxTestIters', 10);
opt=contset(opt,'MaxNumPoints', 10000);
opt=contset(opt,'Increment', 1e-08);
opt=contset(opt,'Eigenvalues', 1);
Since the dynamic system is kind of cascaded, it would be tricky to run it if I send you the .m-files.
Kind regards
Lysandros
Hopf with positive Lyapunov coef. - but followed by stable limit cycles
Sent from sourceforge.net because you indicated interest in https://sourceforge.net/p/matcont/discussion/762214/
To unsubscribe from further messages, please visit https://sourceforge.net/auth/subscriptions/
Yes the nonlinearities are not shown here, because they are rather lengthly and do result from a iterative loop, so symbolic expression of the derivatives is not practical. I suspect the Jacobian Increment has to be adapted in this case.
I would have another question on this dynamic system. After the Hopf point, I arrive at a point, where the cycles grow increasingly and a Limit Point of Cycles is expected (it is codim-1, so one active parameter). Continuation from there cannot proceed due to the step size, but I am not getting the Limit Point of Cycle notice in the command window. I checked the multipliers and one of them is very near to 1.
Should I reverse the direction of continuation at that point? Or is there a setting I need to configure?
Best, lanast
Hard to judge, of course you have at least one multiplier close to 1. If there is a LPC, then there should be a second one.
Sometimes you can just take a point close to the LPC, and use that to start LPC continuation as there is some initial correction anyway.
Best, Hil
From: lanast lanast@users.sourceforge.net
Sent: Wednesday, March 4, 2020 8:50 PM
To: [matcont:discussion]
Subject: [matcont:discussion] Hopf with positive Lyapunov coef. - but followed by stable limit cycles
Yes the nonlinearities are not shown here, because they are rather lengthly and do result from a iterative loop, so symbolic expression of the derivatives is not practical. I suspect the Jacobian Increment has to be adapted in this case.
I would have another question on this dynamic system. After the Hopf point, I arrive at a point, where the cycles grow increasingly and a Limit Point of Cycles is expected (it is codim-1, so one active parameter). Continuation from there cannot proceed due to the step size, but I am not getting the Limit Point of Cycle notice in the command window. I checked the multipliers and one of them is very near to 1.
Should I reverse the direction of continuation at that point? Or is there a setting I need to configure?
Best, lanast
Attachments:
Hopf with positive Lyapunov coef. - but followed by stable limit cycles
Sent from sourceforge.net because you indicated interest in https://sourceforge.net/p/matcont/discussion/762214/
To unsubscribe from further messages, please visit https://sourceforge.net/auth/subscriptions/
Thanks for reminding me, ofcourse the Floquet multiplier in orbit tagential direction has to be equal to 1 . In my system, another one has reached 1.01352 at the very right point of the figure. The long sought-after LPC however didn't appear and the branches do not reverse backwards...It seem like I have to capture the 1-crossing before it happens, right?
Best regards
Lanast