Arc Length, TNB-frames, Curvature, and Torsion 

To execute all commands select "EditExecuteWorksheet". Be warned: This takes about 2-3 minutes to execute on my office desktop. 

Maple worksheet: curvature.mw 

As usual we begin by wiping memory and reloading the relevant packages. We also let Maple know that the variable "t" is a real number for the rest of this worksheet. 

>	restart; with(VectorCalculus): with(plots): assume(t::real);

 

Example: Recall that we can parametrize the ellipse using sines and cosines as follows:  , , (and ) where . Let's define a vector valued function called "" so that the graph of is this ellipse. 

>	r := t -> <2*cos(t)+3,5*sin(t)-2,0>: 'r(t)' = r(t);

 

(1)

 

Next, we define the arc length function for this parametrization and then evalute the arc length function when to compute the arc length of the ellipse. 

Recall that the arc length function is given by  where . 

>	s := t -> int(Norm(diff(r(u),u)),u=0..t): 's(t)' = s(t); s(2*Pi) = evalf(s(2*Pi));

 

 

	(2)

 

Notice that the arc length of an ellipse involves elliptic functions -- surprise! 

Alternatively we could use the VectorCalculus command "ArcLength" (this might be the smarter thing to do rather than reinventing the wheel). 

>	s := t -> ArcLength(r(u),u=0..t): 's(t)' = s(t); s(2*Pi) = evalf(s(2*Pi));

 

 

	(3)

 

Using the "diff", "Norm", and "&x" (cross product) commands, we can construct the unit tangent, unit normal, and binormal functions for this ellipse. 

Recall that T(t) = r'(t) / |r'(t)|,  N(t) = T '(t) / |T '(t)|,  and B(t) = T(t) N(t). 

>	T := t -> simplify(1/Norm(diff(r(t),t))*diff(r(t),t)): 'T(t)' = T(t); N := t -> simplify(1/Norm(diff(T(t),t))*diff(T(t),t)): 'N(t)' = N(t); B := t -> simplify(T(t) &x N(t)): 'B(t)' = B(t);

 

 

 

	(4)

 

VectorCalculus has the built in command "TNBFrame" which does all this for us. 

>	ellipseFrame := TNBFrame(r(u),u); T := t -> simplify(subs(u=t,ellipseFrame[1])): 'T(t)' = T(t); N := t -> simplify(subs(u=t,ellipseFrame[2])): 'N(t)' = N(t); B := t -> simplify(subs(u=t,ellipseFrame[3])): 'B(t)' = B(t);

 

 

 

 

	(5)

 

Note: Our computations yielding the TNB-frame involved VectorCalculus package vectors whereas the TNBFrame command yields "arrays" (old style vectors). 

Next, let's plot the ellipse together with its TNB-frame at the point . 

>

curvePlot := spacecurve(r(t),t=0..2*Pi,color=black,thickness=3): TPlot := arrow(r(Pi/2),T(Pi/2),shape=arrow,color=blue): NPlot := arrow(r(Pi/2),N(Pi/2),shape=arrow,color=green): BPlot := arrow(r(Pi/2),B(Pi/2),shape=arrow,color=red): display({curvePlot,TPlot,NPlot,BPlot},scaling=constrained,orientation=[126,62],viewpoint=circleleft);

 

Now we define the curvature function and plot curvature. Remember that  

>	kappa := t -> simplify(subs(u=t,Norm(diff(T(u),u))/Norm(diff(r(u),u)))): 'kappa(t)' = kappa(t);

 

(6)

 

...Or using the built-in function... 

>	simplify(Curvature(r(t),t)); plot(kappa(t),t=0..2*Pi);

 

The osculating circle at has radius (one over curvature). It's center is located at [1/ units away from the curve in the unit normal direction] and it's parallel to both T(t) and N(t).  

The osculating plane is the plane through the point r(t) with normal vector is given by B(t) or alternatively it's parallel to both T(t) and N(t).  

