I posted yesterday but figured it out; however, a different issue I just detected with the same code arose: namely, why does the solution damp here for an undamped simple harmonic oscillator? I know the exact solution is ##\cos (5\sqrt 2 t)##.
and the function ode45 calls is
global delta alpha beta gamma OMEG delta = 0; % DAMPING alpha = 50; % STIFFNESS beta = 0; % RESTORATION NONLINEARITY OMEG = 1.2; % DRIVING FORCE ANGULAR FREQ gamma = 0; % 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]); plot(t,x(:,2))
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