Using Mathematica To Teach Calculus

Devendra Kapadia

Wolfram Research, Inc.

Goals of this workshop als of this Presentation

During this workshop, we will:

Give an introduction to basic calculus operations using Mathematica

Show some interesting applications to real-world problems

Provide some insight into the technology underlying our calculus functions

Take a trip through the resources available for these functions

Use animations  and interactive graphics for understanding concepts

Convey some of the beauty and utility of  calculus through Mathematica

Overview of Calculus

The Basic Functions of Calculus

Mathematica provides functions for all the standard operations of calculus:

Infinitesimal Calculus

Derivatives			D

Integrals		           Integrate

Limits				Limit

Differential equations	DSolve

Discrete Calculus

Finite and infinite sums	    Sum

Finite and infinite products       Product

Difference equations		     RSolve

Examples from Infinitesimal Calculus

Derivatives

D can be used to find derivatives of elementary functions.

Together[D[(((x^2)+1)/((E^(Sin[x])) (3 x+5))), x]]

The following implicit function has a rather simple derivative.

ImplicitEquation=(((((y(x))/(2^(2/3)))-(x/(2^(2/3)))) Ai((((((y(x))/(2^(2/3)))-(x/(2^(2/3)))))^2)+(1/((2^(1/3)) x)))+Ai'((((((y(x))/(2^(2/3)))-(x/(2^(2/3)))))^2)+(1/((2^(1/3)) x))))/((((y(x))/(2^(2/3)))-(x/(2^(2/3)))) Bi((((((y(x))/(2^(2/3)))-(x/(2^(2/3)))))^2)+(1/((2^(1/3)) x)))+Bi'((((((y(x))/(2^(2/3)))-(x/(2^(2/3)))))^2)+(1/((2^(1/3)) x)))))==1;