Let's plot the curve, the osculating circle at , and the osculating plane at together. 

>

# # rCirc(t,u) sweeps out the osculating circle at t as u varies from 0 to 2*Pi. # rCirc := (t,u) -> simplify(r(t)+1/kappa(t)*N(t)+1/kappa(t)*T(t)*cos(u)+1/kappa(t)*N(t)*sin(u)): 'rCirc(t,u)'=rCirc(t,u); oscCircPlot := spacecurve(rCirc(Pi/2,u),u=0..2*Pi,color=blue,thickness=2): oscPlanePlot := plot3d(r(Pi/2)+T(Pi/2)*t1+N(Pi/2)*t2,t1=-10..10,t2=-5..20,color=yellow): display({curvePlot,oscCircPlot,oscPlanePlot},scaling=constrained,orientation=[61,52],viewpoint=circleleft);

 

Torsion is a measure of how a curve "bends" out of its osculating planes. This curve's torsion is zero because the curve is planar (it doesn't "bend" out of it's osculating plane at all). We can compute torsion using the built in function. 

>	tau := t -> simplify(subs(u=t,Torsion(r(u),u))): 'tau(t)' = tau(t);

 

(7)

 

Example:  Consider the following curve parametrized by 

>	r := t -> <t,t^2,5*sin(t)>: 'r(t)' = r(t);

 

(8)

 

and .   

First, we'll plot this curve. 

>	spacecurve(r(t),t=-2*Pi..2*Pi,numpoints=5000,scaling=constrained,color=blue,thickness=3,axes=boxed,viewpoint=circleleft);

 

Next, let's compute its arc length and plot its curvature. 

Since Maple cannot find an exact arc length, we'll use "evalf" to get a decimal approximation. 

>	"Arc Length" = evalf(ArcLength(r(t),t=-2*Pi..2*Pi)); plot(Curvature(r(t)),t=-2*Pi..2*Pi);

 

Example: Consider the twisted cubic: 

>	r := t -> <t,t^2,t^3>: 'r(t)' = r(t); spacecurve(r(t),t=-1..2,color=blue,thickness=3,axes=boxed,scaling=constrained,orientation=[10,80-5],viewpoint=circleleft);

 

Let's find the TNB-frame of and then plot the curvature and torsion for this curve.
 

>	tnbFrame := TNBFrame(r(t),t): T := tnbFrame[1]: N := tnbFrame[2]: B := tnbFrame[3]: "[T,N,B]" = [T,N,B]; plot(Curvature(r(t),t),t=-1..2,title="Curvature"); plot(Torsion(r(t),t),t=-1..2,title="Torsion");

 

Now let's plot the osculating circle and tangent line for the twisted cubic at the point along with a plot of the twisted cubic itself. 

>

# Plug t=1 into the unit Tangent and unit Normal functions. # The radius of the osculating circle is 1/kappa (the curvature). T1 := subs(t=1,T); N1 := subs(t=1,N); R := 1/subs(t=1,Curvature(r(t),t)); cubicPlot := spacecurve(r(t),t=-1..2,color=blue,thickness=3,axes=boxed,scaling=constrained,orientation=[90,75]): oscCircle := spacecurve(r(1)+R*N1+R*cos(t)*T1+R*sin(t)*N1,t=0..2*Pi,color=green,thickness=2): tangent := spacecurve(r(1)+t*T1,t=-4..6,color=red,thickness=2): display(oscCircle,cubicPlot,tangent,viewpoint=circleleft);

 

Example:  Let's begin by plotting a helix. 

>	spacecurve(<2*cos(t),2*sin(t),t>,t=-2*Pi..2*Pi,axes=boxed,scaling=constrained,viewpoint=circleleft);

 

Now let's try that again using a more "general approach". Instead of punching the parameterization into "spacecure" let's first define it separately so we have it saved for later. 

