Collisions
From ComputingForScientists
Contents 
1. References
 http://introcs.cs.princeton.edu/java/assignments/collisions.html
 http://www.opensourcephysics.org/items/detail.cfm?ID=7375 (page 284)
 https://www.geogebra.org/m/W62QndVz
 http://physics.weber.edu/schroeder/md/
 http://iopscience.iop.org/article/10.1088/01430807/6/3/005/meta
 http://aapt.scitation.org/doi/abs/10.1119/1.1482060
2. Problems
2.1. 1D collision with ground I
A spherical object with mass m with radius r is released from a height h_{o} above the ground with an initial vertical velocity of v_{oy}. When the object strikes the ground, it has an elastic collision.
 Choose values of h_{o}, v_{oy}, and r and draw an approximate sketch of what the curves for the height and velocity of the object as a function of time should look like. Label any key points in terms of the parameters that you chose.
 Write down the equations for the exact solutions for y(t) and v_{y}(t) before the object strikes the ground and then write a program to plot these equations.
 Write down equations for dy / dt and dv_{y} / dt.
 The Forward Euler method requires a parameter Δt to be used. Instead of using trialanderror, use the physics of the problem to come up with a Δt that is reasonable to use.
 Use the Forward Euler method for your equations for dy / dt and dv_{y} / dt to estimate y(t) and v_{y}(t) from the time of release until the object strikes the ground.
 Plot the exact and Forward Euler solutions for y(t) and v_{y}(t) from the time the object is released to the time the object strikes the ground. One plot should show the exact and Forward Euler y(t) and the other plot should show the exact and Forward Euler v_{y}(t).
Your code will be tested using the value of Δt that you chose and then several other values. Make sure that the solution makes sense at the reflection point. A common error is that the exact solution is not correct for arbitrary values of Δt (exact solution should not depend on this parameter  it is only relevant for Euler method solution!).
Answer 

