Cleveland State University

Department of Electrical and Computer Engineering

EEC 510

Linear Systems

function Ex03

% function Ex03.

% Inverted pendulum simulation.

% To see the effects of step size, change dt to 0.1 and tf to 10

% and see what happens.

% Close all figures.

close all;

% Set the step size and simulation length.

dt = 0.001;

tf = 1;

thetainit = 0.1; % initial pendulum angle (radians)

m = 0.2; % pendulum mass

M = 1.0; % cart mass

g = 9.81; % acceleration due to gravity

F = 0.1; % coefficient of viscous friction on cart wheels

l = 1; % length of pendulum

r = 0.02; % radius of pendulum mass

J = m * r * r / 2; % moment of inertia of pendulum mass (cylinder)

% Nonlinear simulation.

theta = thetainit;

thetadot = 0;

thetadotdot = 0;

y = 0;

ydot = 0;

ydotdot = 0;

theta_array = [];

thetadot_array = [];

thetadotdot_array = [];

y_array = [];

ydot_array = [];

ydotdot_array = [];

for t = 0 : dt : tf

% Save the old state values.

theta_old = theta;

thetadot_old = thetadot;

thetadotdot_old = thetadotdot;

y_old = y;

ydot_old = ydot;

ydotdot_old = ydotdot;

sine = sin(theta);

cosine = cos(theta);

% Compute the new state values.

u = 100*theta;

theta = theta + thetadot_old * dt;

thetadot = thetadot + thetadotdot_old * dt;

thetadotdot = (m*g*l*sine - m*l*cosine*(u + ...

m*l*thetadot_old*thetadot_old*sine - F*ydot_old) / (M + m)) / ...

(J + m*l*l - m*m*l*l*cosine*cosine / (M + m));

y = y + ydot_old * dt;

ydot = ydot + ydotdot_old * dt;

ydotdot = (u - m*l*thetadotdot_old*cosine + ...

m*l*thetadot_old*thetadot_old*sine - F*ydot_old) / (M + m);

% Save the state values in arrays for later plotting.

theta_array = [theta_array theta];

thetadot_array = [thetadot_array thetadot];

thetadotdot_array = [thetadotdot_array thetadotdot];

y_array = [y_array y];

ydot_array = [ydot_array ydot];

ydotdot_array = [ydotdot_array ydotdot];

end

% Linearized simulation.

xdot = zeros(4, 1);

x = zeros(4, 1);

x(3) = thetainit;

A = [0 1 0 0;

0 -F/M -m*g/M 0;

0 0 0 1;

0 F/M/l (M+m)*g/M/l 0];

B = [0; 1/M; 0; -1/M/l];

xarray = [];

for t = 0 : dt : tf

u = 100*x(3);

xdot = A * x + B * u;

x = x + xdot * dt;

xarray = [xarray x];

end

% Compare the results of the nonlinear and linearized simulations.

t = [0 : dt : tf];

figure;

plot(t, (180 / pi) * theta_array, '-', t, (180 / pi) * xarray(3, :), '--');

title(['Solid = Nonlinear, Dashed = Linearized']);

xlabel('Time');

ylabel('Angle (radians)');

figure;

plot(t, y_array, '-', t, xarray(1, :), '--');

title(['Solid = Nonlinear, Dashed = Linearized']);

xlabel('Time');

ylabel('Position');

---

Professor Simon's Home Page

Department of Electrical and Computer Engineering

Cleveland State University

---

Last Revised: August 29, 2001