Note: Maple uses the notation "" for the unit vector "", "" for the unit vector "", and "" for the unit vector "". 

>	r := t -> <2*cos(t),2*sin(t),t>: 'r(t)' = r(t); spacecurve(r(t),t=-2*Pi..2*Pi,axes=boxed,scaling=constrained,viewpoint=circleleft);

 

It's easy to differentiate and integrate vector valued functions.  

>	'r(t)' = r(t); Diff(r(t),t) = diff(r(t),t); Int(r(t),t) = int(r(t),t);

 

 

 

	(9)

 

Let's plot the helix together with it's tangent line at . 

The derivative gives the direction of the line, and the parameterization will give a point on the line. The "subs" command allows us to substitute into . 

>	l := t -> r(Pi) + subs(u=Pi,diff(r(u),u))*t: 'l(t)' = l(t); helixPlot := spacecurve(r(t),t=-2*Pi..2*Pi,color=blue): tangentPlot := spacecurve(l(t),t=-3..3,color=red): display({helixPlot,tangentPlot},scaling=constrained,axes=boxed,viewpoint=circleleft);

 

Next, let's compute the arc length function for the helix. [Remember that our parameter ranges from -2π to 2π.] 

Note: Sometimes we have to tell Maple some extra information -- in this problem Maple should know that our dummy variable "" is a real number. 

>	assume(u::real): arcLength := s = int(Norm(diff(r(u),u)),u=-2*Pi..t);

 

(10)

 

It's interesting to note that Maple's built in ArcLength function does not yield a useful formula for arc length here.  

>	simplify(ArcLength(r(u),u=-2*Pi..t));

 

(11)

 

The above formula is far more complicated than the one we originally derived. Oh well, Maple's not perfect. 

Substituting in 2π (the end of the curve) will give us the length of our helix. 

>	simplify(subs(t=2*Pi,arcLength));

 

(12)

 

Now we can make Maple solve for t and then reparametrize our helix in terms of arc length. 

>	reparametrize := t = solve(arcLength,t); rArcLength := simplify(subs(reparametrize,r(t))): 'rArcLength' = rArcLength; spacecurve(rArcLength,s=0..4*sqrt(5)*Pi,scaling=constrained,axes=boxed,viewpoint=circleleft);

 

 

 

 

Warning, solve may be ignoring assumptions on the input variables.

 

Of course, our graph looks no different, since we have the same curve (with a different parametrization).  

Now we can compute the curvature of the helix using the definition: "". 

>	T := diff(rArcLength,s);

 

(13)

 

Notice that we just computed the derivative of and didn't bother dividing by it's length -- we don't need to since is parametrized with respect to arc length! Don't believe me? Let's check the length of "" -- the unit tangent. 

>

Norm(T);

 

(14)

 

Now we compute the curvature. 

>	kappa := Norm(diff(T,s));

 

(15)

 

We have re-discovered that the helix has constant curvature. 

Next, let's make Maple compute the TNB-frame for the helix. We can either do this "step-by-step" or using built in functions. Let's do it "step-by-step" first. 

>	T := t -> diff(r(t),t)/Norm(diff(r(t),t)): 'T(t)' = T(t); N := t -> diff(T(t),t)/Norm(diff(T(t),t)): 'N(t)' = N(t); B := t -> simplify(T(t) &x N(t)): 'B(t)' = B(t);

 

 

 

	(16)

 

Now let's use the built in function "TNBFrame". 

>	TNBFrame(r(t),t); T := t -> subs(u=t,simplify(TNBFrame(r(u),u)[1])): 'T(t)' = T(t); N := t -> subs(u=t,simplify(TNBFrame(r(u),u)[2])): 'N(t)' = N(t); B := t -> subs(u=t,simplify(TNBFrame(r(u),u)[3])): 'B(t)' = B(t);

 

 

 

 

	(17)

 

We can also use a built in command to compute curvature. 

