/* [wxMaxima batch file version 1] [ DO NOT EDIT BY HAND! ]*/
/* [ Created with wxMaxima version 11.08.0 ] */

/* [wxMaxima: title   start ]
Moving in Three Dimensions
   [wxMaxima: title   end   ] */

/* [wxMaxima: comment start ]
Chapter 13 Technology Application Project
   [wxMaxima: comment end   ] */

/* [wxMaxima: section start ]
Introduction
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
OBJECTIVE: Learn to use Maxima to perform calculations to 
analyze motion in the equations that are given parameterically.

Moving in three dimensions is something with which we are 
familiar, but visualizing equations of motion and computations 
involved in analyzing that motion can be cumbersome. This 
module gives a broad overview on the analysis of motion in the 
equations that are given parametrically.
   [wxMaxima: comment end   ] */

/* [wxMaxima: subsect start ]
Technology Guidelines
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
NOTE:  If you have just finished a document, restart Maxima 
before executing a new document.  This can be done by choosing
"Restart" from the Maxima menu.

TO OPEN OR CLOSE CELLS
Click on the arrow at the top of the cell bracket. 

TO STOP AN EXECUTION
Click on STOP button from the toolbar. 

ORDER OF EXECUTION
Execute commands in the order given. Do not skip any Maxima
Input lines within a given document. 

Alternatively, you can execute the entire worksheet by 
selecting the "Evaluate All Cells" command from the "Cell" drop
down menu or simply press Ctrl-r. 

SAVING WORKSHEETS
You can save anytime to any directory you choose, and it is 
wise to save often.

EXPERIENCING MAJOR PROBLEMS
Save if appropriate, and then shut down Maxima and start it up 
again. 
   [wxMaxima: comment end   ] */

/* [wxMaxima: section start ]
Part I:  Parametric Equations of a Curve in
Three-Dimentions
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
First, we load the plots and plottools packages and the define 
x, t, y, and z coordinates for motion parametrically.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
reset()$
kill(all)$
load(draw)$
ratprint:false$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
x(t) := cos(t)$
y(t) := sin(t)$
z(t) := 4-t^2/25$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
Next, we plot the resulting curve in blue. Can you tell in 
which direction you are moving on the curve as t increases?
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
p1: parametric(x(t),y(t),z(t),t,0,10)$
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
If you consider the parametric equation as a vector equation 
for the motion of a particle, the derivative of that vector is 
the velocity vector that we form by differentiating each 
component. The commands below plot both the first and second 
derivative of the vector position function.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
dx: diff(x(t),t,1)$
dy: diff(y(t),t,1)$
dz: diff(z(t),t,1)$
dvector: [dx,dy,dz]$
d2x: diff(x(t),t,2)$
d2y: diff(y(t),t,2)$
d2z: diff(z(t),t,2)$
d2vector: [d2x,d2y,d2z]$
p2: parametric(dx,dy,dz,t,0,10)$
p3: parametric(d2x,d2y,d2z,t,0,10)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("velocity vector: ")$
print(dvector)$
wxdraw3d(
    nticks=100,
    color=red,
    p2,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("acceleration vector: ")$
print(d2vector)$
wxdraw3d(
    nticks=100,
    color=cyan,
    p3,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    color=red,
    p2,
    color=cyan,
    p3,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: subsect start ]
You Try It: Part I
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
Define your own x, y, and z coordinates for motion by changing 
the entries for x, y, z. You may or may not want to change the 
time interval (timeinterval) over which you plot your 
functions.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
load(draw)$
ratprint:false$

x(t) := 10*exp(cos(t))$
y(t) := 10*sin(t^2/15)$
z(t) := t^2/2$
timeinterval: 10$

p1: parametric(x(t),y(t),z(t),t,0,timeinterval)$
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
dx: diff(x(t),t,1)$
dy: diff(y(t),t,1)$
dz: diff(z(t),t,1)$
dvector: [dx,dy,dz]$
d2x: diff(x(t),t,2)$
d2y: diff(y(t),t,2)$
d2z: diff(z(t),t,2)$
d2vector: [d2x,d2y,d2z]$
p2: parametric(dx,dy,dz,t,0,timeinterval)$
p3: parametric(d2x,d2y,d2z,t,0,timeinterval)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("velocity vector: ")$
print(dvector)$
wxdraw3d(
    nticks=100,
    color=red,
    p2,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("acceleration vector: ")$
print(d2vector)$
wxdraw3d(
    nticks=100,
    color=cyan,
    p3,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    color=red,
    p2,
    color=cyan,
    p3,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
The plot of all three together may or may not be instructive.

Compare your results to those of your classmates.
   [wxMaxima: comment end   ] */

