Previous
Contents
Next

Geometry: Tubes around Spacial Curves


Outline of the example:

For an arbitrary curve in three-dimensional space, we use Frenet's formulas to compute a tube around that curve. A tube around a curve is a surface in space; a surface in space is conveniently represented by a vector of three functions

 [ fx(u, v), fy(u, v), fz(u, v)]

where u, v are parameters.

To compute a tube around a curve, we need the normal and the binormal for every curve point. These vectors are obtained by the Frenet formulas:

The tangent:

      g'
t =  ----
     |g'|

The binormal:

      g' × g"
b =  ---------
     |g' × g"|

The normal:

 n = b × t

In these formulas, the cross denotes the cross product between two spacial vectors. Note that g is a function of a parameter; by consequence the vectors t, b, n are also functions of that same parameter.

With these vectors, a tube around the curve is easily defined:

tube(g, r, u, v) :=
    g(u) + r*(cos(v)*n(u) + sin(v)*b(u))

where r is the radius of the tube, u is the parameter of the spacial curve and v is a parameter that runs from 0 to 2*pi.

Note that for the tube we need only the normal and the binormal. Nevertheless we need all the Frenet formulas: The tangent is needed to compute the normal.

We will need the cross product of two vectors a few times, so let us begin with a definition for it:

cross (v1, v2) := determinant(matrix([[1, 0, 0], [0, 1, 0], [0, 0,1]], v1, v2));

Here is an example:

 cross([12, -1, 3], [-9, -2, -2]);
  [8, - 3, - 33]

A quite tricky definition is this one:

First we define a function to rotate the elements of list one position to the left. Here we use our knowledge of Lisp to write:

  (DEFUN rotateLeft
         (Lambda (V) (APPEND (CDR V) (LIST (CAR V))))
  )

Translated into Maxima, this is:

rotateLeft(v1):= append(rest(v1, 1), cons(first(v1), []));

Now comes the trick:

cross (u, v) := rotateLeft( u * rotateLeft(v) - v * rotateLeft(u));

We try that with two general vectors:

 cross([ax, ay, az], [bx, by, bz]);
 [ay bz - az by, az bx - ax bz, ax by - ay bx]

This is in fact the formula that we find in every textbook about analytical geometry.

Now we use Frenet's formulas to compute the tangent, the normal and the binormal for an arbitrary point of a spacial curve:

The tangent:

tangent(fn, x) :=
       diff(fn(x), x) / ((diff(fn(x), x) . diff(fn (x), x))^(1/2))$

The binormal:

binormal (fn, x) := 
         cross (diff(fn(x), x), 
                diff (fn(x), x, 2)
               )  
           /(( cross(diff(fn(x), x), diff (fn(x), x, 2))
               . cross(diff(fn(x), x), diff (fn(x), x, 2))
             ) ^(1/2)
            )$

And the normal:

normal (fn, x) := cross(binormal(fn, x), tangent(fn, x))$

To construct a tube, we use the normal and the binormal vector to drw a circle around a curve point:

 tube(fn, r, t, phi) :=
       fn(t) + r* (cos(phi)* normal(fn, t) + sin(phi)* binormal(fn, t))$

You find these definitions in tubeplot1.mc.

Now let us try these definitions with a very beautiful spacial curve: a helix.

helix(x) := [5*cos(x), 5*sin(x), 3*x];

Maxima answers the definition:

 helix(x) := [5 cos(x), 5 sin(x), 3 x]

The definition is easy to understand:

[5*cos(x), 5*sin(x)]

is a circle with radius 5.
The third coordinate
3*x

translates the points of circle into space. The translation distance for a full revolution of a point on the circle is 2*pi*3.

Now we can try our definitions.

tangent(helix, t);
		      5 sin(t)
(%o14) [- ---------------------------------, 
		     2		  2
	  sqrt(25 sin (t) + 25 cos (t) + 9)

		      5 cos(t)				     3
	  ---------------------------------, ---------------------------------]
		     2		  2		        2	     2
	  sqrt(25 sin (t) + 25 cos (t) + 9)  sqrt(25 sin (t) + 25 cos (t) + 9)
binormal (helix, t);
				 15 sin(t)