>	kappa := simplify(Curvature(r(t),t));

 

(18)

 

Let's plot the TNB-frame for our helix at the point . 

>

helixPlot := spacecurve(r(t),t=-2*Pi..2*Pi,color=black): TPlot := arrow(r(Pi),T(Pi),shape=arrow,color=blue): NPlot := arrow(r(Pi),N(Pi),shape=arrow,color=green): BPlot := arrow(r(Pi),B(Pi),shape=arrow,color=red): display({helixPlot,TPlot,NPlot,BPlot},scaling=constrained,axes=boxed,viewpoint=circleleft);

 

A curve in doesn't neccessarily lie in a plane -- see for example the helix above! Given a particular point on the curve, the osculating plane at that point is the best try at placing the curve in a plane. In other words, near the point in question, the curve 

looks like it lies in the osculating plane. The osculating plane for at is the plane which passes through the point and has normal vector . We will plot the helix together with its osculating point at . 

>	oscPlane := B(Pi).(<x,y,z>-r(Pi))=0; oscPlot := plot3d(solve(oscPlane,z),x=-2..2,y=-2..3): pointPlot := pointplot3d(r(Pi),symbol=sphere,symbolsize=10,color=red): display({oscPlot,pointPlot,helixPlot},scaling=constrained,axes=boxed,viewpoint=circleleft);

 

We know that tangent lines give the best linear approximation of a curve. The osculating circle gives the best approximating circle. That is, the osculating circle at a particular point looks like the curve near that point. The osculating circle lies in the osculating plane, touches the curve at the point in question, and has a raduis of (the reciprical of the curvature at that point). Let's graph the osculating circle for the helix for . 

>	oscCircle := t -> cos(theta)/kappa * T(t) + sin(theta)/kappa * N(t) + (1/kappa * N(t) + r(t)): 'C(t)' = oscCircle(t); oscPlot := spacecurve(oscCircle(Pi),theta=0..2*Pi,color=blue): display({oscPlot,helixPlot,pointPlot},scaling=constrained,axes=boxed,viewpoint=circleleft);

 

The curvature function measures how curved our curve is. Another function called the "torsion" function measures how much the curve "twists" out of its osculating plane. Torsion is defined to be the function where . So to compute torsion directly we would need to parametrize our curve with respect to arc length, compute the curve's TNB-frame, and then differentiate the binormal vector. This is too much work. So in the same way we found a direct formula for curvature, we get a direct formula for torsion: 

Let's compute the torsion of the helix using this formula. 

>	tau := t -> simplify(((diff(r(t),t) &x diff(r(t),t$2)).diff(r(t),t$3))/(Norm(diff(r(t),t) &x diff(r(t),t$2))^2)): 'tau(t)' = tau(t);

 

(19)

 

So we see that the helix not only has constant curvature, but also constant torsion. This means it is "twisted" the "same" way everywhere. Notice if the torsion is 0, then , so the binormal is constant. This in turn implies that the osculating plane is the same at all points so we have a planar curve. In fact, a curve lies in a plane if and only if .  

Note: There is a built in function which computes torsion. 

>	simplify(Torsion(r(t),t));

 

(20)

 

Example: We can parametrize the ellipse using sines and cosines as follows: , , and where . Let's define a vector valued function called "" so that the graph of is this ellipse. Then let's compute its arc length and TNB-frame. 

>

r[1] := t -> <3*cos(t),1*sin(t),0>: 'r[1](t)' = r[1](t); s := t -> int(Norm(diff(r[1](u),u)),u=0..t): 's(t)' = s(t); arcLength := simplify(s(2*Pi))=evalf(s(2*Pi)); ellipseFrame := TNBFrame(r[1](t),t): T[1] := simplify(ellipseFrame[1]); N[1] := simplify(ellipseFrame[2]); B[1] := simplify(ellipseFrame[3]);

 

 

 

 

 

 

	(21)

 