/* [wxMaxima: section start ]
Part II: Equations of Motion 
Velocity -> Position
Velocity -> Acceleration
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
If you are given the velocity, you can differentiate to find 
the acceleration and integrate to find the displacement. Since 
we do three integrations, we need three arbitrary constants to 
get the general solution.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
load(draw)$
ratprint:false$

velocity: [3*sqrt(t+1)/2,%e^(-t),1/(t+1)]$
acceleration: diff(velocity,t,1)$
position: integrate(velocity, t)$
rgeneral: position+[a,b,c]$

print("")$
print("velocity vector:")$
print(velocity)$
print("")$
print("acceleration vector:")$
print(acceleration)$
print("")$
print("general r (position) vector:")$
print(rgeneral)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
s: subst(0,t,rgeneral);
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
a: rhs((solve(s[1]=1,a))[1])$
b: rhs((solve(s[2]=1,b))[1])$
c: rhs((solve(s[3]=0,c))[1])$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
r: position+[a,b,c]$
print("")$
print("r (position) vector:")$
print(r)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: subsect start ]
You Try It: Part II
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
Change the entries in your velocity vector and your initial 
position (velocity and initial), and re-execute the section.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
load(draw)$
ratprint:false$

velocity: [3*sqrt(t+1)/2,%e^(-t),1/(t+1)]$
acceleration: diff(velocity,t,1)$
position: integrate(velocity, t)$
rgeneral: position+[a,b,c]$

print("")$
print("velocity vector:")$
print(velocity)$
print("")$
print("acceleration vector:")$
print(acceleration)$
print("")$
print("general r (position) vector:")$
print(rgeneral)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
s: subst(0,t,rgeneral)$

a: rhs((solve(s[1]=1,a))[1])$
b: rhs((solve(s[2]=1,b))[1])$
c: rhs((solve(s[3]=0,c))[1])$

r: position+[a,b,c]$
print("")$
print("r (position) vector:")$
print(r)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: section start ]
Part III: Computing the Distance Traveled
on a Curved Path
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
First we define the x, y, and z coordinates for motion 
parametrically.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
numer:false$
load(draw)$
ratprint:false$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
r: [3*cos(2*t),t^(10/3)/50,t]$
v: diff(r,t,1)$
speed: sqrt(v.v)$
distance: quad_qags(speed,t,0,10)[1]$
p1: parametric(r[1],r[2],r[3],t,0,10)$
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("The position given in feet is:")$
print(r)$
print("")$
print("The velocity given in feet per minute is:")$
print(v)$
print("")$
print("The speed, ds/dt, given in feet per minute is:")$
print(speed)$
print("")$
print("The distance traveled in feet in 10 minutes is:")$
print(distance)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
In the above example, if you had tried to integrate 
symbolically, instead of numerically, what would have happened?
(The quad_qags( ) function signifies numerical integration.) 
Many integrals arising in the computation of arc length are 
too complicated to evaluate symbolically.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
integrate(speed,t,0,10);
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
numer:true$
print("Speed at 5 minutes in feet per minute is:")$
print(subst(5,t,speed))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
wxdraw2d(
    nticks=100,
    color=red,
    explicit(speed,t,0,10),
    xlabel="t",
    ylabel="speed");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
If t represents time in minutes, what does this say about the 
speed of the hiker? Is it possible?
   [wxMaxima: comment end   ] */

