Calculus III - Mathematica Lab 1

Marcus McGuff
Austin Community College

(Note:  click on the little arrows to “open” each section.)

Surfaces

Cylindrical and Spherical coordinates

Mathematica has commands built in for graphing functions in cylindrical and spherical coordinates:  CylindricalPlot3D and SphericalPlot3D.  These will plot z as a function of r and θ, and ρ as a function of θ and ϕ, respectively.  

The problem comes when you need to graph something like θ==(π/3)or ϕ==(π/6), neither of which are functions of z or ρ.  One possible solution to this is to make the following definitions:

cylindricalToRect[r_,θ_,z_]:={r Cos[θ],r Sin[θ],z}

sphericalToRect[ρ_,θ_,ϕ_]:={ρ Cos[θ]Sin[ϕ],ρ Sin[θ]Sin[ϕ],ρ Cos[ϕ]}

In order to use these to graph something, you must be able to eliminate one of the domain variables.  Then use ParametricPlot3D to graph.  So, for example, to graph the two equations mentioned above, you could execute:

ParametricPlot3D[cylindricalToRect[r,(π/3),z],{r,-5,5},{z,-5,5}]

ParametricPlot3D[sphericalToRect[ρ,θ,(π/6)],{ρ,0,5},{θ,0,2π}]

Notice that SphericalPlot3D[f[θ,ϕ],...] is the same as ParametricPlot3D[sphericalToRect[f[θ,ϕ],θ,ϕ],...] (and similarly for CylindricalPlot3D).  This method allows you to graph a wide variety of graphs in spherical and cylindrical coordinates.

Use the functions developed above to graph the following over domains appropriate to the coordinate system:

1. r=z cos(2θ)

2. ϕ=ρ sin((θ/2))

How does the shape of the surface r=sin(n θ)cos(m z) depend on the parameters m and n?  Clearly state your conclusions, illustrating them with examples.

Hint:  This is in cylindrical coordinates, but you can't easily use CylindricalPlot3D to graph it (because r is a function of θ and z, rather than z being a function of θ and r).

How does the shape of the surface ρ=1+0.2sin(m θ)sin(n ϕ) depend on the parameters m and n?  Clearly state your conclusions, illustrating them with examples.

(You CAN use the SphericalPlot3D function here if you wish.  Be careful to specify ϕ first in your range.)  This particular graphs have been suggested as models for tumors (they are often called "bumpy spheres").

Optimization

Finding critical points using Solve

First, we need to do the following to make life easier:

SetAttributes[Equal,Listable]

This makes it easier to set vector equations equal to each other and solve them.

To find the critical points of a function, you need to find all points where the gradient is equal to the zero vector or undefined.

f[x_,y_]:=(x^3)-2x+3(y^2)-5y

∇f=D[f[x,y],{{x,y}}]

Notice that the gradient in Mathematica is found using the basic derivative function D, but you have to give the list of independent variables in a double curly bracket {{}}.  (Be sure you remember this, or you won’t get what you planned on...)  This is never undefined, so we just need to set each component equal to 0 and solve the system of equations:

critPts=Solve[∇f=={0,0},{x,y}]

There are two critical points for this function.  (Notice the compact way we do this.  You could type in the equations individually, but this is much nicer.)  Also, we give the set of critical points a name (critPts) so we can refer back to them later and use them to calculate with.  You can get Mathematica to pick out a single one of these points using the following notation:

critPts[[1]]

critPts[[2]]

(Yes, you have to use the double square brackets like this.)

Of course, now you have to decide whether each of these is a local minimum, maximum, saddle point, or none of these.  One way to do this would be to graph the function around the critical points and see what it looks like:

Plot3D[f[x,y],{x,-2,2},{y,-2,2}]

Perhaps you can tell what we have here.  If you were to "zoom in" on each critical point more closely (i.e., change your domain displayed), you could probably tell better.  However, it is often easier to look at the contour plot instead:

ContourPlot[f[x,y],{x,-2,2},{y,-2,2},Contours->50]

Notice that this makes it much clearer what is going on.  (It is often useful to use the Contours->50 or some large number so all minima/maxima get completely encircled.)  So, the critical point with the positive x value is a minimum, while the one with the negative x value is a saddle point.  (Be certain that you understand how to tell this from the graph before continuing.)

Let's use the 2 variable second derivative test to check this for sure.  In doing this, I will demonstrate a few tricks you can do with Mathematica to streamline this process (cutting and pasting is a really Bad Idea).

We could do this a couple of different ways.  First, lets find the matrix of 
second derivatives and define D from the second derivative test:

secondDerivTest={secondDMatrix=(({{(∂_(x,x))f[x,y], (∂_(x,y))f[x,y]}, {(∂_(y,x))f[x,y], (∂_(y,y))f[x,y]}})),Det[secondDMatrix]}

The first two sets of expressions give us the rows of the matrix, while the last expression is the value of the determinant.

Now, to plug the critical points in, we could cut and paste them (bad idea), or we could use the replacement operator /.

