## Velocity and Acceleration Vectors in Matlab

In this activity we will introduce Matlab's **quiver** command, which will assist in creating visualizations that help develop a deeper understanding of velocity and acceleration in the plane and three-space.

### Motion on a Line

Recall from your first calculus course (differential calculus) the following idea. Let `x` represent the position of a particle constrained to move along a line (see Figure 1). Because the position of the particle changes with time, `x` is a function of time `t`; i.e., `x=x(t)`.

Now, recall two facts from differential calculus:

- The first derivative of the position with respect to time gives the velocity of the particle on the line as a function of time; i.e.,
`v(t)=x'(t)`.

- The second derivative of the position with respect to time gives the acceleration of the particle on the line as a function of time; i.e.,
`a(t)=x''(t)`.

We will now extend these ideas of velocity and acceleration into the plane.

### Motion in the Plane

Suppose that a particle moves in the plane along the path shown in Figure 2.

**Position: **The *position* of the particle is determined by the *position vector* `vec r(t) = langle x(t), y(t) rangle`. That is, the vector `vec r` gives us the position of the particle as a function of time. The vector `vec r(t)` changes with time, both in magnitude and direction, and as it changes, the tip of the position vector traces out the path of the particle shown in Figure 2.

**Velocity: **As in the case of motion along a line, the first derivative of the position with respect to time yields the velocity. However, in the case of motion in the plane, the position is a vector valued function of time. To take the derivative, we simply take the derivative of each component of the position vector (reasons for this are discussed in your text).

`vec v(t)=vec r'(t)=langle x'(t), y'(t) rangle`

The result should be a vector that determines the *instantaneous* velocity at a particular moment in time. Consequently, the velocity vector should point in the direction of motion at that instant. Therefore, the velocity vector should be *tangent* to the path of the particle at that instant in time, as shown in Figure 3.

**Acceleration: **As in the case of motion along the line, the second derivative of the position with respect to time yields the acceleration. However, in the case of motion in the plane, the position is a vector valued function of time. To take the second derivative, we simply take the second derivative of each component of the position vector (reasons for this are discussed in your text).

`vec a=vec r''(t)=langle x''(t), y''(t) rangle`

Newton's second law states that

`vec F = m vec a.`

This means that the force vector `vec F` is a scalar multiple of the acceleration vector `vec a`, which means that the vectors `vec F` and `vec a` must be parallel. Furthermore, because the mass `m` is a positive scalar quantity, the vectors `vec F` and `vec a` must point in the same direction.

Newton's first law states that "an object in motion tends to stay in motion with the same speed and direction unless acted upon by an external force." Note how the path of the particle bends downward. This cannot happen unless an external force is pushing the particle in that direction. Therefore, the acceleration of the particle must point in that same direction. Thus, the depiction of the acceleration vector in Figure 4 seems reasonable.

An example should reinforce these concepts.

A particle moves on a path determined by the vector-valued function

`vec r(t)=langle 2 cos(t), 3 sin(t) rangle.`

Sketch the path taken by the particle and append velocity and acceleration vectors at several instances of time.

**Solution: **We begin by sketching the path determined by the position vector

`vec r(t)=langle x(t), y(t) rangle = langle 2 cos(t), 3 sin(t) rangle`.

We first initialize some time-values, then compute the position `x(t)=2 cos(t)`, `y=3 sin(t)` at each of these times.

t=linspace(0,2*pi,200); x=2*cos(t); y=3*sin(t);

Create the plot, add a grid, and equalize the axes.

plot(x,y) grid on axis equal

Add appropriate annotations.

xlabel('x-axis') ylabel('y-axis') title('The graph of x = 2 cos(t), y = 3 sin(t)')

The above sequence of commands produces the path of the particle shown in Figure 5.

Now that we have constructed the path, let's see how we can add velocity vectors to our path.

**Adding Velocity Vectors: **The velocity is the first derivative of the position, so:

$\begin{eqnarray} \vec v(t) &=& \frac{d}{dt} \vec r(t)\\ &=& \frac{d}{dt} \langle 2 \cos(t), 3 \sin(t) \rangle\\ &=& \langle -2 \sin(t), 3 \cos(t) \rangle \end{eqnarray}$

Matlab's **quiver** command is used to draw vectors. We first ask for information with the command **help quiver**, part of whose output follows.

QUIVER(X,Y,U,V) plots velocity vectors as arrows with components (u,v) at the points (x,y). QUIVER(U,V,S) or QUIVER(X,Y,U,V,S) automatically scales the arrows to fit within the grid and then stretches them by S. Use S=0 to plot the arrows without the automatic scaling.

We will use the second form of syntax, namely **QUIVER(X,Y,U,V,S)**. A quick explanation of the arguments is in order.

**X**and**Y**are passed the coordinates at which you will plot the tails of your vectors.**U**and**V**are passed the components of the vectors you wish to plot.**S**is a scaling factor. When we sketch*vector fields*in later activities, we will take advantage of this scaling factor. In this activity, we will set**S = 0**, which turns off all scaling of the vectors.

We will now draw a second plot and add velocity vectors at selected times. We begin as above and resketch the ellipse with the following commands.

t=linspace(0,2*pi,200); x=2*cos(t); y=3*sin(t); plot(x,y)

Next, we'll select times at which to plot our velocity vectors and calculate the positions at which the tails of the vectors will be applied.

tt=0:pi/4:2*pi; xx=2*cos(tt); yy=3*sin(tt);

As we saw above, `vec v(t)=langle -2 sin(t), 3 cos(t) rangle.` We will use this formula to calculate the components of the velocity vectors.

vx=-2*sin(tt); vy=3*cos(tt);

Hold the plot, then add the velocity vectors with the **quiver** command. Note that we set the scaling **S = 0**, which turns off all scaling of the vectors drawn by the **quiver** command.

hold on quiver(xx,yy,vx,vy,0)

Add the usual annotations.

xlabel('x-axis') ylabel('y-axis') title('The graph of x = 2 cos(t), y = 3 sin(t)')

The above sequence of commands will produce the plot shown in Figure 6.

In Figure 6, note that some vectors are longer than others; i.e., the magnitude of the vectors (the speed of the particle) changes with position. The speed seems to be at a maximum when the particle crosses the `x`-axis and at a minimum when the particle crosses the `y`-axis.

Let's investigate how we can add acceleration vectors to our path.

**Adding Acceleration Vectors: **The acceleration is the second derivative of position, or the first derivative of the velocity, so:

$\begin{eqnarray} \vec a(t) &=& \frac{d^2}{dt^2} \vec r(t)\\ &=& \frac{d}{dt} \vec v(t)\\ &=& \frac{d}{dt} \langle -2 \sin(t), 3 \cos(t) \rangle\\ &=& \langle -2 \cos(t), -3 \sin(t) \rangle \end{eqnarray}$

We will now draw a third plot which includes both velocity and acceleration vectors at selected times. We repeat our work above in sketching the path and velocity vectors.

t=linspace(0,2*pi,200); x=2*cos(t); y=3*sin(t); plot(x,y) tt=0:pi/4:2*pi; xx=2*cos(tt); yy=3*sin(tt); vx=-2*sin(tt); vy=3*cos(tt); hold on quiver(xx,yy,vx,vy,0)

Next, we calculate the components of the acceleration vectors

ax=-2*cos(tt); ay=-3*sin(tt);

We already "held" the plot, so we need only add the acceleration vectors.

quiver(xx,yy,ax,ay,0)

In this plot, we'll get a little fancy, adding points at the positions where the tails of our velocity and acceleration vectors are applied.

line(xx,yy,... 'linestyle','none',... 'marker','.',... 'markersize',6)

