%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % solves the harmonic oscillator with damping % % % % x' = y... % % y' = -(k/m)*x - (c/m)*y % % % % requires the use of file: harmosc.m % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% clear; global c; global m; global k; c = .3; m = 1; k = 1; % time span of integration tspan = [0 10]; % initial conditions x0 = 1; y0 = 0; v0 = [x0,y0]; % solves ode using ode45 routine [t,vout] = ode45(@harmosc,tspan,v0); % % plots y(t) for xstar=0.0 in an upper subplot % % note: %%paramstring = sprintf( ', {\\itb}=%d, {\\itc}=%d,{\\itx_0}=%d', b, c, xstar); figure(1) clf subplot(3,1,1) plot( t, vout(:,1) ); xlabel( '{\itt}' ); ylabel( '{\itx}({\itt})' ); %title('Prey') axis tight; grid on; subplot(3,1,2) plot( t, vout(:,2) ); xlabel( '{\itt}' ); ylabel( '{\ity}({\itt})' ); %title( 'Predator') axis tight; grid on; subplot(3,1,3) plot( vout(:,1), vout(:,2) ); xlabel( '{\itx}' ); ylabel( '{\ity}' ); title( 'Phase plane') %axis tight; axis equal; grid on; %figure(2) hbead=line(vout(1,1),vout(1,2),'marker','o','markersize',8,'erase','xor'); %axis([-1 1 -1 1]); %axis('square'); for k=2:length(t) set(hbead,'xdata',vout(k,1),'ydata',vout(k,2)); drawnow pause(.1) end