(%o15) [------------------------------------------------------------, 
		    2		 2    2		 2	       2
	sqrt((25 sin (t) + 25 cos (t))  + 225 sin (t) + 225 cos (t))

			   15 cos(t)
- ------------------------------------------------------------, 
	      2		   2    2	   2		 2
  sqrt((25 sin (t) + 25 cos (t))  + 225 sin (t) + 225 cos (t))

			2	     2
		  25 sin (t) + 25 cos (t)
------------------------------------------------------------]
	    2		 2    2		 2	       2
sqrt((25 sin (t) + 25 cos (t))  + 225 sin (t) + 225 cos (t))
(%i13) normal(helix, t);
			  2	       2
(%o16) [- 5 cos(t) (25 sin (t) + 25 cos (t))

	     2		  2		        2	     2	  2
/(sqrt(25 sin (t) + 25 cos (t) + 9) sqrt((25 sin (t) + 25 cos (t))

	  2	        2
 + 225 sin (t) + 225 cos (t))) - 45 cos(t)

	     2		  2		        2	     2	  2
/(sqrt(25 sin (t) + 25 cos (t) + 9) sqrt((25 sin (t) + 25 cos (t))

	  2	        2			  2	       2
 + 225 sin (t) + 225 cos (t))), - 5 sin(t) (25 sin (t) + 25 cos (t))

	     2		  2		        2	     2	  2
/(sqrt(25 sin (t) + 25 cos (t) + 9) sqrt((25 sin (t) + 25 cos (t))

	  2	        2
 + 225 sin (t) + 225 cos (t))) - 45 sin(t)

	     2		  2		        2	     2	  2
/(sqrt(25 sin (t) + 25 cos (t) + 9) sqrt((25 sin (t) + 25 cos (t))

	  2	        2
 + 225 sin (t) + 225 cos (t))), 0]

When we try

tube(helix, 1, t, rho);

We get a page full of unsimplified expressions.

The obvious problem is that these expressions were not simplified. This is not the right way to do symbolic computations. We need to care for simplification and, as a general rule, we should simplify an expression as soon as possible.

For this example, simplification can be easily achieved with the trigsimp function:

trigsimp(tangent(helix, t));
			    5 sin(t)  5 cos(t)	   3
      			 [- --------, --------, --------]
			    sqrt(34)  sqrt(34)  sqrt(34)
trigsimp(binormal(helix, t));
			3 sin(t)    3 cos(t)	 5
(%o18) 		       [--------, - --------, --------]
			sqrt(34)    sqrt(34)  sqrt(34)
trigsimp(normal(helix, t));
  [- cos(t), - sin(t), 0]

In this example, trigonometric simplification was both helpful and necessary. Let us examine a different spacial curve before we draw our conclusions about improved definitions for the functions tangent, binormal, normal.

w(x) := [2*t^2 - t, t^2, 2*t^2 + 4*t];
				  2	  2     2
(%o8) 		      w(x) := [2 t  - t, t , 2 t  + 4 t]
tg: tangent(w, t);
		     4 t - 1
(%o9) [------------------------------------, 
		     2		  2	 2
       sqrt((4 t + 4)  + (4 t - 1)  + 4 t )

		    2 t				        4 t + 4
    ------------------------------------, ------------------------------------]
		  2	       2      2		        2	     2	    2
    sqrt((4 t + 4)  + (4 t - 1)  + 4 t )  sqrt((4 t + 4)  + (4 t - 1)  + 4 t )

It is obvious that trigonometric simplification will not help here, but it is equally obvious that all three denominators can be simplified. We try radcan, which is often useful to rewrite radicals in canonical form:

radcan(tg);
		4 t - 1			   2 t
(%o10) [-----------------------, -----------------------, 
		 2			  2
	sqrt(36 t  + 24 t + 17)	 sqrt(36 t  + 24 t + 17)

							       4 t + 4
						       -----------------------]
								2
						       sqrt(36 t  + 24 t + 17)

This looks a lot better.

For the binormal and the normal we obtain:

bn: binormal(w, t);
							    2
(%o11) [(8 t - 2 (4 t + 4))/sqrt((4 (4 t + 4) - 4 (4 t - 1))

		      2			     2
 + (8 t - 2 (4 t + 4))  + (2 (4 t - 1) - 8 t) ), 

							    2
(4 (4 t + 4) - 4 (4 t - 1))/sqrt((4 (4 t + 4) - 4 (4 t - 1))

		      2			     2
 + (8 t - 2 (4 t + 4))  + (2 (4 t - 1) - 8 t) ), 

						    2			   2
(2 (4 t - 1) - 8 t)/sqrt((4 (4 t + 4) - 4 (4 t - 1))  + (8 t - 2 (4 t + 4))

		      2
 + (2 (4 t - 1) - 8 t) )]