Some explanatory comments are in order:

- The
**line**command is used to draw a line using the data in the vectors**xx**and**yy**. - Note the three dots (ellipsis) at the end of some lines. This is Matlab's
*line continuation*character, which allows the user to split up very long lines over several lines. - Next come several "key-value" pairs. The first sets the key 'linestyle' to a value of 'none'. This turns off connecting each pair of consecutive points with a line segment.
- The next key 'marker' is set to the value '.', which results in a dot being used to plot each point.
- Finally, the key 'markersize' is set to a value of 6, which sets the size of the marker (the dot) to 6 points.

Add the usual annotations.

xlabel('x-axis') ylabel('y-axis') title('The graph of x = 2 cos(t), y = 3 sin(t)')

The above sequence of commands will produce the plot shown in Figure 7.

Note how all of the acceleration vectors point inward towards the origin. This example should bring to mind the physical situation of a planet orbiting the sun. The force of gravitation always points toward the sun, so the planet is always accelerating toward the sun.

### Motion in Space

The technique shown in our planar examples works equally well in three-space. Let's just jump into an example.

The position of a particle is given by the vector:

`vec r(t)=langle cos(t), sin(t), t//(8 pi) rangle`

**Solution: **It is a simple matter to sketch the path in three-space.

t=linspace(0,8*pi,200); x=cos(t); y=sin(t); z=t/(8*pi); plot3(x,y,z)

Add a grid and box for depth, equalize the axes, then orient the view to the standard orientation used in class.

grid on box on axis equal view([130,30])

The above sequence of commands will produce the helix shown in Figure 8.

Select some times at which to plot the velocity and acceleration vectors, then compute the position of the tails of the vectors at each time element.

tt=0:pi/2:8*pi; xx=cos(tt); yy=sin(tt); zz=tt/(8*pi);

**Adding Velocity Vectors: **The velocity is the first derivative of position, so:

$\begin{eqnarray} \vec v(t) &=&\frac{d}{dt} \vec r(t)\\ &=&\frac{d}{dt} \langle \cos(t), \sin(t), t/(8\pi) \rangle\\ &=& \langle -\sin(t), \cos(t), 1/(8\pi)\rangle \end{eqnarray}$

Use this result to compute the components of the velocity vectors at each time element.

vx=-sin(tt); vy=cos(tt); vz=ones(size(tt))/(8*pi);

The computation of **vz** requires a bit of explanation. We need each entry of **vz** to equal `1//(8 pi)`. So we first construct a vector of all ones that is the same size as the vector **tt**, then divide each entry by `8 pi`.

Hold the plot, then plot the vectors with the **quiver3** command.

hold on quiver3(xx,yy,zz,vx,vy,vz,0)

The resulting velocity vectors are shown in Figure 9.

**Adding Acceleration Vectors: **The acceleration is the second derivative of the position, or equivalently, the first derivative of the velocity.

$\begin{eqnarray} \vec a(t) &=&\frac{d^2}{dt^2} \vec r(t)\\ &=&\frac{d}{dt} \vec v(t)\\ &=&\frac{d}{dt} \langle -\sin(t), \cos(t), 1/(8\pi) \rangle\\ &=& \langle -\cos(t), -\sin(t), 0\rangle \end{eqnarray}$

Compute the components of the acceleration vectors.

ax=-cos(tt); ay=-sin(tt); az=zeros(size(tt));

The command **az=zeros(size(tt))** creates a vector having he same size as the vector **tt** and fills it with zeros. Add the acceleration vectors to the plot.

quiver3(xx,yy,zz,ax,ay,az,0);

Add some appropriate annotations.

xlabel('x-axis') ylabel('y-axis') zlabel('z-axis') title('The graph of x = cos(t), y = sin(t), z=t/8\pi')

The result of the above sequence of commands is captured in Figure 10.

### Constant Speed Curves

