MATLAB: Bouncing ball’s height and velocity!

animationMATLABplot

I want to make a graph when the ball drop from the height(20m). I insert what I wanted! I want to know what I do wrong.
clear all
h0 = 20;      
v = 0;      
g = 10;         
t=0;
dt = 0.01;     
rho = 0.75;     
tau = 0.10 ;    
hmax = h0   ;   
h = h0;
hstop = 0.01;   
freefall = 1; 
t_last = -sqrt(2*h0/g); 
vmax = sqrt(2 * hmax * g);
H = [];
T = [];
while(hmax > hstop)
    if(freefall==1)
        hnew = h + v*dt - 0.5*g*dt*dt;
    if(hnew<0)
         t = t_last + 2*vmax;
         freefall = 0;
         t_last = t + tau;
         h = 0;
    else
      t = t + dt;
      v = v - g*dt;
      h = hnew;
    end
    else
      t = t + tau;
      vmax = vmax * rho;
      v = vmax;
      freefall = 1;
      h = 0;
    end
    hmax = 0.5*vmax*vmax/g;
  H.append(h);
  T.append(t);
end 
  
%% Simulation
plot(time, height, 'r.', 'MarkerSize', 50);
axis([-2, 20, 0 25]); grid;
xlabel('ball position X [m]')
ylabel('ball position Y [m]')
title('TaengTaeng Ball')
drawnow

Best Answer

Well there was a lot that was wrong with that code. Here is the fix:
% Demo to simulate a bouncing ball.
clc;    % Clear the command window.
fprintf('Beginning to run %s.m.\n', mfilename);
close all;  % Close all figures (except those of imtool.)
clear;  % Erase all existing variables. Or clearvars if you want.
workspace;  % Make sure the workspace panel is showing.
format long g;
format compact;
fontSize = 20;
h0 = 20;	% Initial drop height in meters.

v = 0;		% Initial y velocity in m/sec.
g = 9.8;	% Gravitational acceleration in m/s^2;
t = 0;		% Initial time when dropped
dt = 0.01;	% Delta time in seconds.
rho = 0.75;	% Velocity reduction factor.  Velocity reduces this much after a bounce.
peakHeight = h0;	% Initial drop height in meters.
h = h0;		% Instantaneous height.
hstop = 0.01;	% Height at which if the peak height after a bounce is less than this, stop the simulation.
% Preallocate arrays for time and height.  Make them plenty large - we will crop to the final size later.
T = 0 : dt : 1000;
H = zeros(1, length(T));
% Setup a failsafe.  Don't do more than this number of iterations or else we might have an infinite loop.  This will prevent that.
maxIterations = 100000;		
loopCounter = 1;
while (peakHeight > hstop) && (loopCounter < maxIterations)
	% Compute new height.
	hNew = h + v * dt - 0.5 * g * dt ^ 2;
	% 	fprintf('After iteration %d, time %f, hmax = %f.\n', loopCounter, T(loopCounter), hNew);
	if(hNew<0)
		% Ball hit the ground.
		
		% Find index of last time h was 0
		lastBounceIndex = find(H(1 : loopCounter-1) == 0, 1, 'last');
		if isempty(lastBounceIndex)
			% If it hasn't bounced yet, start looking from the beginning.
			lastBounceIndex = 1;
		end
		% Compute the greatest height since the last bounce, or the initial release.
		[peakHeight, index] = max(H(lastBounceIndex : end)); % Record height
		% Find time when it was at that height.
		tMax = T(index + lastBounceIndex - 1);
		plot(tMax, peakHeight, 'b+', 'MarkerSize', 18, 'LineWidth', 2);
		hold on;
		fprintf('After iteration %d, time %f, hmax = %f.\n', loopCounter, tMax, peakHeight);
		
		% Reflect it up.  For example, if at this time,
		% the ball was going to be at -4 (with no ground in the way)
		% Now, after bouncing, it would be at +4 above the ground.
		h = 0; %abs(hNew);
		v = -v * rho;
		
	else
		% Ball is falling or rising.
		v = v - g*dt;
		h = hNew;
	end
	H(loopCounter) = h;
	loopCounter = loopCounter + 1;
end
% Crop times.
T = T(1 : loopCounter - 1);
H = H(1 : loopCounter - 1);
% Plot the trajectory.
plot(T, H, 'r.', 'MarkerSize', 5);
grid on;
xlabel('Time [seconds]', 'FontSize', fontSize)
ylabel('Ball Position Y [m]', 'FontSize', fontSize)
title('Bouncing Ball', 'FontSize', fontSize)
% Maximize window
g = gcf;
g.WindowState = 'maximized';
fprintf('Done running %s.m.\n', mfilename);
Related Question