secondDerivTest /.{y->(5/6),x->(Sqrt[(2/3)])}

This replaces x with (Sqrt[(2/3)]) and y with (5/6).  These are called replacement rules. However, if you notice, our list of critical points is already composed of a list of these replacement rules:

critPts

So, we could do the following:

secondDerivTest /.critPts

The first set of values corresponds to the values at the first critical point and the second set to the second point.  So, at the first point, D is negative (-12 (Sqrt[6])), which means it is a saddle point.  At the second point, D is positive (12 (Sqrt[6])) and (Subscript[f, x,x]) is positive (2 (Sqrt[6])), so the point is a local minimum.

Finding critical points using FindRoot/FindMinimum

Unfortunately, Solve doesn't always work to find the critical points (if your function is exponential or trigonometric, for example).  In this case, there are two options.

Use FindRoot

The problem here is that, since it is based on Newton's method, you need an initial guess.  This means two things:  First, you need to look at a graph of the function (contour lines are good here) to make a reasonable guess.  Second, you can only realistically find critical points in a particular region (not true if Solve works; you can find all critical points then).  After finding each critical point, you still have to classify it (min/max/saddle point/neither).

Use FindMinimum/FindMaximum

This command is kind of like the FindRoot command, but you don’t have to give it an initial guess); it then automatically finds a nearby minimum of the function (you CAN give it an initial guess, in order to find other minima/maxima).  This has the advantage that you know your answer will be a minimum/maximum.  Of course, it won't help you find saddle points.

Use NMinimize/NMaximize

This is similar to FindMinimum/FindMaximum, excet it looks for an absolute min/max, not local ones.  (Minimize/Maximize does the same thing symbolicaly.)

Example

Consider the function:

f[x_,y_]:=(2y-x+(E^(1-(x^2)-(y^2))))/(1+(x^2)+(y^2))

∇f=D[f[x,y],{{x,y}}]

Solve[∇f=={0,0},{x,y}]

In other words, forget it.  So, let's get a quick look at this near the origin:

ContourPlot[f[x,y],{x,-2,2},{y,-2,2},Contours->20]

It appears there is a minimum around (0.8,-1.5) and a maximum around (0,0.2) or so.  Let's find them:

FindRoot[∇f,{x,0.8},{y,-1.5}]

FindRoot[∇f,{x,0},{y,0.2}]

You can use the graph (or the second derivative test) to see what each of these is.  If there were a saddle point, you could find it this way as well.

Or, using the FindMinimum command:

FindMinimum[f[x,y],{x,0},{y,0}]

Notice that you use this command on f(x,y), not on ∇f  (and it gives you the value of the function at the point as an added bonus).  Also, since there is only one minimum nearby, you don't have to be very accurate with your initial guess.

In order to find the maximum:

FindMaximum[f[x,y],{x,0},{y,0}]

You would not be able to find saddle points reliably this way, however.  (Sometimes it will find them, but not always.)

Lagrange multipliers

Fortunately, the method of Lagrange multipliers can be relatively easy to implement using Mathematica.  For example, if you need to optimize f(x,y)=3(x^2)(y^2) subject to the contraint that (x^2)+(y^2)=4, you could set it up as follows:

f[x_,y_]:=3(x^2)(y^2)

∇f=D[f[x,y],{{x,y}}]

g[x_,y_]:=(x^2)+(y^2)

∇g=D[g[x,y],{{x,y}}]

system={∇f==λ ∇g,g[x,y]==4}//Flatten

The //Flatten command at the end is necessary here for the Solve command to work properly with this.  It strips out the extra {}'s you would get otherwise.  Now you can use the Solve command to solve the system (even though you don't care what λ equals, you should still solve for it):

solutionPts=Solve[system,{x,y,λ}]

Notice that, again, we can give a name to the entire set of solutions to this sytem.  This way, we can work with all of them at one time, or pick a single one out to work with. 

Now, we need to find which of these is the minimum and which is the maximum:

f[x,y]/.solutionPts

(This plugs all of the solution points into f and displays the results.)  So, the minimum value is 0 (which occurs at the first 4 solution points) and the maximum value is 12 (which occurs at the last 4 solution points).

Optimization Problems:

For each of the following, find the critical points, show them graphically on a contour plot (and 3-D if you wish), and classify them using the second derivative test as minima/maxima/saddle points/none of these.  If possible, use Solve to find all critical points.  If Solve doesn't work (or you don't get an answer back from it in a few minutes of waiting), then use FindRoot or FindMinimum to find all critical points on the region -2≤x≤2 and -2≤y≤2.

1.

f[x_,y_]:=((y^2)-y+x y+Log[(x^4)+(y^2)+1])(E^(-((x^2)+(y^2))))

2.

f[x_,y_]:=5y-(y^3)-2((x )^2)y+3x

3. Find the minimum and maximum values of the function f subject to the given restriction.

f[x_,y_,z_,w_]:=4(y^2)-5(x^2)+2(w^2)-3(z^2)

subject to (x^2)+(y^2)+3(z^2)+4(w^2)=10