Next, let's plot the ellipse together with its TNB-frame at the point . 

>

ellipsePlot := spacecurve(r[1](t),t=0..2*Pi,color=black): TPlot := arrow(r[1](Pi/2),subs(t=Pi/2,T[1]),shape=arrow,color=blue): NPlot := arrow(r[1](Pi/2),subs(t=Pi/2,N[1]),shape=arrow,color=green): BPlot := arrow(r[1](Pi/2),subs(t=Pi/2,B[1]),shape=arrow,color=red): display({ellipsePlot,TPlot,NPlot,BPlot},scaling=constrained,axes=boxed,viewpoint=circleleft);

 

Now let's graph the curvature function for . 

>	plot(Curvature(r[1](t),t),t=0..2*Pi);

 

Next, let's compute the torsion function for the ellipse. Notice that is it zero. So the derivative of the binormal should be the zero vector. So we'll also compute that derivative. 

>	Torsion(r[1](t),t);

 

(22)

 

>	diff(B[1],t);

 

(23)

 

The torsion is zero since the ellipse's osculating plane is the same for every value of -- this means our curve is not "twisting" out of the osculating plane. In other words, the torsion is zero because our curve lies in the xy-plane. In fact, it is true that any curve lying in a plane will have zero torsion. 

Example: Consider the twisted cubic given by the following parametric equations:  , , and where . Let's compute its arc length, curvature, and TNB-frame. 

>	r := t -> <t,t^2,t^3>: 'r(t)' = r(t); s := t -> int(Norm(diff(r(u),u)),u=0..t): 's(t)' = s(t); arcLength := simplify(s(2))=evalf(s(2)); ellipseFrame := TNBFrame(r(t),t): T := simplify(ellipseFrame[1]); N := simplify(ellipseFrame[2]); B := simplify(ellipseFrame[3]);

 

 

 

 

 

 

	(24)

 

Next, let's graph the curvature and the torsion in the same plot window. From this we can see that the twisted cubic gets less curved and less twisted as . 

>	plot({Curvature(r(t),t),Torsion(r(t),t)},t=0..2);

 

Now let's graph the twisted cubic along with its osculating circle and TNB-frame at . 

>

cubicFrame := TNBFrame(r(t),t): T1 := subs(t=1,cubicFrame[1]): N1 := subs(t=1,cubicFrame[2]): B1 := subs(t=1,cubicFrame[3]): kappa1 := subs(t=1,Curvature(r(t),t)): oscCircle := cos(theta)/kappa1*T1 + sin(theta)/kappa1*N1 + 1/kappa1*N1+r(1): cubicPlot := spacecurve(r(t),t=0..2,color=black): TPlot := arrow(r(1),T1,shape=arrow,color=blue): NPlot := arrow(r(1),N1,shape=arrow,color=green): BPlot := arrow(r(1),B1,shape=arrow,color=red): oscPlot := spacecurve(oscCircle,theta=0..2*Pi,color=grey): display({cubicPlot,TPlot,NPlot,BPlot,oscPlot},scaling=constrained,axes=boxed,viewpoint=circleleft);

 

Finally, let's graph the twisted cubic along with its osculating plane at . 

>	oscPlane := B1.(<x,y,z>-r(1))=0; oscPlanePlot := plot3d(solve(oscPlane,z),x=0..2,y=0..4): cubicPlot := spacecurve(r(t),t=0..2,color=black,thickness=3): display({oscPlanePlot,cubicPlot},axes=boxed,orientation=[-150,70],viewpoint=circleleft);

 

Example: Consider the curve parametrized by 

>	r := t -> <(2*cos(7*t)+5)*cos(t),(2*cos(7*t)+5)*sin(t),2*sin(7*t)>: 'r(t)' = r(t);

 

(25)

 

where  

First, we'll plot this curve. 

