% Modeling the velocity of a skydiver from free-fall through deployment
%    of the parachute.  We take
%    v' = -32.2 + 0.00104 v^2 (before parachute)
%    v' = -32.2 - 2.01 v (after)
% and v(0) = 0.
%
% by Justin Dunmyre
% modified by Gavin LaRose, 8/21/12

% (c)2012 Regents of the University of Michigan
% This work is licensed under the
% Creative Commons Attribution-NonCommercial-ShareAlike 3.0 Unported License.
% To view a copy of this license, visit
% http://creativecommons.org/licenses/by-nc-sa/3.0/
% or send a letter to
% Creative Commons,
% 444 Castro Street, Suite 900,
% Mountain View, CA, 94041, USA.

function skydiver()

%-------------------------------------------------------------------------------
% configuration variables

% the time at which to open the parachute
    timeOpenParachute=15;
% the end time to which to graph the solution
    totalTime=30;

% the color of the solution for the velocity pre-deployment    
    closedColor='c';
% the color of the solution for the velocity post-deployment    
    openColor='g';

%-------------------------------------------------------------------------------
% the solutions to the differential equations

% the solution for the velocity pre-deployment
    function v=ClosedExactSol(t)
        v = 176*(1 - exp(.366*t))./(1 + exp(.366*t));
    end

% the solution for the velocity post-deployment
    function v=OpenExactSol(t)
       % Initial Velocity based on closed parachute solution, adjusted
       initVelocity=ClosedExactSol(timeOpenParachute);
       v = -16 + (initVelocity+16)*exp(-2.01*(t-timeOpenParachute));
    end

% the times at which we graph the velocities and corresponding
%    values of the velocities
    tClosed=linspace(0,timeOpenParachute,100);
    tOpen=linspace(timeOpenParachute,totalTime,100);
    vClosed=ClosedExactSol(tClosed);
    vOpen=OpenExactSol(tOpen);
    
% graph these
    scrsz = get(0,'ScreenSize'); % Get the screen size
    figh  = figure('Position',[1 scrsz(4)/4 scrsz(3)/4 scrsz(4)/2]);
    axes('Parent', figh, 'Color', [0.6 0.6 0.6]);
    hold on;
    axis([0 totalTime min(vOpen)-20 0]);
    xlabel('t','FontSize',20);
    ylabel('v','FontSize',20);
    title('Skydiver velocity vs t','FontSize',20);
    
    eqntext = annotation( figh, 'textbox', [0.21 0.66 0.33 0.092], ...
			  'String', ...
			  {'v\prime=-32.2 + 0.00104 v^2', ...
			   '(before parachute)'}, ...
			  'FontSize', 20, 'Color',closedColor, ...
			  'FitBoxToText', 'on' );
    eqntext = annotation( figh, 'textbox', [0.57 0.32 0.28 0.092], ...
			  'String', ...
			  {'v\prime=-32.2 - 2.01 v', ...
			   '(after parachute)'}, ...
			  'FontSize', 20, 'Color', openColor,...
			  'FitBoxToText','on' );
    
    plot(tClosed, vClosed, [closedColor '-'], 'LineWidth',2);
    plot(tOpen, vOpen, [openColor '-'], 'LineWidth',2);
    
end

% all done
%-------------------------------------------------------------------------------
