bouncing ball with Lite HyEQ Solver

Consider a model for a bouncing ball system where {x_1} is the distance of the ball regarding the floor and {x_2} is its speed. Now consider the hybrid model that describes the bouncing ball written as

Example_1_20x

where {\gamma >0} is the gravity constant and {\lambda \in [0,1)} is the restitution coefficient. Consider the ball to be bouncing on a floor at zero height. The constants for the bouncing ball system are {\gamma = 9.81} and {\lambda=0.8}.
The following procedure is used to simulate this example in the Lite HyEQ Solver:

  • Inside the MATLAB script run.m, initial conditions, simulation horizons, a rule for jumps, ode solver options, and a step size coefficient are defined. The function HyEQsolver.m is called in order to run the simulation, and a script for plotting solutions is included.
  • Then the MATLAB functions f.m, C.m, g.m, D.m are edited according to the data given above.
  •  Finally, the simulation is run by clicking the run button in run.m or by calling run.m in the MATLAB command window.

The simulation results are shown in the following images, and the example code for each of the MATLAB files run.m, f.m, C.m, g.m, and D.m is given below.

fig1BouncingBallfig2BouncingBall

Simulation code


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Matlab M-file       Project: HyEQ Toolbox  @ Hybrid Dynamics and Control
% Lab, http://www.u.arizona.edu/~sricardo/index.php?n=Main.Software
%
% Filename: run.m
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Bouncing ball example - Lite simulator version
%% initial conditions
clear all
close all
clc

x1_0 = 1;
x2_0 = 0;
x0 = [x1_0;x2_0];
% simulation horizon
TSPAN=[0 10];
JSPAN = [0 20];
% rule for jumps
% rule = 1 -> priority for jumps
% rule = 2 -> priority for flows
rule = 1;
options = odeset('RelTol',1e-6,'MaxStep',.1);
% constants
global gamma lambda;
    gamma = 9.81; %# gravity constant
    lambda = 0.8; % restitution coefficient

%% simulate
[t,j,x] = HyEQsolver(@f,@g,@C,@D,x0,TSPAN,JSPAN,rule,options);

%% plot solution
figure(1) % position
clf
subplot(2,1,1),plotHarc(t,j,x(:,1))
grid on
ylabel('x1')
subplot(2,1,2),plotjumps(t,j,x(:,1))
grid on
ylabel('x1')
figure(2) % velocity
clf
subplot(2,1,1),plotflows(t,j,x(:,2))
grid on
ylabel('x2')
subplot(2,1,2),plotjumps(t,j,x(:,2))
grid on
ylabel('x2')
% plot hybrid arc
plotHybridArc(t,j,x)
xlabel('j')
ylabel('t')
zlabel('x1')

Flow map


function xdot = f(x)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Matlab M-file       Project: HyEQ Toolbox  @ Hybrid Dynamics and Control
% Lab, http://www.u.arizona.edu/~sricardo/index.php?n=Main.Software
%
% Filename: f.m
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Description: Flow map
%
% Version: 1.0
% Required files: -
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% constants
global gamma;

% state
x1 = x(1);
x2 = x(2);
% differential equations
xdot = [x2 ; -gamma];
end

Flow set


function [value] = C(x)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Matlab M-file       Project: HyEQ Toolbox  @ Hybrid Dynamics and Control
% Lab, http://www.u.arizona.edu/~sricardo/index.php?n=Main.Software
%
% Filename: C.m
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Description: Flow set
%
% Version: 1.0
% Required files: -
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

x1 = x(1);
if x1 >= 0% flow condition
    value = 1;% report flow
else
    value = 0;% do not report flow
end
end

Jump map


function xplus = g(x)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Matlab M-file       Project: HyEQ Toolbox  @ Hybrid Dynamics and Control
% Lab, http://www.u.arizona.edu/~sricardo/index.php?n=Main.Software
%
% Filename: g.m
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Description: Jump map
%
% Version: 1.0
% Required files: -
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% constants
global lambda;

% state
x1 = x(1);
x2 = x(2);
xplus = [-x1 ; -lambda*x2];
end

Jump set


function inside = D(x)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Matlab M-file       Project: HyEQ Toolbox  @ Hybrid Dynamics and Control
% Lab, http://www.u.arizona.edu/~sricardo/index.php?n=Main.Software
%
% Filename: D.m
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Description: Jump set
%
% Version: 1.0
% Required files: -
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

x1 = x(1);
x2 = x(2);
if (x1 <= 0 && x2 <= 0)% jump condition
    inside = 1;% report jump
else
    inside = 0;% do not report jump
end
end

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s