Sage Quickstart for Multivariable Calculus¶
This Sage quickstart tutorial was developed for the MAA PREP Workshop “Sage: Using Open-Source Mathematics Software with Undergraduates” (funding provided by NSF DUE 0817071). It is licensed under the Creative Commons Attribution-ShareAlike 3.0 license (CC BY-SA).
Because much related material was covered in the calculus tutorial, this quickstart is designed to:
Give concrete examples of some things not covered there, without huge amounts of commentary, and
Remind the reader of a few of the things in that tutorial.
Vector Calculus¶
In Sage, vectors are primarily linear algebra objects, but they are slowly becoming simultaneously analytic continuous functions.
Dot Product, Cross Product¶
sage: v = vector(RR, [1.2, 3.5, 4.6])
sage: w = vector(RR, [1.7,-2.3,5.2])
sage: v*w
17.9100000000000
sage: v.cross_product(w)
(28.7800000000000, 1.58000000000000, -8.71000000000000)
Lines and Planes¶
Intersect \(x-2y-7z=6\) with
\(\frac{x-3}{2}=\frac{y+4}{-3}=\frac{z-1}{1}\). One way to plot the
equation of a line in three dimensions is with parametric_plot3d.
sage: # designed with intersection at t = 2, i.e. (7, -10, 3)
sage: var('t, x, y')
(t, x, y)
sage: line = parametric_plot3d([2*t+3, -3*t-4, t+1], (t, 0, 4),color='red')
sage: plane = plot3d((1/5)*(-12+x-2*y), (x, 4, 10), (y, -13,-7), opacity=0.5)
sage: intersect=point3d([7,-10,3],color='black',size=30)
sage: line+plane+intersect
Graphics3d Object
Vector-Valued Functions¶
We can make vector-valued functions and do the usual analysis with them.
sage: var('t')
t
sage: r=vector((2*t-4, t^2, (1/4)*t^3))
sage: r
(2*t - 4, t^2, 1/4*t^3)
sage: r(t=5)
(6, 25, 125/4)
The following makes the derivative also a vector-valued expression.
sage: velocity = r.diff(t) # velocity(t) = list(r.diff(t)) also would work
sage: velocity
(2, 2*t, 3/4*t^2)
Currently, this expression does not function as a function, so we need to substitute explicitly.
sage: velocity(t=1)
(2, 2, 3/4)
sage: T=velocity/velocity.norm()
sage: T(t=1).n()
(0.683486126173409, 0.683486126173409, 0.256307297315028)
Here we compute the arclength between \(t=0\) and \(t=1\) by
integrating the normalized derivative. As pointed out in the
calculus tutorial, the syntax for
numerical_integral is slightly nonstandard – we just put in the
endpoints, not the variable.
sage: arc_length = numerical_integral(velocity.norm(), 0,1)
sage: arc_length
(2.3169847271197814, 2.572369791753217e-14)
We can also plot vector fields, even in three dimensions.
sage: x,y,z=var('x y z')
sage: plot_vector_field3d((x*cos(z),-y*cos(z),sin(z)), (x,0,pi), (y,0,pi), (z,0,pi),colors=['red','green','blue'])
Graphics3d Object
If we know a little vector calculus, we can also do line integrals. Here, based on an example by Ben Woodruff of BYU-Idaho, we compute a number of quantities in a physical three-dimensional setting.
Try to read through the entire code. We make an auxiliary function that will calculate \(\int (\text{integrand})\, dt\). We’ll use our formulas relating various differentials to \(dt\) to easily specify the integrands.
sage: density(x,y,z)=x^2+z
sage: r=vector([3*cos(t), 3*sin(t),4*t])
sage: tstart=0
sage: tend=2*pi
sage: t = var('t')
sage: t_range=(t,tstart,tend)
sage: def line_integral(integrand):
....: return RR(numerical_integral((integrand).subs(x=r[0], y=r[1],z=r[2]), tstart, tend)[0])
sage: _ = var('x,y,z,t')
sage: r_prime = diff(r,t)
sage: ds=diff(r,t).norm()
sage: s=line_integral(ds)
sage: centroid_x=line_integral(x*ds)/s
sage: centroid_y=line_integral(y*ds)/s
sage: centroid_z=line_integral(z*ds)/s
sage: dm=density(x,y,z)*ds
sage: m=line_integral(dm)
sage: avg_density = m/s
sage: moment_about_yz_plane=line_integral(x*dm)
sage: moment_about_xz_plane=line_integral(y*dm)
sage: moment_about_xy_plane=line_integral(z*dm)
sage: center_mass_x = moment_about_yz_plane/m
sage: center_mass_y = moment_about_xz_plane/m
sage: center_mass_z = moment_about_xy_plane/m
sage: Ix=line_integral((y^2+z^2)*dm)
sage: Iy=line_integral((x^2+z^2)*dm)
sage: Iz=line_integral((x^2+y^2)*dm)
sage: Rx = sqrt(Ix/m)
sage: Ry = sqrt(Iy/m)
sage: Rz = sqrt(Iz/m)
Finally, we can display everything in a nice table.
Recall that we use the r"stuff" syntax to indicate “raw” strings
so that backslashes from LaTeX won’t cause trouble.
sage: html.table([
....: [r"Density $\delta(x,y)$", density],
....: [r"Curve $\vec r(t)$",r],
....: [r"$t$ range", t_range],
....: [r"$\vec r'(t)$", r_prime],
....: [r"$ds$, a little bit of arclength", ds],
....: [r"$s$ - arclength", s],
....: [r"Centroid (constant density) $\left(\frac{1}{m}\int x\,ds,\frac{1}{m}\int y\,ds, \frac{1}{m}\int z\,ds\right)$", (centroid_x,centroid_y,centroid_z)],
....: [r"$dm=\delta ds$ - a little bit of mass", dm],
....: [r"$m=\int \delta ds$ - mass", m],
....: [r"average density $\frac{1}{m}\int ds$" , avg_density.n()],
....: [r"$M_{yz}=\int x dm$ - moment about $yz$ plane", moment_about_yz_plane],
....: [r"$M_{xz}=\int y dm$ - moment about $xz$ plane", moment_about_xz_plane],
....: [r"$M_{xy}=\int z dm$ - moment about $xy$ plane", moment_about_xy_plane],
....: [r"Center of mass $\left(\frac1m \int xdm, \frac1m \int ydm, \frac1m \int z dm\right)$", (center_mass_x, center_mass_y, center_mass_z)],
....: [r"$I_x = \int (y^2+z^2) dm$", Ix],[r"$I_y=\int (x^2+z^2) dm$", Iy],[mp(r"$I_z=\int (x^2+y^2)dm$"), Iz],
....: [r"$R_x=\sqrt{I_x/m}$", Rx],[mp(r"$R_y=\sqrt{I_y/m}"), Ry],[mp(r"$R_z=\sqrt{I_z/m}"),Rz]
....: ])
Functions of Several Variables¶
This connects directly to other issues of multivariable functions.
How to view these was mostly addressed in the various plotting tutorials. Here is a reminder of what can be done.
sage: # import matplotlib.cm; matplotlib.cm.datad.keys()
sage: # 'Spectral', 'summer', 'blues'
sage: g(x,y)=e^-x*sin(y)
sage: contour_plot(g, (x, -2, 2), (y, -4*pi, 4*pi), cmap = 'Blues', contours=10, colorbar=True)
Graphics object consisting of 1 graphics primitive
Partial Differentiation¶
The following exercise is from Hass, Weir, and Thomas, University Calculus, Exercise 12.7.35. This function has a local minimum at \((4,-2)\).
sage: f(x, y) = x^2 + x*y + y^2 - 6*x + 2
Quiz: Why did we not need to declare the variables in this case?
sage: fx(x,y)= f.diff(x)
sage: fy(x,y) = f.diff(y)
sage: fx; fy
(x, y) |--> 2*x + y - 6
(x, y) |--> x + 2*y
sage: f.gradient()
(x, y) |--> (2*x + y - 6, x + 2*y)
sage: solve([fx==0, fy==0], (x, y))
[[x == 4, y == -2]]
sage: H = f.hessian()
sage: H(x,y)
[2 1]
[1 2]
And of course if the Hessian has positive determinant and \(f_{xx}\) is positive, we have a local minimum.
sage: html("$f_{xx}=%s$"%H(4,-2)[0,0])
sage: html("$D=%s$"%H(4,-2).det())
Notice how we were able to use many things we’ve done up to now to solve this.
Matrices
Symbolic functions
Solving
Differential calculus
Special formatting commands
And, below, plotting!
sage: plot3d(f,(x,-5,5),(y,-5,5))+point((4,-2,f(4,-2)),color='red',size=20)
Graphics3d Object
Multiple Integrals and More¶
Naturally, there is lots more that one can do.
sage: f(x,y)=x^2*y
sage: # integrate in the order dy dx
sage: f(x,y).integrate(y,0,4*x).integrate(x,0,3)
1944/5
sage: # another way to integrate, and in the opposite order too
sage: integrate( integrate(f(x,y), (x, y/4, 3)), (y, 0, 12) )
1944/5
sage: var('u v')
(u, v)
sage: surface = plot3d(f(x,y), (x, 0, 3.2), (y, 0, 12.3), color = 'blue', opacity=0.3)
sage: domain = parametric_plot3d([3*u, 4*(3*u)*v,0], (u, 0, 1), (v, 0,1), color = 'green', opacity = 0.75)
sage: image = parametric_plot3d([3*u, 4*(3*u)*v, f(3*u, 12*u*v)], (u, 0, 1), (v, 0,1), color = 'green', opacity = 1.00)
sage: surface+domain+image
Graphics3d Object
Quiz: why did we need to declare variables this time?