% ### EXprojectileMOD.m ### 2020.01.28 C. Bergevin % [REF: ex.4.3.2 from Fowles & Cassidy 2005] % Purpose: Modified version of EXprojectile.m to show "optimal" angle for % projectile flights w/ drag is less than 45 degrees % ---- Notes % o v0= 143.2 mph ~ 64 m/s clear % ------------------------------------------ P.thetaA= linspace(36,42,21); % launch angle [degrees] P.g= 9.8; % grav. const. [m^2/s] {9.8} P.drag= 1; % boolean to incl. drag: 0=no drag, 1=drag {1} P.v0= 64; % launch velocity [m/s] {64} %P.theta= 39; % launch angle [degrees] {45} P.D= 0.073; % diameter of object [m] {0.073} P.m= 0.145; % mass of object [kg] {0.145} P.coord0= [0 0]; % initial [x z] coords [m] {[0 0]} P.tLim= [0 15]; % time limits of integration [s] P.tRez= 300; % # of (interp.) time points for integration interval {300?} % ------------------------------------------ % --- derived params. if (P.drag==0), P.gamma= 0; % determine assoc. const. from input params. else P.gamma= 0.15*P.D^2/P.m; end P.y0(1)= P.coord0(1); P.y0(3)= P.coord0(2); % initial horiz. and vert. positions % --- set up fig. plus color-coding scheme colormap(jet(numel(P.thetaA))) jetcustom = jet(numel(P.thetaA)); figure(1); clf; hold on; grid on; % --- for nn=1:numel(P.thetaA) P.theta= P.thetaA(nn); P.theta= pi*P.theta/180; % convert launch angle to rads P.y0(2)= P.v0*cos(P.theta); % initial horiz. velocity P.y0(4)= P.v0*sin(P.theta); % initial vert. velocity % --- use built-in solver ode45 to numerically integrate [t vals] = ode45('PROJECTILEfunction',linspace(P.tLim(1),P.tLim(2),P.tRez),P.y0,[],P); indxG= find(vals(:,3)<0,1); % find when object hits the ground % --- rename vars. (excluding those in the ground!) & plot x= vals(1:indxG,1); xdot= vals(1:indxG,2); z= vals(1:indxG,3); zdot= vals(1:indxG,4); plot(x,z,'-','LineWidth',1,'Color',jetcustom(nn,:)); end xlabel('x [m]'); ylabel('z [m]'); hC=colorbar; caxis([min(P.thetaA) max(P.thetaA)]); ylabel(hC, 'Launch angle [deg]') title('Projectile range w/ quadratic drag for multiple launch angles')