/* [wxMaxima: subsect start ]
You Try It: Part III
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
Choose your own position function (r) and the interval over 
which you wish to integrate.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
numer:false$
load(draw)$
ratprint:false$

r: [2*t^(.3),cos(t-1),exp(-t*t)]$
v: diff(r,t,1)$
speed: sqrt(v.v)$
distance: quad_qags(speed,t,0,10)[1]$
p1: parametric(r[1],r[2],r[3],t,0,10)$
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    xlabel="x",
    ylabel="y",
    zlabel="z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("The position given in feet is:")$
print(r)$
print("")$
print("The velocity given in feet per minute is:")$
print(v)$
print("")$
print("The speed, ds/dt, given in feet per minute is:")$
print(speed)$
print("")$
print("The distance traveled in feet in 10 minutes is:")$
print(distance)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: section start ]
Part IV: Computing Curvature and Torsion
for a Space Curve
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
Maxima simplifies the process of finding curvature and torsion.
The computations below use formulas directly from your text, 
and we graphically explore the interpretations of the 
curvature and torsion functions. 

The vect package define the operator ~ as the cross product.
So v~a finds the cross product of vectors v and a.

Note that the package lapack may take some time to load, 
especially if this is the first time loading this package in 
Maxima.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
numer:false$
load(draw)$
load(vect)$
ratprint:false$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
mag(v) := (sqrt(v.v))$
r: [10*sin(t),10*exp(-t),t^(10/3)/100]$
v: diff(r,t,1)$
a: diff(v,t,1)$
speed: mag(v)$
utan: [v[1]/speed,v[2]/speed,v[3]/speed]$
crossva: express(v~a)$
transcrossva: transpose(crossva)$
curvature: mag(crossva)/speed^3$
m2: matrix(v,a,diff(a,t,1))$
det2: determinant(m2)$
prod2: crossva.transcrossva$
torsion: det2/prod2$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("The position vector is:")$
print(r)$
print("")$
print("The velocity vector is:")$
print(v)$
print("")$
print("The acceleration vector is:")$
print(a)$
print("")$
print("The speed is:")$
print(speed)$
print("")$
print("The unit tangent vector is:")$
print(utan)$
print("")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("The curvature is:")$
print(curvature)$
print("")$
print("The torsion is:")$
print(torsion)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
p1: parametric(r[1],r[2],r[3],t,0,10)$
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    xlabel="x",
    ylabel="y",
    zlabel="Z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
p2: explicit(speed,t,0,10)$
wxdraw2d(
    nticks=200,
    color=red,
    p2,
    xlabel="t",
    ylabel="speed");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
p3: explicit(curvature,t,0,10)$
wxdraw2d(
    nticks=200,
    color=green,
    p3,
    xlabel="t",
    ylabel="curvature");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
p4: explicit(torsion,t,0,10)$
wxdraw2d(
    nticks=200,
    color=blue,
    p4,
    xlabel="t",
    yrange=[-1,10],
    ylabel="torsion");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
wxdraw2d(
    nticks=200,
    color=red,
    p2,
    color=green,
    p3,
    color=blue,
    p4,
    yrange=[-10,20],
    xlabel="t");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
Look at your space curve and see if you can identify the points
on your space curve at which the curvature or torsion spike. 
What has caused these events? Are they related at all to the 
speed?

Recall from the definition of curvature that relates it to the 
cross product of the velocity and acceleration that the 
curvature is smallest when the velocity and acceleration are 
in the same direction. The curvature is largest when the 
velocity and acceleration are perpendicular to one another.

Looking at the torsion as the dot product of the derivative of 
the unit binormal and the normal vectors. Torsion will peak 
when those two vectors are in the same direction and be close 
to 0 when those two vectors are perpendicular.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    xlabel="x",
    ylabel="y",
    zlabel="Z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
You can see the TNB frame for a curve in action in Part VI of 
this module.
   [wxMaxima: comment end   ] */