Write down the equations for the exact solutions for y(t) and v_{y}(t) before the object strikes the ground and then write a program to plot these equations. For constant acceleration, the following kinematic equations apply y(t) = y_{o} + v_{oy}t + a_{y} / 2t^{2} v_{y}(t) = v_{oy} + a_{y}t where for the given problem, a_{y} = − g and y_{o} = h. These are exact solutions for the position and velocity for a particle the experiences constant acceleration. They are derived from the equations dv_{y} / dt = − g and dy / dt = v_{y}; see an introductory physics textbook for the derivation. To get the exact value of time for when the object strikes the ground, use the kinematic equation y(t) = y_{o} + v_{oy}t + (a_{y} / 2)t^{2} and set y(t) = r and solve for t. Assuming v_{oy} = 0, the result is
clear r = 1; % m ho = 10; % m voy = 0; % m/s g = 9.8; % m/s^2 % thit only correct when voy = 0. Otherwise need to solve quadratic equation. thit = sqrt(2*(hor)/g); %t = [0:0.1:thit]; % Don't use this. Use linspace to ensure last value in t is thit t = linspace(0,thit,100); y = ho + voy*t + (g/2)*t.^2; v = voy + (g)*t; % Or, using a loop instead of vectorization for i = 1:length(t) y(i) = ho + voy*t(i) + (g/2)*t(i)^2; v(i) = voy + (g)*t(i); end figure(1);clf; grid on; hold on; plot(t,y,'k','LineWidth',2); plot(t,v,'g','LineWidth',2); xlabel('time [s]'); legend('y [m]','v [m/s]'); Write down equations for dy / dt and dv_{y} / dt. By definition, dy / dt = v_{y} From Newton's second law, dv_{y} / dt = − g The Forward Euler method requires a parameter Δt to be used. Instead of using trialanderror, use the physics of the problem to come up with a Δt that is reasonable to use. We want many steps to be taken using the Forward Euler method as the object falls, so we could use the time it takes the object to hit the ground and divide by 100, so that
Thinking more about the problem, you may realize that we want to have small steps taken just before the object hits the ground. In that case, we want Δt to be a fraction of the time it takes the object to travel a distance r. Just before the object hits the ground, it has a velocity of and so it takes for the object to travel a distance r, and so a better starting point for Δt is
As usual, we will try larger and smaller values of Δt and consider how the numerical solution changes and then use the a value of Δt for which the solution is both accurate and can be computed in a reasonable amount of time. Use the Forward Euler method for your equations for dy / dt and dv_{y} / dt to estimate y(t) and v_{y}(t) before the object strikes the ground. The equations dv_{y} / dt = − g and dy / dt = v_{y} do not tell us y(t) and v_{y}(t) directly. Using an analytic method (continuous integration), one can get the exact solutions, which correspond to the kinematic equations given above. If one does not know the analytic method, one can make the (Forward Euler) approximations
and then get an approximate (or numerical) solution by iterating Δy / Δt = v_{y} Δv_{y} / Δt = − g In code, the first two iterations are clear r = 1; % m ho = 10; % m voy = 0; % m/s g = 9.8 % m/s^2 % Initial conditions y(1) = ho; vy(1) = voy; % Forward Euler time step dt = (r/sqrt(2*g*y(1)))/100; y(2) = y(1) + vy(1)*dt vy(2) = vy(1) + (g)*dt y(3) = y(2) + vy(2)*dt vy(3) = vy(2) + (g)*dt To get value of y(t) and v_{y}(t) to plot for the Forward Euler (approximate) solution, modify the code from the previous question to compute more time steps. A The full solution is: % dx/dt = vx (by definition) % m*dvx/dt = mg (Newton's second law) clear; ho = 10; % m vyo = 0; % m/s (wall is at x = 0, puck starts at x = L) r = 1; % m g = 9.8; % m/s^2 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Exact % Before collision, exact (assuming vyo = 0). thit = sqrt(2*(hor)/g); % Above quation only correct if vyo is not zero, otherwise need % to solve quadratic equation for t: hor = vyo*t  g*t^2/2. te = linspace(0,thit,100); ye = ho + vyo*te  0.5*g*te.^2; vye = vyo  g*te; %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % Forward Euler % m*dvx/dt = mg => vx(i+1) = vx(i)  dt*g % dx/dt = v => x(i+1) = x(i) + dt*vx(i) y(1) = ho; vy(1) = vyo; t(1) = 0; dt = r/sqrt(2*g*ho)/10; i = 1; while (y(i) > r) fprintf('x(%d) = %.2f; vy(%d) = %.2f\n',i,y(i),i,vy(i)); vy(i+1) = vy(i)  g*dt; y(i+1) = y(i) + dt*vy(i); t(i+1) = t(i) + dt; i = i+1; end figure(1);clf plot(te,ye,'k.','MarkerSize',15); hold on;grid on; plot(t,y,'g.','MarkerSize',10); plot([t(1),t(end)],[r,r],'k'); % Plot ground line. legend(' Exact',' Forward Euler','x = r line'); xlabel('t [s]'); ylabel('x [m/s]'); figure(2);clf plot(te,vye,'k.','MarkerSize',15); hold on;grid on; plot(t,vy,'g.','MarkerSize',10); xlabel('t [s]'); ylabel('v_x [m/s]'); legend(' Exact',' Forward Euler'); 
Student Question 

In this code, I made it so the mass lost energy with each collision. Everything looks fine for a few bounces, but then the solution flatlines. Why? Hint: Zoom in on the plot where the trajectory becomes flat and look at what is happening. r = 1; a = 9.8; y(1) = 100; v(i) = 3; g = 9.8 dt = r/(100*sqrt(2*g*y(1))); t(1) = 0; for i = 1:200000 v(i+1) = v(i) +dt*a; y(i+1) = y(i) +v(i)*dt; t(i+1) = t(i) +dt; if (y(i+1) <= r) v(i+1) = 1*v(i)*0.8; end end plot(t,y) 
2.2. 1D collision with ground II
Build on your answer to the previous problem and
 Write down the exact equations for y(t) and v_{y}(t) after the object strikes the ground for the first time and before it strikes the ground for the second time.
 Use the Forward Euler method to estimate y(t) and v_{y}(t) from the time of release until it strikes the ground for the second time.
 Plot the exact and Forward Euler solutions for y(t) and v_{y}(t) from the time the object is released to the time that it hits the ground for the second time. One plot should show the exact and Forward Euler y(t) and the other plot should show the exact and Forward Euler v_{y}(t).
Hint 

For the exact equations, you'll use the kinematic equations with different initial conditions for y and v_{y} that correspond to their values just after the first bounce. For the Forward Euler method, you will need to make two changes:

Partial Answer 

Modify the Forward Euler loop from the previous answer to be i = 1; bounce = 0; while (bounce < 3) fprintf('x(%d) = %.2f; vy(%d) = %.2f, bounces = %d\n',... i,y(i),i,vy(i),bounce); vy(i+1) = vy(i)  g*dt; y(i+1) = y(i) + dt*vy(i); t(i+1) = t(i) + dt; if (y(i+1) < r) vy(i+1) = vy(i+1); bounce = bounce + 1; end i = i+1; end 
2.3. 1D collision with wall
This problem is similar to the previous problem except that the object does not experience the force of gravity.
An hockey puck of mass m, radius r, and velocity v_{ox} is moving along the xaxis towards a wall with which it has an elastic collision. Initially it is a distance L from the wall. There is no friction between the puck and the ice.
 Choose values of L, r, m, and v_{ox} and draw an approximate sketch of what the curves for the position and velocity of the puck as a function of time should look like from its time of release until it returns to its starting point. Label any key points.
 Write down the exact equations for y(t) and v_{y}(t) before the puck strikes the wall.
 Write down equations for dx / dt and dv_{x} / dt of the puck.
 Use the Forward Euler method and the equations you wrote down for dx / dt and dv_{x} / dt to estimate x(t) and v_{x}(t) before the puck strikes the wall.
 Plot the exact and Forward Euler solutions for x(t) and v_{x}(t) before the puck strikes the wall.
 Write down the exact equations for x(t) and v_{x}(t) after the puck strikes the wall.
 Use the Forward Euler method to estimate x(t) and v_{x}(t) after the puck strikes the wall.
 Plot the exact and Forward Euler solutions for x(t) and v_{x}(t) after the puck strikes the wall.
Basic Solution 

% An hockey puck of mass m, radius r, and velocity vox is moving along the % xaxis towards a wall. Initially it is a distance L from the wall. % There is no friction between the puck and the ice. % The solution here is written in a basic way to illustrate the key points. % After writing code in a simple way, one should look at what is done % and look for generalizations and simplifications. % dx/dt = vx (by definition) % m*dvx/dt = 0 (Newton's second law) clear; L = 10; % m vxo = 1; % m/s (wall is at x = 0, xo = L) r = 1; % m % Forward Euler % vx(i+1) = vx(i) + dt*0 => vx(i+1) = vx(i) % x(i+1) = x(i) + dt*vx(i) x(1) = L; vx(1) = vxo; t(1) = 0; % dt is 1/100 of time it takes puck to move its radius. dt = r/abs(vxo)/100; %%%%% Before collision i = 1; % Use while loop to iterate until puck center is at r or % past r. while (x(i) >= r) vx(i+1) = vx(i); x(i+1) = x(i) + dt*vx(i); t(i+1) = t(i) + dt; i = i+1; end % Exact solution % dvx/dt = 0 => vx = const % constant must be vxo % % dx/dt = vx = vxo => dx = vx*dt => % xxo = vxo*(tto). Let to = 0 and xo = L. % x = xo + vxo*t = L + vxo*t % Hit time is computed from exact equation. thit = abs((Lr)/vxo); te = linspace(0,thit,100); xe = L + vxo*te; % Make ve an array with same number of elements as x. % Could also do vxe = vxo + 0*te vxe = vxo*ones(1,length(te)); figure(1);clf plot(te,xe,'k','LineWidth',5); hold on;grid on; plot(t,x,'g.','MarkerSize',10); legend(' Exact',' Forward Euler'); xlabel('t [s]'); ylabel('x [m/s]'); figure(2);clf plot(te,vxe,'k','LineWidth',5); hold on;grid on; plot(t,vx,'g.','MarkerSize',10); xlabel('t [s]'); ylabel('v_x [m/s]'); legend(' Exact',' Forward Euler'); %%%%% After collision i = 1; % Use while loop to iterate until puck center is at r or % past r. % Note that the first time value after collision will be the last time value % before collision. When we combine arrays, there will be a duplicate. t2(1) = t(end); vx2(1) = vx(end); % Implement the collision x2(1) = x(end); while (x2(i) <= L) vx2(i+1) = vx2(i); x2(i+1) = x2(i)+dt*vx2(i); t2(i+1) = t2(i) + dt; i = i+1; end vxo2 = vxo; xo2 = r; te2 = te(end) + linspace(0,thit,100); % Reuse t array for computing exact solution. xe2 = xo2 + vxo2*(te2te(end)); % Make ve an array with same number of elements as x. vxe2 = vxo2*ones(1,length(te2)); %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% figure(3);clf plot([te,te2],[xe,xe2],'k','LineWidth',5); hold on;grid on; plot([t,t2],[x,x2],'g.','MarkerSize',5); legend(' Exact',' Forward Euler'); xlabel('t [s]'); ylabel('x [m/s]'); figure(4);clf plot([te,te2],[vxe,vxe2],'k','LineWidth',5); hold on;grid on; plot([t,t2],[vx,vx2],'g.','MarkerSize',5); xlabel('t [s]'); ylabel('v_x [m/s]'); legend(' Exact',' Forward Euler'); 
Improved Solution 

clear; L = 10; % m vxo = 1; % m/s (wall is at x = 0, puck starts at x = L) r = 1; % m % Forward Euler % m*dvx/dt = 0 => vx(i+1) = vx(i) + dt*0 => vx(i+1) = vx(i) % dx/dt = v => x(i+1) = x(i) + dt*vx(i) x(1) = L; vx(1) = vxo; t(1) = 0; % dt is 1/10 of time it takes puck to move its radius. dt = r/abs(vxo)/10; i = 1; while (x(i) <= L) % Run between start and when it returns near start pos. fprintf('x(%d) = %.2f; vx(%d) = %.2f\n',i,x(i),i,vx(i)); vx(i+1) = vx(i); x(i+1) = x(i) + dt*vx(i); t(i+1) = t(i) + dt; if (x(i+1) <= r) vx(i+1) = vx(i); % Implement the collision end i = i+1; end % Exact solution % Hit time is computed from exact equation. thit = abs((Lr)/vxo); % Before collision te1 = linspace(0,thit,100); xe1 = L + vxo*te1; vxe1 = vxo*ones(1,length(te1)); % After collision te2 = te1(end) + linspace(0,thit,100); xe2 = r + (vxo)*(te2te1(end)); vxe2 = (vxo)*ones(1,length(te2)); te = [te1,te2]; xe = [xe1,xe2]; vxe = [vxe1,vxe2]; % Remove duplicate values due to how calculation was done. [te,I] = unique(te); xe = xe(I); vxe = vxe(I); % A clearer way to compute exact solution  use time values from % Forward Euler (xc = x_clearer). This will make xc and x have % same time values and same length, thus allowing comparison by doing % xxc, for example. A drawback is that we may not have a time % value that corresponds to the exact time of the collision (if dt was % chosen so that thit does not appear in array t). To see this, divide % dt by 9.1 instead of 10.0 and zoom in on the collision point on the x(t) % plot. for i = 1:length(t) if (t(i) < thit) xc(i) = L + vxo*t(i); else xc(i) = r + (vxo)*(t(i)thit); end end % Use the zoom tool to see on plot how the Forward Euler solution shows % the puck having a position that is < r. This is still a problem that % needs to be addressed. figure(1);clf plot(te,xe,'k.','MarkerSize',15); hold on;grid on; plot(t,x,'g.','MarkerSize',10); plot(t,xc,'r.','MarkerSize',5); legend(' Exact',' Forward Euler','Alternative Exact'); xlabel('t [s]'); ylabel('x [m/s]'); figure(2);clf plot(te,vxe,'k.','MarkerSize',15); hold on;grid on; plot(t,vx,'g.','MarkerSize',10); xlabel('t [s]'); ylabel('v_x [m/s]'); legend(' Exact',' Forward Euler'); 
2.4. 1D Collision with walls
You are in a room with a ceiling height of 2 meters. You throw a rubber ball with a radius of 0.1 meters straight downward with a velocity of 10 m/s. The ball is released from your hand at a height of 0.5 meters. The ball has elastic collisions with the floor and ceiling.
Use the Forward Euler method to compute the height and velocity of the ball as a function of time. The last time plotted should be at about the time the ball hits the floor for the third time.
When your program is executed, two figures should appear. The first should show y(t) and the second should show v(t).
Answer 

clear r = 0.1; % m vo = 10.0; % m/s ho = 0.5; % m h = 2; g = 9.8; % Equation of motion mdv/dt = mg % Forward Euler: v(i+1) = v(i)  g*Dt; % Definition of velocity dy/dt = v % Forward Euler: y(i+1) = y(i) + v(i)*Dt; t(1) = 0; v(1) = vo; y(1) = ho; dt = abs(r/vo)/10; for i = 1:1000 v(i+1) = v(i)  g*dt; y(i+1) = y(i) + v(i)*dt; t(i+1) = t(i) + dt; if (y(i+1) < r) v(i+1) = v(i); end if (y(i+1) > h) v(i+1) = v(i); end end figure(1);clf plot(t,y); hold on; 
2.5. 1D collision between two masses
Repeat problem 2. but replace the wall with a puck of mass M and radius R that is initially at rest. Select your own values for the parameters.
2.6. 2D collision with ground
Repeat problem 1. except give the object an initial nonzero horizontal velocity of v_{ox}
2.7. 2D collision with wall
An hockey puck of mass m, radius r, and velocity v_{o} is moving towards a wall at an angle of 45^{o}. Initially it is a distance L from the wall. There is no friction between the puck and the ice and the collision is elastic. Compute and plot the exact and Forward Euler solutions for x, y, v_{x}, and v_{y} of the the puck from the initial time until the puck is a distance L from the wall.
2.8. 2D collision between two masses
An hockey puck of mass m, radius r, and velocity v_{o} is moving towards the origin at an angle of 45^{o}. A second stationary puck of mass M and radius R is at the origin. Initially the moving puck is a distance L from the origin. There is no friction between the puck and the ice and the collision is elastic. Compute and plot the exact and Forward Euler solutions for x, y, v_{x}, and v_{y} of both pucks from the initial time until the pucks are separated by a distance L.