Solving the Duffing Equation with ODE45

[member 428835]

Hi PF

I am trying to replicate the phase portrait at the bottom right of this page: https://en.wikipedia.org/wiki/Duffing_equation

What I have is this:

global delta alpha beta gamma OMEG
delta   = 0.3;  % DAMPING
alpha   = -1;   % STIFFNESS
beta    = 1;    % RESTORATION NONLINEARITY
OMEG    = 1.2;  % DRIVING FORCE ANGULAR FREQ
gamma   = 0.2;  % FORCING

Fs = 100;       % Sampling frequency                    
T  = 1/Fs;      % Sampling period       
L  = 10000;     % Length of signal

[t, x]  = ode45(@duffing,(0:L-1)*T,[0 1]);

[~, idx] = min( abs( t-OMEG/(2*pi)*38 ) );

plot(x(1:idx,2),x(1:idx,1))
title('phase space')
xlabel('position')
ylabel('velocity')

and the function ode45 calls is

function xdot = duffing(t,x)
global delta alpha beta gamma OMEG
xdot(1) = -delta*x(1) - alpha*x(2) - beta*x(2)^3 + gamma*cos(OMEG*t);
xdot(2) = x(1);
xdot    = xdot';
end

For these parameter values I should get a cyclic phase portrait, but I don't. Can you see why?

 

[member 428835]

Sheesh don't know how to delete this but I just realized they must be going from a long time out, not t=0 as shown, since their IC doesn't align to the phase portrait.