Finding the range of the robot.

The robot will throw a ball which is attached to the end of an arm of length 0.25 m rotating around the origin. We want to experiment with the rotation speed of the robot to figure out what the maximum range and height of the throw are. 

Our basic setup is that while on the arm, the ball follows the path

xArm[t_] := (1/4)*Cos[k*t];
yArm[t_] := (1/4)*Sin[k*t];

and once in the air, the ball follows the path (constant speed motion in x, constant acceleration motion in y):

g = -9.8;

xAir[t_] := a*t + b;
yAir[t_] := (g/2)*t^2 + c*t + d;

At the time of release (t0), the position and velocity of these two curves are supposed to agree. We can use this information to solve for a, b, c, and d.

absol = Solve[{xArm[t0] == xAir[t0],D[xArm[t],t]== D[xAir[t],t] /. {t -> t0}},{a,b}]

cdsol =Solve[{yArm[t0] == yAir[t0],D[yArm[t],t]== D[yAir[t],t] /. {t -> t0}},{c,d}]

We now use this information to write the path in the air very concretely in terms of k and t0 as xFlight and yFlight. These functions will describe the path of the ball as a function of t0 and k.

xFlight[t_] := xAir[t] /. absol[[1]];
yFlight[t_] := yAir[t] /. cdsol[[1]];

We can see that these formulae look pretty unpleasant:

xFlight[t]
yFlight[t]

But they actually look a lot better when you write them in terms of (t - t0):

Simplify[xFlight[t]]

The computer isn’t smart enough to do this on its own for y, but we can write it out ourselves:

yFlightNew[t_] := (g/2)*(t - t0)^2 + (1/4)*k*(t - t0)*Cos[k t0] + (1/4)*Sin[k t0];

yFlightNew[t]

We can check that this formula is the same as the old one by subtracting them and simplifying:

Simplify[yFlightNew[t] - yFlight[t]]

Now these equations are still fairly easy to work with using the computer! We want to now find the maximum height reached by the ball in terms of k and t0. This is a classic max-min problem: we want to find the maximum of yFlight on the interval [t0,infty). To do that, we need to differentiate yFlight and find the point where the derivative is zero.

yflightcriticalpt =t /. Solve[D[yFlight[t],t] == 0,t][[1]]

Now we have to be a little careful! This is the time at which the projectile reaches the maximum height IF this time is > t0. We can write down this function in Mathematica using an “If” statement:

yMax[T0_,K_] := (If[yflightcriticalpt >T0,yFlight[yflightcriticalpt],yFlight[T0]] /. {t0 -> T0,k->K}) /. {t0 -> T0,k->K};

yMax[3,10]

It is now time for some graphics, to make sure that all our formulas look like they are making sense:

ballPath[t_,T0_,K_] := If[t > T0,{xFlight[t],yFlight[t]},{xArm[t],yArm[t]}] /. {t0 -> T0,k -> K};

Launcher[r_,theta_] := {
Graphics[Rectangle[{-0.01,-1.5 r},{0.01,0.01}]],
ParametricPlot[{r*Cos[u],r*Sin[u]},{u,0,2*Pi}]
 };

We can now see the ball and the launcher and the maximum height on a single plot:

Manipulate[Show[ParametricPlot[ballPath[t,LaunchAngle/myK,myK],{t,0,maxT}],Launcher[0.22,T0],ParametricPlot[{t,-0.35},{t,-2,2}],Plot[yMax[LaunchAngle/myK,myK],{t,-2,2}],
PlotRange->{{-0.4,1.5},{-0.5,1}},AxesOrigin->{0,0}],{maxT,0.01,4},{myK,1,30},{LaunchAngle,0,2*Pi}]

Let’s get more abstract. Suppose we study the maximum value as a function of K alone. We’ll plot the maximum height for each t0 value between 0 and 2 Pi/K.

Manipulate[Plot[yMax[LaunchAngle/myK,myK],{LaunchAngle,0,2*Pi}],{{myK,20},1,25}]

This is a really interesting function! We can see that the maximum occurs close to zero across the entire range. Now what is the region of the plane that we can target? The simplest way to answer that question is to plot all of the trajectories for a given K.

Manipulate[Show[Table[ParametricPlot[ballPath[t,LaunchAngle/myK,myK],{t,0,10}],{LaunchAngle,0,2*Pi,2*Pi/40}],Launcher[0.22,T0],ParametricPlot[{t,-0.35},{t,-3,3}],
PlotRange->{{-2,2},{-0.5,2}},AxesOrigin->{0,0},ImageSize->600],{{myK,20},1,30}]