>	spacecurve(r(t),t=0..2*Pi,scaling=constrained,numpoints=5000,color=blue,thickness=3,axes=boxed,viewpoint=circleleft);

 

Next, let's compute its arc length and plot its curvature. 

Since Maple cannot find an exact arc length, we'll use "evalf" to get a decimal approximation. 

>	"Arc Length" = evalf(ArcLength(r(t),t=0..2*Pi)); plot(Curvature(r(t)),t=0..2*Pi);

 

The following commands compute the TNB-frame of our curve for values of   ranging from 0 to 2π. Then the second group of commands plots these vectors along with the curve. [I will plot sets of TNB-vectors. If you would like to plot more or less, just modify the value of n.] 

>

TNBf := TNBFrame(r(t)): n := 40:  # The number of triples of vectors to plot. T := [seq(0,i=1..n)]: N := T: B := T:  # T,N,B are initially arrays of zeros. # As i = 1, 2, ..., n, we substitute in t = 2Pi/n, 4Pi/n, ... 2n Pi/n = 2Pi. for i from 1 to n do   T[i] := simplify(subs(t=2.0*Pi*(i/n),TNBf[1]));   N[i] := simplify(subs(t=2.0*Pi*(i/n),TNBf[2]));   B[i] := simplify(subs(t=2.0*Pi*(i/n),TNBf[3])); end do: crv := spacecurve(r(t),t=0..2*Pi,color=gray,thickness=3): Tplot := [seq(0,i=1..n)]: Nplot := Tplot: Bplot := Tplot: for i from 1 to n do   Tplot[i] := arrow(r(2.0*Pi*(i/n)),T[i],shape=arrow,color=green):   Nplot[i] := arrow(r(2.0*Pi*(i/n)),N[i],shape=arrow,color=red):   Bplot[i] := arrow(r(2.0*Pi*(i/n)),B[i],shape=arrow,color=blue): end do: display({crv,seq(op([Tplot[i],Nplot[i],Bplot[i]]),i=1..n)},scaling=constrained,axes=boxed,viewpoint=circleleft);

 

>

TNBf := TNBFrame(r(t)): n := 100:  # The number of triples of vectors to plot. T := [seq(0,i=1..n)]: N := T: B := T:  # T,N,B are initially arrays of zeros. # As i = 1, 2, ..., n, we substitute in t = 2Pi/n, 4Pi/n, ... 2n Pi/n = 2Pi. for i from 1 to n do   T[i] := simplify(subs(t=2.0*Pi*(i/n),TNBf[1]));   N[i] := simplify(subs(t=2.0*Pi*(i/n),TNBf[2]));   B[i] := simplify(subs(t=2.0*Pi*(i/n),TNBf[3])); end do: crv := spacecurve(r(t),t=0..2*Pi,color=gray,thickness=3,numpoints=5000): Tplot := [seq(0,i=1..n)]: Nplot := Tplot: Bplot := Tplot: myFrame := Tplot: for i from 1 to n do   Tplot[i] := arrow(r(2.0*Pi*(i/n)),T[i],shape=arrow,color=green):   Nplot[i] := arrow(r(2.0*Pi*(i/n)),N[i],shape=arrow,color=red):   Bplot[i] := arrow(r(2.0*Pi*(i/n)),B[i],shape=arrow,color=blue):   myFrame[i] := display({crv,Tplot[i],Nplot[i],Bplot[i]}): end do: display(seq(myFrame[i],i=1..n),insequence=true,scaling=constrained,axes=boxed,viewpoint=circleright);

 

Example: Consider the following vector valued function: (where ) 

>	r := t -> <2*sqrt(t)*cos(t),sqrt(t)*sin(t),t>: 'r(t)' = r(t); spacecurve(r(t),t=0..4*Pi,color=blue,thickness=3,axes=boxed,viewpoint=circleleft);

 

First, let's compute r'(t), |r'(t)|, and the curve's arc length (using "evalf" to get a decimal approximation). 