Simplify[Solve[D[ImplicitEquation, x],y'[x]]]

Integrals

Integrate can compute indefinite and definite integrals. Conditions are returned wherever it is necessary.

∫x sin(x)\[DifferentialD]x

(∫_0^∞)(E^(a x))\[DifferentialD]x

Limits

This computes the limit of a sequence as n approaches infinity.

Limit[((2 (n^2)+5)/(3 (n^2)+1)),n->∞]

Left-hand and right-hand limits can be found using the Direction option for Limit.

f[x_]:=Sign[x]

Plot[f[x], {x,-1,1}, PlotStyle -> {Thickness[0.02],Magenta}, Exclusions -> {0}]

Limit[f[x], x-> 0]

Limit[f[x], x-> 0, Direction -> 1]

Differential Equations

This first-order differential equation is solved in terms of exponentials and trigonometric functions.

DifferentialEquation=y'(x) + 2 y(x)==sin(x);

Solution=DSolve[DifferentialEquation,y,x]

The solution can be verified and plotted as follows.

Simplify[DifferentialEquation/.Solution[[1]]]

Plot[y(x)/.Solution[[1]]/.{(Subscript[c, 1])->4},{x,0,13},PlotStyle->{Thickness[0.01],Red},Background->Yellow]

Examples from Discrete Calculus

Sums

The sum of the numbers in any row of Pascal’s triangle is a power of 2.

Manipulate[Pane[Text[Style[Grid[CellularAutomaton[{#[[1]]+#[[3]]&,{},1},{{1},0},12]/. 0->"",Spacings->{-2/3,0},ItemSize->2,ItemStyle->{Automatic,Join[Table[Black,{n+1}],Table[White,{13-n}]]}],15]],{500,275},Alignment->{Left,Center}],{{n,12,"rows"},1,12,1},AutorunSequencing->Automatic, ControllerLinking->True]

(∑_(k=0)^n)(({{n}, {k}}))

Products

This product defines the Factorial function.

(∏_(k=1)^n)k

The following product is evaluated in terms of special functions.

(∏_(k=1)^n)(((3 k-1))^(3 k-1)) (E^(-k))

Difference Equations

This difference equation generates the Fibonacci sequence.

RSolve[{a(n)==a(n-2)+a(n-1),a(1)==1,a(2)==1},a,n]

Table[a(n)/.%[[1]],{n,1,10}]

Symbolic Integration with Mathematica

Scope of the Integrate Function

Integrate can compute:

A wide variety of  indefinite integrals

Definite integrals over finite and infinite intervals

Multiple integrals in an arbitrary number of dimensions

Integrals of piecewise functions and distributions

Indefinite Integrals

The indefinite integral of this rational function has a rational part and a transcendental logarithmic part.

∫((28 (x^4)+70 (x^3)+83 (x^2)+49 x+11)/((7 x+4) (((2 (x^2)+3 x+1))^2)))\[DifferentialD]x

We can use D to verify the result.

D[%, x]//Together

The following elliptic integral has an algebraic part and a transcendental part involving EllipticF.

∫((z^4)/(((((z^4)-1))^2) (Sqrt[(z^4)+1])))\[DifferentialD]z

Collect[%, EllipticF[a_, b_]]

The result can be verified using D.

D[%, z]//Simplify//TraditionalForm

An Improper Definite Integral

The following improper logarithmic integral has a simple closed form answer.

(∫_0^1)(((log^3)(x))/(1-x))\[DifferentialD]x

N[%]

NIntegrate[(((Log[x])^3)/((1-x))), {x,0,1}]

Multiple Integrals

The volume of the unit ball in four dimensions can be found using Boole.

(∫_(-∞)^∞)(∫_(-∞)^∞)(∫_(-∞)^∞)(∫_(-∞)^∞)Boole[(∑_(i=1)^4)((x(i))^2)≤1]\[DifferentialD]x(1)\[DifferentialD]x(2)\[DifferentialD]x(3)\[DifferentialD]x(4)

N[%]

NIntegrate[Boole[(∑_(i=1)^4)((x[i])^2)<1],{x[1],-1,1},{x[2],-1,1},{x[3],-1,1},{x[4],-1,1}]

Exact solutions of differential equations

Scope of DSolve

DSolve can find exact solutions for the following:

Ordinary differential equations (linear, nonlinear, single, system)

Partial  differential equations  (general first-order PDEs, some higher order equations)

Differential-algebraic equations (linear constant coefficient systems)

The Logistic Equation

The logistic equation appears in the study of population dynamics.

We can solve this equation and plot the solution for different values of the parameters as follows.

sol=Quiet[DSolve[{y'(x)==y(x) (1-b y(x)),y(0)==a},y,x]]

Plot[Evaluate[y(x)/.sol/.{b->(1/27)}/.{{a->(1/13)},{a->(1/2)},{a->(1/6)},{a->(1/28)}}],{x,0,18},Background->Yellow,PlotStyle->Thickness[0.01]]

The Aging Spring

This ODE is a model for an aging spring. It has an exponential function in its coefficient and is solved by transforming it to Bessel’s equation.

AgingSpring=x''(t)+4 (E^(-(t/4))) x(t) == 0;

sol=DSolve[{AgingSpring,x(0)==1,x'(0)==2},x,t]

Plot[Evaluate[x(t)/.sol[[1]]],{t,0,30},PlotStyle->{Thickness[0.0075],Red},ImageSize->400, Background -> Yellow]

A Cornu Spiral

The following equations describe a plane curve through the origin whose curvature is equal to the parameter value 's' at every point.

PlaneCurve={x'(s)==cos(t(s)),y'(s)==sin(t(s)),t'(s)==s,x(0)==0,y(0)==0,t(0)==0};

The solution of the ODEs contains Fresnel functions.

sol=DSolve[PlaneCurve,{x,y,t},s]

A graph of the curve clearly shows that the curvature becomes large as s approaches infinity.

ParametricPlot[Evaluate[{x(s),y(s)}/.sol[[1]]],{s,-20,20},PlotStyle->Red,PlotPoints->30, Background -> Yellow]

A Linear PDE with Gaussian Function Solutions

The general solution of this PDE contains an arbitrary function C[1].

sol=DSolve[y ((∂u(x,y))/(∂x))+x ((∂u(x,y))/(∂y))==-x y u(x,y),u,{x,y}]

We can obtain Gaussian and other solutions by making a particular choice of C[1].

particularsolution=u(x,y)/.sol[[1]]/.(Subscript[c, 1])[a_]:>ExpandAll[(E^(-k a))]

Animate[Plot3D[Evaluate[particularsolution/.{k->a}],{x,-6,6},{y,-6,6},PlotRange->All,ColorFunction->(Hue[0.9 (1-#3)]&),PlotPoints->50,Mesh->False,Background->Cyan],{a,0,(1/2)},DefaultDuration->10,AnimationRepetitions->1]

A Linear Constant Coefficient DAE

This example has  a differential equation, an algebraic constraint and an initial condition.

DifferentialEquation=x'(t)-y(t)==cos(t);

AlgebraicConstraint=x(t)+y(t)==1;

InitialCondition=x((π/2))==3;

DAE={DifferentialEquation,AlgebraicConstraint,InitialCondition};

sol=DSolve[DAE,{x,y},t]

Plot[{Tooltip[x(t)/.sol,Style["x[t]",Red,18]],Tooltip[y(t)/.sol,Style["y[t]",Blue,18]],Tooltip[x(t)+y(t)/.sol,Style["x[t]+ y[t]",Magenta,18]]},{t,3,12},PlotStyle->{{Thickness[0.01],Red},{Thickness[0.01],Blue},{Thickness[0.01],Dashing[{0.05}],Magenta}},Background->Yellow,Filling->{1->{2}}]

How to create a soliton with Mathematica

Scott Russell and Solitary Waves

({{In 1834, a young Scottish engineer named John Scott Russell (1808-1882) observed 
a solitary wave on the Union Canal and conducted experiments on the properties of 
these waves. , }})

In 1895, Korteweg and de Vries derived the equation for weakly nonlinear shallow water waves such as the solitary wave observed by Scott Russell.

Solving the KdV Equation with DSolve

Special solutions of the Korteweg - deVries (KdV) equation can be found using DSolve.

KdV=((∂u(x,t))/(∂t))+(((∂^3)u(x,t))/(∂(x^3)))+6 u(x,t)((∂u(x,t))/(∂x)) ==0;

sol=DSolve[KdV,u,{x,t}]

Animating the Soliton

We can get a soliton of the type observed by Russell by choosing the parameters carefully, as in the animation below.

Animate[Plot[Evaluate[u(x,a)/.sol[[1]]/.{(Subscript[c, 1])->1,(Subscript[c, 2])->-4,(Subscript[c, 3])->(1/4)}],{x,-7,20},PlotStyle->{Thickness[0.007],Blue},PlotRange->All,Frame->True,Background->Orange,FillingStyle->GrayLevel[0.6],Filling->Axis,ImageSize->900],{a,0.1,3,0.02},DefaultDuration->20,AnimationRepetitions->1]

Minimal surfaces

Definition of a Minimal Surface

Differential geometry is the application of concepts from differential calculus to the study of curves, surfaces and higher-dimensional objects.

An important class of surfaces studied in differential geometry are minimal surfaces.

The term minimal refers to the fact that they have minimal surface area for given boundary conditions.

The Helicoid and the Catenoid

The helicoid is a basic example of a minimal surface.

ParametricPlot3D[{sinh(v) sin(u),-sinh(v) cos(u),u},{u,0,4 π},{v,-2,2},Mesh->False,ColorFunction->"Rainbow"]

Another well-known example of a minimal surface is the catenoid.

ParametricPlot3D[{cosh(v) cos(u),cosh(v) sin(u),v},{u,0,4 π},{v,-2,2},Mesh->False,ColorFunction->"Rainbow"]

Here is an isometric deformation between this pair of surfaces.

Animate[ParametricPlot3D[cos(t) {sinh(v) sin(u),-sinh(v) cos(u),u}+{cosh(v) cos(u),cosh(v) sin(u),v} sin(t),{u,0,6 π},{v,-2,2},ColorFunction->"Rainbow",PerformanceGoal->"Quality",Mesh->False,PlotStyle->Directive[Orange,Opacity[1]],Boxed->False,Axes->False],{t,0,(π/2),(π/50)},AnimationRepetitions->1]

The Weierstrass Representation

In 1866, Karl Weierstrass gave a method for generating minimal surfaces, starting with a pair of functions f and g.

For example start with the following:

f=1;

g=z;

Compute the following integrals.

Integrals={∫f (1-(g^2))\[DifferentialD]z,∫I f ((g^2)+1)\[DifferentialD]z,∫2 f g\[DifferentialD]z}

Clear[f,g, r]

Transform to polar coordinates and obtain a surface.

z(r_,ϕ_)=FullSimplify[Re(Integrals/.{z->r (E^(I ϕ))}),r>0∧Im(ϕ)==0]

The above surface is minimal.

Labeled[ParametricPlot3D[Evaluate[z(r,ϕ)],{r,0,3},{ϕ,0,2 π},Boxed->False,Axes->False,PlotRange->All,PlotPoints->40,ColorFunction->"BlueGreenYellow",MeshStyle->Red],Text[Style["Enneper's Surface",Large,Bold,Red]]]

Costa's Minimal Surface

In 1984, Costa introduced an interesting minimal surface which has no self-intersections. This surface can be represented using Weierstrassian elliptic functions.

CostaMinimalSurface(u_,v_):=With[{g=WeierstrassInvariants[N[{(1/2),(I/2)}]],e1=WeierstrassP[(1/2),WeierstrassInvariants[N[{(1/2),(I/2)}]]]},Module[{wz1=WeierstrassZeta[u+I v,g],wz2=WeierstrassZeta[u+I v-(1/2),g],wz3=WeierstrassZeta[u+I v-(I/2),g],wp=WeierstrassP[u+I v,g]},{(1/2) Re(π u+((π^2)/(4 e1))-wz1+((π (wz2-wz3))/(2 e1))),(1/2) Re(π v+((π^2)/(4 e1))-I wz1-((I π (wz2-wz3))/(2 e1))),(1/4) (Sqrt[2 π]) log(\[LeftBracketingBar]((wp-e1)/(e1+wp))\[RightBracketingBar])}]]

ParametricPlot3D[Evaluate[CostaMinimalSurface(u,v)],{u,0.0001,1-0.0001},{v,0.0001,1-0.0001},Mesh->None,Boxed->False,Axes->False,PlotStyle->Directive[Yellow,Specularity[White,20]],SphericalRegion->False, Background -> Cyan]

The Painlevé Equations

The Notion of a Critical Point

The solution of the following ODE has a branch point which moves as the value of C[1] changes.

DSolve[y'(x) + ((y(x))^3)==0,y(x),x][[1]]

This is an example of a differential equation with a movable critical point.

Painlevé's Classification of Nonlinear Second-Order ODEs

Around 1900, the French mathematician Paul Painlevé showed that  ODEs of type

y''(x)==R(x,y(x),y'(x))

with no movable critical points, where R is  a rational function, can  be reduced to fifty basic  equations.

These fifty equations are called the Painlevé equations and are divided into two groups :

Forty - four of these equations can either be solved explicitly (soluble) or have first integrals (semi - transcendental).

The remaining six equations require new special functions and are called Painlevé transcendents.

The  Painlevé transcendents are labelled as Painlevé  I - VI.

Soluble Painlevé Equations

Exact solutions for several Painlevé  equations can be found using DSolve.

For example, Painlevé  Equation 3  has a solution is given in terms of Weierstrassian elliptic functions.

DSolve[y''(x)==6 ((y(x))^2)-(1/24),y(x),x]

Rational Solutions of Painlevé II

The second Painleve transcendent is given by

PainleveII=y''(x)==2 ((y(x))^3)+x y(x)+α;

We can find special rational solutions of this equation using the family of Vorob'ev-Yablonskii polynomials Q[n, z], given by

Clear[Q]

Q(0,x_):=1

Q(1,x_):=x

Q(n_,x_):=Q(n,x)=Together[((4 (((((∂Q(n-1,x))/(∂x))))^2)+x ((Q(n-1,x))^2)-4 Q(n-1,x) (((∂^2)Q(n-1,x))/(∂(x^2))))/(Q(n-2,x)))]

Q(6,x)

The roots of these polynomials lie on approximate equilateral triangles.

Animate[sol=NRoots[Together[Q(n,x)]==0,x];ListPlot[{Re(x),Im(x)}/.{ToRules[sol]},PlotStyle->{PointSize[0.035],Red},Background->Cyan,AspectRatio->1],{n,6,18,1},DefaultDuration->10, AnimationRepetitions -> 1]

The rational solutions of Painlevé II can be constructed using the Vorob' ev-Yablonskii polynomials, by taking logarithms, differentiation, etc.

Notebook on "The Painlevé equations and Mathematica"

An elementary exposition of the Painlevé equations using Mathematica 6.0  is  available in notebook form.

The Evaluation of Integrals

Methods of Integration

Indefinite integrals are often evaluated using the Risch algorithm.

Two important techniques for evaluating definite integrals are to use the corresponding indefinite integral or to apply  the Mellin transform.

Special methods are used for improper integrals, logarithmic integrals, piecewise integrals, etc.

An Exponential Indefinite Integral

Consider the following exponential indefinite integral.

∫(((E^x) ((E^x) x-8 (E^(2 x)) x-15 (E^(3 x)) x-3 (E^(4 x)) x+3 (E^(5 x)) x+4 x-5 (E^x)-(E^(2 x))+(E^(3 x))-3))/((-3+(E^x)) (((1+(E^x)))^2) x))\[DifferentialD]x

We replace exponentials by powers of a new variable t (substitution).

(((E^x) ((E^x) x-8 (E^(2 x)) x-15 (E^(3 x)) x-3 (E^(4 x)) x+3 (E^(5 x)) x+4 x-5 (E^x)-(E^(2 x))+(E^(3 x))-3))/((-3+(E^x)) (((1+(E^x)))^2) x))/.{(E^(a_. x)):>(t^a)}

Next, we apply Apart as in the method of partial fractions.

Apart[%, t]

The Risch algorithm works by making a similar replacement and then computes the polynomial, fractional, logarithmic and transcendental part of the integral:

∫(((E^x) ((E^x) x-8 (E^(2 x)) x-15 (E^(3 x)) x-3 (E^(4 x)) x+3 (E^(5 x)) x+4 x-5 (E^x)-(E^(2 x))+(E^(3 x))-3))/((-3+(E^x)) (((1+(E^x)))^2) x))\[DifferentialD]x

D[%, x]//Together//TraditionalForm

The Internals of DSolve

Structure of DSolve Code

Many important algorithms in DSolve are written using the Mathematica programming language.

The design is modular, with separate functions for solving each subproblem.

The modular design helps in automatic selection of the best algorithm for each problem.

Also, it is easy to add new algorithms whenever they become available.

Solving a second-order ODE

The general solution of

DSolve[y''(x) + y(x)==tan(x),y,x]

is found using

DSolve`DSolveLinearHomogeneousWithConstantCoefficients[{1,0,1},3,x]

followed by

DSolve`DSolveVariationOfParameters[ {Cos[x], Sin[x]}, 2, Tan[x], x, {1, 0, 1}]

Evaluation of Infinite Series

Examples of Infinite Series

The following series defines the exponential function.

(∑_(n=0)^∞)((x^n)/(n!))

The harmonic series is divergent.

(∑_(n=1)^∞)(1/n)

Steps of Summing an Infinite Series

Infinite series are summed in the following steps:

The convergence of the series is checked using the Ratio Test, Root Test or other methods.

Depending on the summand, the answer is found in terms of elementary or special functions.

(User visible step) A simplified form of the answer or a divergence message is returned.

Evaluation of an Infinite Binomial Series

Our aim is to evaluate the following infinite series,

(∑_(n=0)^∞)(((e^(-n)) (({{2 n}, {n}})) ((((2/11)))^(n+1)))/(n+3))

The nth term of the series is given by:

t(n_):=(((E^(-n)) (({{2 n}, {n}})) ((((2/11)))^n))/(n+3))

We find the ratio of adjacent terms and simplify it.

ratio=Factor[FullSimplify[((t(n+1))/(t(n)))]]

Now rewrite the ratio.

ratio = (((n+3) (n+(1/2)))/(n+4))*(1/(n+1))*(8/(11 E));

Use the ratio to build  a hypergeometric function.

HypergeometricPFQ[{3,(1/2)},{4},(8/(11 E))] t[0]

Together[%]

This agrees with the answer given by Sum.

(∑_(n=0)^∞)(((E^(-n)) (({{2 n}, {n}})) ((((2/11)))^n))/(n+3))

N[%]

Finally, we use NSum to verify the answer.

NSum[(((({{2 n}, {n}})) ((((2/11)))^n))/((E^n) (n+3))),{n,0,∞}]

Resources for Calculus functions

Function Pages

Every function in Mathematica has its own homepage called a function page.

The function page has the definition of the function, further information about it and a collection of examples on its use.

Let us look at the Integrate function page:

?Integrate

Interesting connections between infinitesimal and discrete calculus can be found here:

?RSolve

Advanced Documentation for DSolve

The DSolve Advanced Documentation was introduced in Mathematica 5.1 and provides a comprehensive reference for all aspects of DSolve.

It includes a chapter called "Working with DSolve" with useful tips and techniques.

There is a link to the Advanced Documentation from the DSolve function page.

?DSolve

The Wolfram Demonstrations Project

With Mathematica 6, users have access to hundreds of useful demonstration notebooks for classroom instruction.

You can become a part of the Mathematica community by contributing interesting demonstrations.

Summary

DSolve, Integrate and other functions are well-equipped with techniques for solving calculus problems of all types.

The ability to create animations and interactive graphics should be very useful for understanding concepts in calculus.

We believe that the new symbolic, graphical and other features of Mathematica 6 make it an excellent tool for education at all levels.