In Figure 10, note that all of the acceleration vectors are orthogonal (perpendicular) to the velocity vectors. This was not the case in the elliptical example of Figure 7. There is a reason why the acceleration vectors are orthogonal to the velocity vectors in the current example, the reason being that the position vector `vec r(t)=langle cos(t), sin(t), t//(8 pi) rangle` traces what is know as *constant speed curve*.

Because `vec r(t)=langle cos(t), sin(t), t//(8 pi) rangle`, we saw that the velocity is given by

`vec v(t)=langle -sin(t), cos(t), 1//(8 pi) rangle`.

Speed is the magnitude of velocity, so the speed of the particle traveling on the helix in Figure 10 is

`|vec v(t)|=sqrt(sin^2(t)+cos^2(t)+(1//(8 pi))^2)=sqrt(1+1/(64 pi^2))`.

Thus, the speed is constant. Note that the velocity is not constant, because the direction of the velocity vectors in Figure 10 is always changing. But the magnitude of the vectors does not change; i.e., the speed is constant.

Now, what happens when the speed is constant; i.e., what happens when `|vec v|=c`, where `c` is a constant? Well, in this case, the square of the speed will be `|vec v|^2=c^2`, so the square of the speed is also constant. Recall that `|vec v|^2=vec v cdot vec v`, so we can write

`vec v cdot vec v = c^2`.

Differentiate both sides of this last result with respect to `t`. Remember that `vec v=vec v(t)`; i.e., `vec v` is a function of time. Also, remember that `c` is constant. Note how we use the product rule for dot products in the following calculation.

$\begin{eqnarray} \frac{d}{dt}(\vec v \cdot \vec v)&=&\frac{d}{dt}(c^2)\\ \vec v\cdot\frac{d\vec v}{dt}+\frac{d\vec v}{dt}\cdot \vec v &=& 0\\ 2\vec v\cdot\frac{d\vec v}{dt}&=&0 \end{eqnarray}$

Dividing both sides of the last result by 2 yields the result we need.

`vec v cdot frac(d vec v)(dt)=0.`

Because the derivative of the velocity is acceleration,

`vec v cdot vec a=0.`

That is, in the case of constant speed curves, the velocity and the acceleration vectors must be orthogonal to each other at each point of the curve. This is evident in the case of the constant speed curve shown in Figure 10.

### Matlab Files

Although the following file features advanced use of Matlab, we include it here for those interested in discovering how we generated the images for this activity. You can download the Matlab file at the following link. Download the file to a directory or folder on your system.

The file **velacc.m** is designed to be run in "cell mode." Open the file **velacc.m** in the Matlab editor, then enable cell mode from the **Cell Menu**. After that, use the entries on the **Cell Menu** or the icons on the toolbar to execute the code in the cells provided in the file. There are options for executing both single and multiple cells. After executing a cell, examine the contents of your folder and note that a PNG file was generated by executing the cell.

### Exercises

In each of the following exercises, perform each of the following tasks:

- Sketch the path traced by the tip of the given position vector `vec r(t)`.
- Provide pencil and paper calculations to determine the velocity and acceleration vectors.
- Plot velocity vectors at several positions along the path.
- Plot acceleration vectors at the same positions used in plotting velocity vectors along the path.
- Publish your results to HTML print the results. Hand in the printed results in class.

- `vec r(t)= langle t^2,t^3 rangle` in the plane over the time interval `-1 le t le 1`.
- `vec r(t)= langle 2 cos(t), 2 sin(t) rangle` in the plane over the time interval `0 le t le 2 pi`.
- `vec r(t)= langle 2 cos(2t), 2 sin(2t), t//(2 pi) rangle` in three-space over the time interval `0 le t le 2 pi`.
- `vec r(t)= langle e^(-t) cos(t), e^(-t) sin(t), e^(-t) rangle` in three-space over the time interval `0 le t le 2 pi`.