/* [wxMaxima: subsect start ]
You Try It: Part IV
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
Redefine your position function (r) and compute. It should be 
noted that the computations here are not trivial, so more 
complicated functions could lead to long waits.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
numer:false$
load(draw)$
load(vect)$
ratprint:false$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
mag(v) := (sqrt(v.v))$
r: [2*t^.3,cos(t-1),%e^(-t*t)]$
v: diff(r,t,1)$
a: diff(v,t,1)$
speed: trigsimp(mag(v))$
utan: [v[1]/speed,v[2]/speed,v[3]/speed]$
crossva: express(v~a)$
transcrossva: transpose(crossva)$
curvature: trigsimp(mag(crossva)/speed^3)$
m2: matrix(v,a,diff(a,t,1))$
det2: determinant(m2)$
prod2: crossva.transcrossva$
torsion: trigsimp(det2/prod2)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("The position vector is:")$
print(r)$
print("")$
print("The velocity vector is:")$
print(v)$
print("")$
print("The acceleration vector is:")$
print(a)$
print("")$
print("The speed is:")$
print(speed)$
print("")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("The unit tangent vector is:")$
print(utan)$
print("")$
print("The curvature is:")$
print(curvature)$
print("")$
print("The torsion is:")$
print(torsion)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
p1: parametric(r[1],r[2],r[3],t,0,10)$
wxdraw3d(
    nticks=100,
    color=blue,
    p1,
    xlabel="x",
    ylabel="y",
    zlabel="Z")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
p2: explicit(speed,t,0,10)$
wxdraw2d(
    nticks=200,
    color=red,
    p2,
    xlabel="t",
    ylabel="speed");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
p3: explicit(curvature,t,0,10)$
wxdraw2d(
    nticks=200,
    color=green,
    p3,
    xlabel="t",
    ylabel="curvature");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
p4: explicit(torsion,t,0,10)$
wxdraw2d(
    nticks=200,
    color=blue,
    p4,
    xlabel="t",
    yrange=[-1,10],
    ylabel="torsion");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
wxdraw2d(
    nticks=200,
    color=red,
    p2,
    color=green,
    p3,
    color=blue,
    p4,
    yrange=[-3,3],
    xrange=[0,10],
    xlabel="t");
/* [wxMaxima: input   end   ] */

/* [wxMaxima: section start ]
Part V:  TNB Frame - Computation and 
Visualization
   [wxMaxima: section end   ] */