Now we get to the group work problems for this week:

1. Find the t0 and t values which make the ball land on the table (y = -0.35) at the location x = 1 when k = 20. This means that you have to solve the system of equations:

xFlight[t] = 1
yFlight[t] = -0.35

or 

(1/4) (Cos[k t0]+k (-t+t0) Sin[k t0]) = -1

-4.9` (((t-t0))^2)+(1/4) k (t-t0) Cos[k t0]+(1/4) Sin[k t0] = -0.35

for t0 and t, given that k = 20.

Solution:

First, we solve for the t value when y is equal to -0.35 in terms of t0 using the quadratic formula:

tSolns = Solve[{yFlight[t] == -0.35} /. {k -> 20},{t}]

Two solutions for t in terms of t0? Of course, this makes sense mathematically-- we’re solving a quadratic, after all. But it doesn’t make much sense physically. After all, when we throw the ball, there should be only one time where it touches the ground after we release it (the other time must be less than t, and so correspond to some solution to the y quadratic before the actual release of the ball!). 

To see this, we can plot both of these functions along with t0 as functions of t0. Here’s the plot of the first solution (t0 is a line on the plot). What should be the range for t0? Well, we’re interested in the range of angles 3*Pi/2 to 2*Pi. But we are covering k=20 radians of angle per second, so this is the range of times 3*Pi/40 to 2*Pi/20.

Plot[{t /. tSolns[[1]],t0},{t0,3*Pi/40,2*Pi/20},Filling->Axis]

That’s not going to work: we are ``striking the ground’’ at a time t (the lower curve) before we released the ball! Let’s try the other solution:

Plot[{t /. tSolns[[2]],t0},{t0,3*Pi/40,2*Pi/20},Filling->Axis]

Much better. The curve is above the line, so the ball is striking the ground after it’s thrown. This means we should go with this solution from here on out.

So we plug this solution into our xFlight function to get the function which describes where the ball strikes the floor in terms of t0.

xStrike[T0_] := xFlight[t] /. tSolns[[2]] /. {t0 -> T0,k->20};

xStrike[x] // InputForm

Now we can plot the strike position in terms of the release time and see where it intersects the line y = 1. We see that we should have two solutions-- one around 0.24 and another around 0.31.

Plot[{xStrike[t0],1},{t0,3*Pi/40,2*Pi/20},Filling->Axis]

So let’s try to solve the equations:

Solve[xStrike[t0] == 1.0,{t0}]

Hmm. No dice. Well, we can try a numerical method:

FindRoot[xStrike[t0] == 1.0,{t0,0.24}]

FindRoot[xStrike[t0] == 1.0,{t0,0.31}]

Wolfram|Alpha will do this with the command “find roots of (Cos[20*x] + 20*x*Sin[20*x])/4 - 0.0510204081632653*Sin[20*x]*
  (98.*x + 50.*Cos[20.*x] + 1.4142135623730951*
    Sqrt[968. + 625.*Cos[40.*x] + 245.*Sin[20.*x]]) - 1”

These solutions should do nicely for us: Let’s check on our original plot that they do really correspond to two trajectories which strike (1,-0.35).

Show[{ParametricPlot[ballPath[t,0.2401,20],{t,0,10}],ParametricPlot[ballPath[t,0.307105,20],{t,0,10}],Launcher[0.22,3],ParametricPlot[{t,-0.35},{t,-3,3}]},
PlotRange->{{-0.75,1.5},{-0.5,1.5}},AxesOrigin->{0,0},ImageSize->600]

2. Now suppose that our robot can launch the ball at this time t0 (accurate to within 1/1000 of a second). Estimate the resulting error in the x position of the ball when it hits the table. There are a few steps to this problem.

a. Write the x position of the ball when it hits the table in terms of t0 alone. You’ll need to use the yFlight equation

-4.9` (((t-t0))^2)+(1/4) k (t-t0) Cos[k t0]+(1/4) Sin[k t0] = -0.35

to figure out when (at what value of t) the ball hits the table (in terms of t0) and then plug that time into the xFlight equation 

xFlight[t] = (1/4) (Cos[k t0]+k (-t+t0) Sin[k t0])

in order to figure out where the ball hits the table. Your answer should be a function

xStrike[t0]

b. Now do an error estimate for this new function xStrike of t0, using t0 as the variable, the value of t0 you got in part 1 for the point (that is, the x0 value in the examples we did in class) and 1/1000 for the error (Δ t0). As we did before, you’ll have to differentiate xStrike[t0] and evaluate at the given value of t0.

