An example of a nonlinear dynamical system with numerous application in engineering is the single well-forced Duffing-Van der Pol oscillator:
d^2y_1/dt^2-/Mu*(1-y_1^2)*dy_1/dt+y_1^3=f*cos(/Omega*t) (1')
where /Mu, f, and /Omega are parameters, y_1 is the variable.
For convenience, the nonlinear second-order function (1') should be transformed to its corresponding state functions with the setting x=y_1, y=dy_1/dt, z=/Omega*t. In this case, an autonomous system is derived, which contains a set of first-order ordinary differential equations, descibed as follows:
/dot(x)=y; (1)
/dot(y)=/Mu*(1-x^2)*y-x^3+f*cos(z) (2)
/dot(z)=/Omega. (3)
1 Phase plane and Phase space
Matlab Codes:
function xp=VanDerPol(t,x)
% Duffing-Van Der Pol equations expressed as a first order system
% (see equations 3-5). Parameter values (mu, f and omega)
% can be chosen to demonstrate various trajectories (see text).
% IN:
% t: time (not used, but necessary for ODE solver)
% x: Input vector
% OUT:
% xp: dx/dt
% Current parameters demonstrate chaotic behavior of the oscillator.
mu=0.2; f=1.0; omega=0.94;
%There are other three sets of parameters should be varify, just do it
%without any hesitate,because different behaviors can be generated!
%(1)mu=0.0;f=0.0;(2)mu=0.2; f=0.0; In these two cases, we only have a phase
%plane with respect to varibals of x(1) and x(2), means phase plane x-y.
%(3)mu=0.2,f=1.0;omega=0.9;(4) Current parameters. In these two cases, we
%can have eigher phase plane with respect two any two varibals of
%x(1),x(2),x(3) or phase space of the above three varibals via differet
%visual angles, however, we should bound the values of variable x(3) for
%its unbounded solution evoles with time t incresing to infinity, like
%sin(x3) or cos(x3),...
% define equations
xp(1)=x(2);
xp(2)=mu*(1-x(1)^2)*x(2)-x(1)^3+f*cos(x(3));
xp(3)=omega;
% transpose into column vector for ODE solver
xp=xp';
% File: VanDerPolSolv.m
% Numerical integration of Duffing-Van Der Pol Oscillator.
% Ordinary differential equation solver using in-built
% Matlab function (ode45).
% Uses: VanDerPol.m
% Set initial conditions for the three differential equations
x0=[1;0;0];
% Set integration time
timePoints=0:0.01:10000;
%timePoints=5000:0.3:7500; (Not smooth enough, as shown in Fig.1)
% Solve equations
options=odeset('RelTol',1e-10,'AbsTol',1e-12);
[t,x]=ode45('VanDerPol',timePoints,x0);
2 Poincare section
% File: SurfaceOfSection.m
% Two dimension Poincare surface of section on plane x3=0
% (data assumed to be stored in matrix 'xnew' generated from PSplot.m).
% Start with empty surface
clear Surface
OldDistance=0;
k=1;
for i=1:size(xnew)*[1;0]
NewDistance=xnew(i,3);
if (NewDistance>=0 && OldDistance<0)
% Add new point to the surface
TotalDistance=NewDistance-OldDistance;
Surface(k,:)=xnew(i-1,:)-(OldDistance/TotalDistance)*(xnew(i,:)-xnew(i-1,:));
k=k+1;
end
OldDistance=NewDistance;
end
% Generate 2D plot
plot(Surface(:,1),Surface(:,2),'.','markersize',1);
xlabel('x1');
ylabel('x2');
title('Poincare Surface of Section');