/* [wxMaxima: comment start ]
The steps here resemble the ones above, except that we choose 
a simpler function and extend the formulas to the unit normal 
and binormal vectors.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
numer:false$
load(draw)$
load(vect)$
ratprint:false$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
mag(v) := (sqrt(v.v))$
r: [sin(2*t),cos(2*t),t^2/10]$
v: diff(r,t,1)$
a: diff(v,t,1)$
speed: trigsimp(mag(v))$
utan: [v[1]/speed,v[2]/speed,v[3]/speed]$
crossva: express(v~a)$
transcrossva: transpose(crossva)$
curvature: trigsimp(mag(crossva)/speed^3)$
m2: matrix(v,a,diff(a,t,1))$
det2: determinant(m2)$
prod2: crossva.transcrossva$
torsion: trigsimp(det2/prod2)$
top: trigsimp(diff(utan,t))$
bottom: trigsimp(radcan(mag(top)))$
un: trigsimp(radcan(top/bottom))$
ubn: trigsimp(radcan(express(utan~un)))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("The position vector is:")$
print(r)$
print("")$
print("The velocity vector is:")$
print(v)$
print("")$
print("The acceleration vector is:")$
print(a)$
print("")$
print("The speed is:")$
print(speed)$
print("")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("The unit tangent vector is:")$
print(utan)$
print("")$
print("The curvature is:")$
print(curvature)$
print("")$
print("The torsion is:")$
print(torsion)$
print("")$
print("The unit normal is:")$
print(un)$
print("")$
print("The unit binormal is:")$
print(ubn)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
Let's find the tangential and normal components of the 
accelerations.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
at: trigsimp(a*utan).[1,1,1]$
an: ratsimp(trigreduce(ratsimp((a*un)).[1,1,1]))$
print("")$
print("The tangential component of acceleration is ", at)$
print("The normal component of acceleration is ", an)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
You can get a visual perspective of the tangential and normal 
components of the acceleration in two dimensions by exploring 
the Java applet, "Tangent and Normal Vectors."

If you compare these components found by the dot products to 
the formulas given in your book, you will see that they agree, 
once the absolute value sign is removed from the curvature 
formula.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
print("The tangential component of acceleration is ", diff(speed,t,1))$
print("The normal component of acceleration is ", radcan(curvature*speed^2))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
The equality of the tangential components can be verified 
easily if you notice that the first two terms of the result 
found previously by dotting the acceleration vector into the 
unit tangent vector subtract out with one another. The equality
of the normal components is not as obvious. You can see the 
advantage in using the formulas used in the latter part, due to
the complexity of the unit normal vector. 

The following set of commands will plot the unit tangent, unit 
normal, and unit binormal vectors as you move along the curve. 
The resulting animated image will be saved in the wxMaxima 
folder with the filename "13_m3_1.gif".
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
p1: parametric(r[1],r[2],r[3],t,0,2*%pi)$
numer: true$
P: []$

for i: 0 thru 6.4 step .4 do block(
    p2: vector(subst(i,t,r),subst(i,t,utan)),
    p3: vector(subst(i,t,r),subst(i,t,un)),
    p4: vector(subst(i,t,r),subst(i,t,ubn)),
    p5: points([
        subst(i,t,r),
        subst(i,t,r+0.2*ubn),
        subst(i,t,r+0.2*utan+0.2*ubn),
        subst(i,t,r+0.2*utan),
        subst(i,t,r),
        subst(i,t,r+0.2*un),
        subst(i,t,r+0.2*utan+0.2*un),
        subst(i,t,r+0.2*utan),
        subst(i,t,r+0.2*utan+0.2*ubn),
        subst(i,t,r+0.2*utan+0.2*ubn+0.2*un),
        subst(i,t,r+0.2*utan+0.2*un),
        subst(i,t,r+0.2*utan+0.2*ubn+0.2*un),
        subst(i,t,r+0.2*ubn+0.2*un),
        subst(i,t,r+0.2*un),
        subst(i,t,r+0.2*ubn+0.2*un),
        subst(i,t,r+0.2*ubn)]),
                    
    P: append(P,[gr3d(
        nticks=100,
        color=black,
        line_width=2,
        p1,
        color=red,
        line_width=3,
        p2,
        color=green,
        p3,
        color=blue,
        p4,
        points_joined=true,
        color=black,
        line_width=1,
        point_type=-1,
        p5,
        xrange=[-1.5,1.5],
        yrange=[-1.5,1.5],
        zrange=[0,4],
        title="Tangent (Red), Unit Normal (Green), Unit Binormal (Blue)")]))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
draw(
    terminal=animated_gif,
    delay=10,
    file_name="13_m3_1",
    P)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: subsect start ]
You Try It: Part V
   [wxMaxima: subsect end   ] */