c. What size basket can we hit with this accuracy? Suppose the error in the release time were 1/100 of a second because we bought cheaper parts. How would that affect our accuracy? Suppose the error was reduced to 1/10,000 of a second because we bought better parts? How would that affect our accuracy?

Solution:

We already know our xStrike[t0] function:

xStrike[t0]

So taking the derivative isn’t that hard:

xsD = D[xStrike[t0],t0]

For the first error estimate, we take t0 = 0.2401 and Δt0 = 1/1000. We use our derivative to compute (d/d(t0)) xStrike[t0] at t0 = 0.2401 and get

error1 = (xsD /. {t0 -> 0.2401}) * (1/1000)

That is, an error on the order of 0.07 m, or about 7 cm or 4 inches. Ouch! What about the second release time? Will we get better or worse error? We compute again with the same Δt0 and this time take t0 = 0.307105. We get

error2 = (xsD /. {t0 -> 0.307105}) * (1/1000)

Nope! Even worse! This time, we are off by about 10 cm or 6 inches. We plot these on the picture to see how bad it is:

Show[{ParametricPlot[ballPath[t,0.2401 - 1/1000,20],{t,0,10}],ParametricPlot[ballPath[t,0.2401 + 1/1000,20],{t,0,10}],ParametricPlot[ballPath[t,0.307105 - 1/1000,20],{t,0,10}],ParametricPlot[ballPath[t,0.307105 + 1/1000,20],{t,0,10}],Launcher[0.22,3],ParametricPlot[{t,-0.35},{t,-3,3}]},
PlotRange->{{-0.75,1.5},{-0.5,1.5}},AxesOrigin->{0,0},ImageSize->600]

Of course, improving our release time accuracy by a factor of 10 to an error of 1/10,000 reduces the error in strike position by a factor of 10, to about 1 cm. This is tolerable:

Show[{ParametricPlot[ballPath[t,0.2401 - 1/10000,20],{t,0,10}],ParametricPlot[ballPath[t,0.2401 + 1/10000,20],{t,0,10}],ParametricPlot[ballPath[t,0.307105 - 1/10000,20],{t,0,10}],ParametricPlot[ballPath[t,0.307105 + 1/10000,20],{t,0,10}],Launcher[0.22,3],ParametricPlot[{t,-0.35},{t,-3,3}],
ParametricPlot[{1,t},{t,-3,3}]},
PlotRange->{{-0.75,1.5},{-0.5,1.5}},AxesOrigin->{0,0},ImageSize->600]

And a release time accuracy of 1/100 second just gives comical results as it increases the error by a factor to 10 to about 1 m! At this point, we might as well just release the ball randomly!

Show[{ParametricPlot[ballPath[t,0.2401 - 1/100,20],{t,0,10}],ParametricPlot[ballPath[t,0.2401 + 1/100,20],{t,0,10}],ParametricPlot[ballPath[t,0.307105 - 1/100,20],{t,0,10}],ParametricPlot[ballPath[t,0.307105 + 1/100,20],{t,0,10}],Launcher[0.22,3],ParametricPlot[{t,-0.35},{t,-3,3}],
ParametricPlot[{1,t},{t,-3,3}]},
PlotRange->{{-0.75,1.5},{-0.5,1.5}},AxesOrigin->{0,0},ImageSize->600]

Extra Credit Problems.

I. I am fascinated by the question of the “envelope” curve in the picture

as a function of k. Is it really a parabola? If so, what’s the formula for the parabola in terms of k? 

a. For a given x value, the height of this function f(x) is the largest y of any trajectory of the ball (that is, it is the largest y value for that x as we vary t0 and keep k fixed). Looking at the xFlight and yFlight functions again, we see that the question here is to solve

xFlight[t] = x or (1/4) (Cos[k t0]+k (-t+t0) Sin[k t0]) = x

for t (in terms of t0, k, and x) and then plug that into yFlight[t] or -4.9` (((t-t0))^2)+(1/4) k (t-t0) Cos[k t0]+(1/4) Sin[k t0] to get a function yFlight[x] in terms of the variables t0 and k. 

Now our function is defined by 

f[x] = Max_t0 yFlight[x]

We can get this by taking the derivative of yFlight[x] with respect to t0 and finding the maximum value.

b. Once you get the f[x], we have to ask whether it is quadratic in x. I have no idea whether this is true, but it might take some fairly serious algebra in order to prove it. 

II. Suppose that in the classwork problem, there was error in the rotation speed k as well as in the time of release t0. Supposing we release the ball at exactly the correct time t0 (calculated in part 1) by the rotation speed k varies by 1% from the assumed value of k = 20. Estimate the resulting error in the x position when the ball strikes the table.