>	diff(r(t),t); Norm(diff(r(t),t)); s := t -> int(Norm(diff(r(u),u)),u=0..t): 's(t)' = s(t); s(4*Pi) = evalf(s(4*Pi));

 

 

 

 

	(26)

 

...OR... 

>	evalf(ArcLength(r(t),t=0..4*Pi));

 

(27)

 

The following commands compute the unit tangents, normals, and binormals for this curve and then plots them. 

>

TNBf := TNBFrame(r(t)): n := 30:  # The number of triples of vectors to plot. T := [seq(0,i=1..n)]: N := T: B := T:  # T,N,B are initially arrays of zeros. # As i = 1, 2, ..., n, we substitute in t = 4Pi/n, 8Pi/n, ... 4n Pi/n = 4Pi. for i from 1 to n do   T[i] := simplify(subs(t=4.0*Pi*(i/n),TNBf[1]));   N[i] := simplify(subs(t=4.0*Pi*(i/n),TNBf[2]));   B[i] := simplify(subs(t=4.0*Pi*(i/n),TNBf[3])); end do: crv := spacecurve(r(t),t=0..4*Pi,color=gray,thickness=3): Tplot := [seq(0,i=1..n)]: Nplot := Tplot: Bplot := Tplot: for i from 1 to n do   Tplot[i] := arrow(r(4.0*Pi*(i/n)),T[i],shape=arrow,color=green):   Nplot[i] := arrow(r(4.0*Pi*(i/n)),N[i],shape=arrow,color=red):   Bplot[i] := arrow(r(4.0*Pi*(i/n)),B[i],shape=arrow,color=blue): end do: display({crv,seq(op([Tplot[i],Nplot[i],Bplot[i]]),i=1..n)},scaling=constrained,axes=boxed,viewpoint=circleleft);

 

Now let's find and plot the curvature of r(t).  

>	kappa := simplify(Curvature(r(t))); plot(kappa,t=0..4*Pi);

 

Example: Consider the graph of for . 

>	f := x -> ln(x): 'f(x)' = f(x); plot(f(x),x=0.1..5);

 

Let's compute the arc length of for . The most natural parametrization for this graph is and where . In vector valued form this is . 

>	Int(sqrt(1+(diff(f(x),x))^2),x=1..5) = ArcLength(<x,f(x)>,x=1..5); "Arc Length" = evalf(ArcLength(<x,f(x)>,x=1..5));

 

 

	(28)

 

Next, let's find the curvature of and plot it for . We can use Maple's "Curvature" function or our special formula for the curvature of a graph of a function. Both of these answers are equivalent, but the special formula yields a simpler answer. 

>	"Built in curvature function" = Curvature(<x,f(x)>,x); "Special 2D formula for curvature" = abs(diff(f(x),x,x))/(1+(diff(f(x),x))^2)^(3/2); kappa := Curvature(<x,f(x)>,x): plot(kappa, x=0.1..5);

 

Where is the curvature of maximized?  We should look for critical points (where the derivative of the curvature function is 0). 

>	simplify(diff(kappa,x));

 

(29)

 

Interestingly, if we attempt to solve the unsimplified form of , Maple fails to find any solutions. However, if we simplify and then solve, Maple finds 2 critical points (we are interested in the positive one). 

>	solve(diff(kappa,x)=0);

 

(30)

 

>	solve(simplify(diff(kappa,x)=0));

 

(31)

 

Alternatively, we can use "fsolve" to numerically solve (looking for an approximate answer). However, "fsolve" (using the default guess "x=0") picks the "wrong" solution. We need to give it some help and suggest that the solution we want is close to . 

>	fsolve(simplify(diff(kappa,x))=0);

 

(32)

 

>	fsolve(diff(kappa,x)=0,x=1);

 

(33)

 

The curvature is maximized (on the interval ) at about (the exact answer is /).