/* [wxMaxima: comment start ]
The computational code is written so that all you have to do is 
change the initial vector ( r ) function and then re-execute the 
commands. We suggest that you choose some of the homework 
problems from your text for the TNB frame, since the computation 
for the unit normal and unit binormal can get you very bogged 
down, even in Maple. You may also wish to change the domain over 
which you extend your function in the second section.

If you check out the following example, expect to wait a few 
minutes for the computations.
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
kill(all)$
numer:false$
load(draw)$
load(vect)$
ratprint:false$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
mag(v) := (sqrt(v.v))$
r: [2*t,cos(t-1),%e^(-t)]$
v: diff(r,t,1)$
a: diff(v,t,1)$
speed: trigsimp(mag(v))$
utan: [v[1]/speed,v[2]/speed,v[3]/speed]$
crossva: express(v~a)$
transcrossva: transpose(crossva)$
curvature: trigsimp(mag(crossva)/speed^3)$
m2: matrix(v,a,diff(a,t,1))$
det2: determinant(m2)$
prod2: crossva.transcrossva$
torsion: trigsimp(det2/prod2)$
top: trigsimp(diff(utan,t))$
bottom: trigsimp(radcan(mag(top)))$
un: trigsimp(radcan(top/bottom))$
ubn: trigsimp(radcan(express(utan~un)))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("")$
print("The position vector is:")$
print(r)$
print("")$
print("The velocity vector is:")$
print(v)$
print("")$
print("The acceleration vector is:")$
print(a)$
print("")$
print("The speed is:")$
print(speed)$
print("")$
print("The unit tangent vector is:")$
print(utan)$
print("")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("The curvature is:")$
print(curvature)$
print("")$
print("The torsion is:")$
print(torsion)$
print("")$
print("The unit normal is:")$
print(un)$
print("")$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
print("The unit binormal is:")$
print(ubn)$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: comment start ]
The following set of commands will plot the unit tangent, unit 
normal, and unit binormal vectors as you move along the curve. 
The resulting animated image will be saved in the wxMaxima 
folder with the filename "13_m3_2.gif".
   [wxMaxima: comment end   ] */

/* [wxMaxima: input   start ] */
p1: parametric(r[1],r[2],r[3],t,0,10)$
numer: true$
P: []$

for i: 0 thru 10 step .5 do block(
    p2: vector(subst(i,t,r),subst(i,t,utan)),
    p3: vector(subst(i,t,r),subst(i,t,un)),
    p4: vector(subst(i,t,r),subst(i,t,ubn)),
    p5: points([
        subst(i,t,r),
        subst(i,t,r+0.2*ubn),
        subst(i,t,r+0.2*utan+0.2*ubn),
        subst(i,t,r+0.2*utan),
        subst(i,t,r),
        subst(i,t,r+0.2*un),
        subst(i,t,r+0.2*utan+0.2*un),
        subst(i,t,r+0.2*utan),
        subst(i,t,r+0.2*utan+0.2*ubn),
        subst(i,t,r+0.2*utan+0.2*ubn+0.2*un),
        subst(i,t,r+0.2*utan+0.2*un),
        subst(i,t,r+0.2*utan+0.2*ubn+0.2*un),
        subst(i,t,r+0.2*ubn+0.2*un),
        subst(i,t,r+0.2*un),
        subst(i,t,r+0.2*ubn+0.2*un),
        subst(i,t,r+0.2*ubn)]),
                    
    P: append(P,[gr3d(
        nticks=100,
        color=black,
        line_width=2,
        p1,
        color=red,
        line_width=3,
        p2,
        color=green,
        p3,
        color=blue,
        p4,
        points_joined=true,
        color=black,
        line_width=1,
        point_type=-1,
        p5,
        xrange=[-1,21],
        yrange=[-2,2],
        zrange=[-1,2],
        title="Tangent (Red), Unit Normal (Green), Unit Binormal (Blue)")]))$
/* [wxMaxima: input   end   ] */

/* [wxMaxima: input   start ] */
draw(
    terminal=animated_gif,
    delay=10,
    file_name="13_m3_2",
    P)$
/* [wxMaxima: input   end   ] */

/* Maxima can't load/batch files which end with a comment! */
"Created with wxMaxima"$