Information for Maple assignment 2
Each student will get individual data for the assignment.
Below is helpful background information for the second Maple
assignment. I repeat what I wrote in the last set of instructions:
You are allowed, and indeed encouraged, to copy and modify the commands discussed here.
Packages
For this assignment you will need commands from the plots and
VectorCalculus packages. Load these packages using the
with command, as discussed in Packages section of the Background Information for the first Maple
assignment.
Plotting Curves
We will again be using the spacecurve command from the plots package. For these instructions, we will be using as an example the curve given by:
r(t) = < 5 cos(2t) + 5 sin(t), 2 cos(2t) + 4 cos(t), sin(3t) + 4 sin(t)>
or if you prefer the parametric notation:
x=5 cos(2t) + 5 sin(t)
y=2 cos(2t) + 4 cos(t)
z= sin(3t) + 4 sin(t)
I'm beginning with this curve because it is fairly simple. More complicated curves will appear later. This curve is a closed
curve (i.e. a loop) on the interval [0,2π] because the formulas are all
combinations of sines and cosines with integers multiplying the
parameter, t.
Since we will be doing a lot with this curve, it's wise to store it to a variable. There are several approaches to doing this, but the best for our purposes is storing it as a "live" function, r(t)
r := t -> < 5*cos(2*t) + 5*sin(t), 2*cos(2*t) + 4*cos(t), sin(3*t) + 4*sin(t)>
This should allow us to call r(t) simply by typing r(t), but also evaluate r at a value, say t=2, by typing in r(2). (Another option would have been to store the right hand expression directly to r, omitting the "t ->" portion. In this case we would need to use the subs command to evaluate r at a specific value.)
Let's start by just graphing our curve.
spacecurve(r(t),t=0..2*Pi);
We can use just r(t) because we've already stored the equation to r(t). Notice that we used the capitalization Pi to get the constant 3.141... (pi and PI do not give this constant). The result is the picture shown below.
This picture is not very useful. The rainbow coloration does not add anything and will make the curve fade away if we print out the picture in black and white. It is not clear whether that crossing in the middle is actually a crossing or just an illusion due to the viewpoint. Plus there are no axes to give reference so do not know how far to the left/right, up/down, front/back the curve stretches. Let's add some options to this command to make the picture better. You enter options in the command after the range, separated by commas.
- The thickness option will make the curve a little thicker and easier to see. The default is thickess=1, but let's make thickess=2 to beef up our curve slightly.
- We can make the whole curve a single color, say black. We do this via the option color=black. Current versions of Maple have hundreds of pre-set colors ranging from black to beige to SteelBlue. Feel free to experiment, but remember the limitations of your printer for the end result. (For the record Maple is Canadian, so it also accepts the "colour" option.)
- We will add the standard coordinate axes with the option axes=normal. Other choices for axes include axes=none (the default), axes=boxed, and axes=frame.
- We can label the x-, y-, and z-axis with the option labels=[x,y,z].
- We can make the curve smoother using the numpoints option. Maple plots curves by evaluating r(t) at many different values for t, and then connecting the dots. A higher numpoints setting makes more dots, thus smoothing the curve but lengthening the computation time. Some experimentation is often necessary to strike a balance. For this assignment numpoints=150 works well. For contrast, try numpoints=15 (chunky curve, short computation) and numpoints=1500 (smoother curve, longer computation).
- You can use the mouse to adjust the viewpoint, but you can also use the orientation command. Often the orientation command is not necessary because of the mouse, but it's good to know it's there if you need it.
- You can look up the list of all available options for three-dimensional plotting by looking up help(plot3d/options).
Putting these options into the command we get
spacecurve(r(t),t=0..2*Pi, thickness=2, color=black, axes=normal, labels=[x,y,z], numpoints=150);
Certainly the view shown is much nicer, and has some readable
quantitative information (on the axes). It is still a bit
deceptive, however, at the center. The picture seems to show the curve intersecting itself. In
fact, this particular curve does not have any
self-intersections. This is what's called a simple closed
curve. If you right-click on a Maple
picture, the Manipulator button allows you to rotate or
scale or pan. So you can get many different views. Here
is a view which shows the same curve from a different angle, and zoomed out farther. Apparently this curve does not intersect
itself!
Consider the following picture of the curve, s(t)=< sin(2t), cos(2t), sin(t) >
This shares the a figure-8 shape with our original picture of our curve r(t). This curve actually has an intersection, however, as illustrated by the animation below (created by the animate command - hover your cursor over the animate for the Maple code.)
Since you cannot print an animation, you should include several views of your curve to convince the reader of any intersections.
Warning: There are curves which do not intersect themselves,
but have the irritating (interesting?) property that every two dimensional
view shows intersections! See the bottom of this page for examples. For example, q(t)=<sin(t),cos(3t),sin(2t)> (animated below) has this property. In this case you would include several views of your curve to convince the reader that there are no intersections.
Watch as intersections appear and disappear. Compare this to a curve with an actual intersection, such as s(t)=<sin(2t), cos(2t), sin(t)> above where the intersection does not disappear. More examples of such curious curves appear at the bottom of this page.
Analyzing Curves
We can try various views to get the size of the curve. Two
images are shown below (here we used the "axes=frame" option instead of "axes=normal"). The image on the left has the z-axis sticking directly out from the
viewing plane, so we are looking down from above. Alternately, you can see this as projecting (i.e. squashing) the curve into the xy-plane. The image on the right has the y-axis sticking directly out from the viewing plane. Therefore it shows a picture of the curve squashed
into the xz-plane.
From these pictures I can get a rough idea of the dimensions of this
curve. Apparently the curve fits inside a box with -10≤x≤6 and
-3≤y≤6 and -3.6≤z≤3.6, as illustrated below. You may draw such markings by hand.
How long is this curve? This is called the arc length, which is calculated by integrating the speed of r(t) for the desired range of t (i.e. ∫ab |r'(t)| dt ). We need to differentiate r(t), and then find its length or magnitude (this is the speed). The diff command will differentiate r(t). The VectorCalculus package contains a command called Norm which computes the norm of a vector. We will use the variable speed to represent |r'(t)|.
speed:=Norm( diff(r(t),t) );
speed:=((-10*sin(2*t)+5*cos(t))^2+(-4*sin(2*t)-4*sin(t))^2+(3*cos(3*t)+4*cos(t))^2)^(1/2)
We are interested in the length of one loop of the curve. For our curve "one loop" is given by the range 0 ≤ t ≤ 2π, so we need to integrate speed from t=0 to
t=2*Pi to get the length. We can integrate using the int command, so may naturally try:
int(speed, t=0..2*Pi);
Notice that when you use int here, however, Maple usually returns the original integral, not evaluated, after about 12 seconds' thought. When Maple returns the original command, it means that Maple has given up trying to evaluate whatever you asked. Usually this is because you typed in something that the program does not recognize, but in this case it is because Maple cannot solve this integral symbolically. This is not surprising -- most functions defined by formulas can't be antidifferentiated in terms of familiar functions.
Maple needs to be "convinced" to evaluate this integral
approximately. You can do this via:
evalf(Int(speed,t=0..2*Pi));
56.13189562
Please note the capital I in the
integrate command. This is an important difference which you can read
about on a help page if you want. The numerical computation of the
integral took less than a hundredth (.01) of a second. This curve
is about 56 units long.
Warning: Do not confuse the Norm command in Maple with Norm from This Old House and The New Yankee Workshop, nor with Norm from Cheers. One is a command, one is a carpenter, and one is a carouser.
The VectorCalculus package also contains the Curvature command,
which computes the curvature of a curve.
The curvature is often a very messy formula:
Curvature(r(t),t);
((-1/2/(-200*sin(t)*cos(t)^2+144*cos(t)^6-584*cos(t)^4+498*cos(t)^2+64*cos(t)+
16-64*cos(t)^3)^(3/2)*(-10*sin(2*t)+5*cos(t))*(-200*cos(t)^3+400*sin(t)^2*cos(
t)-864*cos(t)^5*sin(t)+2336*cos(t)^3*sin(t)-996*sin(t)*cos(t)-64*sin(t)+192*
sin(t)*cos(t)^2)+1/(-200*sin(t)*cos(t)^2+144*cos(t)^6-584*cos(t)^4+498*cos(t)^
2+64*cos(t)+16-64*cos(t)^3)^(1/2)*(-20*cos(2*t)-5*sin(t)))^2+(-1/2/(-200*sin(t
)*cos(t)^2+144*cos(t)^6-584*cos(t)^4+498*cos(t)^2+64*cos(t)+16-64*cos(t)^3)^(3
/2)*(-4*sin(2*t)-4*sin(t))*(-200*cos(t)^3+400*sin(t)^2*cos(t)-864*cos(t)^5*sin
(t)+2336*cos(t)^3*sin(t)-996*sin(t)*cos(t)-64*sin(t)+192*sin(t)*cos(t)^2)+1/(-\
200*sin(t)*cos(t)^2+144*cos(t)^6-584*cos(t)^4+498*cos(t)^2+64*cos(t)+16-64*cos
(t)^3)^(1/2)*(-8*cos(2*t)-4*cos(t)))^2+(-1/2/(-200*sin(t)*cos(t)^2+144*cos(t)^
6-584*cos(t)^4+498*cos(t)^2+64*cos(t)+16-64*cos(t)^3)^(3/2)*(3*cos(3*t)+4*cos(
t))*(-200*cos(t)^3+400*sin(t)^2*cos(t)-864*cos(t)^5*sin(t)+2336*cos(t)^3*sin(t
)-996*sin(t)*cos(t)-64*sin(t)+192*sin(t)*cos(t)^2)+1/(-200*sin(t)*cos(t)^2+144
*cos(t)^6-584*cos(t)^4+498*cos(t)^2+64*cos(t)+16-64*cos(t)^3)^(1/2)*(-9*sin(3*
t)-4*sin(t)))^2)^(1/2)/(-200*sin(t)*cos(t)^2+144*cos(t)^6-584*cos(t)^4+498*cos
(t)^2+64*cos(t)+16-64*cos(t)^3)^(1/2)
Ugh. This is certainly an example of why people use : instead of
;. A command ending with a colon will not display its output. Knowledge
of the specific formula of the curvature is not likely to be useful, so we can safely omit it
(use Ctrl+Del to delete a line from a Maple document).
We can compute the curvature at, for example, t=2 using the subs command.
evalf( subs(t=2, Curvature(r(t),t) ));
.2586332282
The subs command substitutes 2 for t in the result
of the Curvature(r(t),t) command, and then evalf
asks for a floating point (i.e. decimal) approximation. Pay attention to the "nesting" of the commands. Just as in algebra, Maple works outward from the inside, first calculating Curvature(r(t),t), then substituting 2 for t, then evaluating the expression to get a decimal approximation.
The curvature expression for r(t) is a (complicated) real-valued function of t, so we can plot it like any other one-variable real-valued function.
plot(Curvature(r(t),t),t=0..2*Pi,thickness=2);
The height of this function at t=t0 gives the curvature of r(t) at the point r(t0). We can use this graph to locate (approximately) the curviest (most
curvy?) point on the curve. If you position your mouse cursor on the
graph, and left click, the cursor position relative to the graph is
shown on the left side of the toolbar at the top of the window. In this case, the coordinates
that I see for the marked maximum are about (4.69, 1.81). So the highest
curvature in this curve is about 1.81 and this occurs when the
parameter t is about 4.69. Maple does have commands that find the minima and maxima, but they are often very slow, especially for functions as complicated as our Curvature(r(t),t). For our purposes these approximations should suffice.
Where (in space) is the curviest point? That is, what are the coordinates for the point on r(t) corresponding to the value of t that gives the highest curvature? We must substitute the value of the
parameter, t=4.69, into the curve function. We can either use subs(t=4.69, r(t)) or more simply r(4.69). Both should return a vector approximately < -9.994,-2.086,-3.001 >, which means the most curvy point has coordinates approximately (-9.994,-2.086,-3.001).
To graph a point, in our case the most curvy point, we can use the pointplot3d command (in the plots package). For example, try
pointplot3d([-9.994,-2.086,-3.001], symbol=box, symbolsize=20, color=blue);
You will get a thoroughly unhelpful picture (although it's more interesting than if you do not include any of the options), since this point is only meaningful relative to the curve r(t). Use the display command
(discussed in the
Background for the first assignment)
to plot both the curve and the point together.
POINT:=pointplot3d([-9.994,-2.086,-3.001], symbol=box, symbolsize=20, color=blue):
CURVE:=spacecurve(r(t),t=0..2*Pi,thickness=2, color=black, axes=normal, labels=[x,y,z], numpoints=150):
display({MCP,CURVE});
This is a bit confusing, since the little blue box
is supposed to be the most curvy point, but the curve does not look very curvy at that point compared to other parts of the curve. But then I adjusted the
image by adding the scaling=constrained option and the
result is the more convincing picture below. So, again:
pictures are wonderful but they can be deceptive unless you stay
alert!
Round-off Warning: Round-off and estimation errors from earlier steps may cause the point you plot to be slightly "off" the actual point of maximum curvature. If this happens, you may want to go back to previous steps and add an extra one or two decimal digits of precision.
I chose the example above to analyze in detail because it is
fairly simple. You will get randomly generated curves which may be a bit more wild. You have been
given the curvature of your curve at t=2, so that you can
check initially that you've entered the functions correctly.
Weird Pictures
Some strange things can happen. Let me show you a few of them.
Below is a plot of curvature for a curve which isn't very nice.
The curvature has an enormous thin peak. Experiments seem to show that
such behavior is common.
Below is the curve which has the curvature graph shown above.
On the left is an unconstrained view, and on the right, a
constrained view (observe the scale of the y-axis in the first picture).
The very high curvature on the right side of the curve corresponds to the
single high peak for the curvature graph. The other bends of the curve also
make the curvature increase and decrease, but these are drowned out in the curvature graph above
due to this single extreme bend.
Below are three pictures of another curve. I emphasize that these are
all pictures of the same curve. This curve does not
intersect itself, but every two dimensional view does seem to show an
intersection. So detecting intersections by "just looking" may be a
bit difficult. Some examination of the images is necessary.
Here are two pictures of yet another closed curve. This closed curve
does have two self-intersections, each of which can appear to be an illusion. It can be difficult to convince people of these intersections without much numerical work (e.g. showing that r(t0)=r(t1) for some 0 ≤ t0 < t1 < 2 π) or many carefully chosen
pictures. If you have suggestions, please let me know.
I am very curious to learn if this homework assignment will lead to
any bizarre pictures and numbers. Please tell me if you think you have
found some.
I thank Ms. L. Pudwell for several useful conversations about this
material.
Note on torsion: After experiments with some students, I am very
glad that torsion (roughly the twistiness of the curve) was not mentioned in this assignment! The torsion
computations are quite elaborate, taking large amounts of
memory and computational time. I think we could have frozen huge
chunks of university computer resources with any general requests for
torsion. I know that I caused my computers at home and at school to be
very unhappy by asking for the torsion of some students' curves.
Maintained by
greenfie@math.rutgers.edu and last modified 8/3/2009 by Andrew Baxter.