radcan(bn);
			  4	      10	    1
(%o12) 		   [- ----------, ----------, - ----------]
		      3 sqrt(13)  3 sqrt(13)    3 sqrt(13)
nm: normal(w, t);
(%o13) [(4 t + 4) (4 (4 t + 4) - 4 (4 t - 1))

		2	     2	    2				       2
/(sqrt((4 t + 4)  + (4 t - 1)  + 4 t ) sqrt((4 (4 t + 4) - 4 (4 t - 1))

		      2			     2
 + (8 t - 2 (4 t + 4))  + (2 (4 t - 1) - 8 t) ))

					  2	       2      2
 - 2 t (2 (4 t - 1) - 8 t)/(sqrt((4 t + 4)  + (4 t - 1)  + 4 t )

				 2		        2
 sqrt((4 (4 t + 4) - 4 (4 t - 1))  + (8 t - 2 (4 t + 4))

		      2
 + (2 (4 t - 1) - 8 t) )), (4 t - 1) (2 (4 t - 1) - 8 t)

		2	     2	    2				       2
/(sqrt((4 t + 4)  + (4 t - 1)  + 4 t ) sqrt((4 (4 t + 4) - 4 (4 t - 1))

		      2			     2
 + (8 t - 2 (4 t + 4))  + (2 (4 t - 1) - 8 t) ))

						2	     2	    2
 - (4 t + 4) (8 t - 2 (4 t + 4))/(sqrt((4 t + 4)  + (4 t - 1)  + 4 t )

				 2		        2
 sqrt((4 (4 t + 4) - 4 (4 t - 1))  + (8 t - 2 (4 t + 4))

		      2
 + (2 (4 t - 1) - 8 t) )), 2 t (8 t - 2 (4 t + 4))

		2	     2	    2				       2
/(sqrt((4 t + 4)  + (4 t - 1)  + 4 t ) sqrt((4 (4 t + 4) - 4 (4 t - 1))

		      2			     2
 + (8 t - 2 (4 t + 4))  + (2 (4 t - 1) - 8 t) ))

 - (4 t - 1) (4 (4 t + 4) - 4 (4 t - 1))

		2	     2	    2				       2
/(sqrt((4 t + 4)  + (4 t - 1)  + 4 t ) sqrt((4 (4 t + 4) - 4 (4 t - 1))

		      2			     2
 + (8 t - 2 (4 t + 4))  + (2 (4 t - 1) - 8 t) ))]
radcan(nm);
					       2
	(42 sqrt(13) t + 40 sqrt(13)) sqrt(36 t  + 24 t + 17)
(%o14) [-----------------------------------------------------, 
			      2
			1404 t  + 936 t + 663

				       2
(12 sqrt(13) t + 17 sqrt(13)) sqrt(36 t  + 24 t + 17)
-----------------------------------------------------, 
		      2
		1404 t  + 936 t + 663

					 2
  (48 sqrt(13) t - 10 sqrt(13)) sqrt(36 t  + 24 t + 17)
- -----------------------------------------------------]
			2
		  1404 t  + 936 t + 663

Now let us try to improve the definitions of tangent, normal, binormal:

We may assume that for expressions that contain trigonometric functions both the dot product and the cross product will frequently yield expressions that contain powers of trigonometric functions. It is therefore a good idea to apply trigsimp to these results.

Here are the new definitions:

tangent(fn, x) :=
       diff(fn(x), x) / radcan (trigsimp((diff(fn(x), x) . diff(fn (x), x)))^(1/2))$
binormal (fn, x) := 
         trigsimp(cross (diff(fn(x), x), 
                         diff (fn(x), x, 2)
                         )
                    ) / 
            radcan
                 (trigsimp(
                     ( cross(diff(fn(x), x), diff (fn(x), x, 2))
                        . cross(diff(fn(x), x), diff (fn(x), x, 2))
                      ) ^(1/2)
               )       )$
normal (fn, x) := radcan (trigsimp(cross(binormal(fn, x), tangent(fn, x))))$

You find these definitions in tubeplot2.mc.

With these definitions, we obtain reasonably well simplified results.

What we have Learned:



Previous
Contents
Next