Divergence and Curl and Their Geometric Interpretations

Setting up the combined potential due to three point potentials as an example to visualize gradients, divergences, and curls

Simple (2  D)  (1/r) potential :   ( gravity,  
          electrostatic potential have this dependency)

potential[x_ ,  y_,  xo_ ,  yo_] :=  -1/
        Sqrt[(x - xo)^2  +  (y - yo)^2]

A field source located a distance 1 south of the origin

HoleSouth[x_,  y_]  :=  
      potential[x, y,  Cos[3  Pi/2],  Sin[3 Pi/2]]

Sources located distance 1 at 30° and 150°:

HoleNorthWest[x_ ,  y_]  :=  
      potential[x, y,  Cos[Pi/6],  Sin[ Pi/6]]

HoleNorthEast[x_ ,  y_]  :=  
      potential[x, y,  Cos[ 5 Pi/6],  Sin[5 Pi/6]]

Function that returns the two dimensional (x,y) gradient field of any function declared a function of two arguments:

gradfield[
        scalarfunction_]  :=  {D[scalarfunction[x, y], x] // Simplify,  
        D[scalarfunction[x, y], y] // Simplify}

Generalizing the function to any arguments:

gradfield[scalarfunction_,  x_ ,  
        y_]  :=  {D[scalarfunction[x, y], x] // Simplify,  
        D[scalarfunction[x, y], y] // Simplify}

The sum of three potentials:

ThreeHolePotential[x_,  y_]  :=  
      HoleSouth[x, y]  +  HoleNorthWest[x, y]  +  
        HoleNorthEast[x, y]

f(x,y) visualization of the scalar potential:

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

Contour visualization of the three-hole potential

ContourPlot[ThreeHolePotential[x, y], {x, -2, 2}, {y, -2, 2}, 
      PlotPoints -> 40, 
      ColorFunction -> (Hue[1 - #*0.66] &)]

Gradient field of three-hole potential

gradthreehole  =  gradfield[ThreeHolePotential]

<< Graphics`PlotField`

PlotVectorField[gradthreehole, {x, -2, 2}, {y, -2, 2}, 
      ScaleFactor -> 0.2, ColorFunction -> (Hue[1 - #*0.66] &), 
      PlotPoints -> 21]

Function that takes a two-dimensional vector function of (x,y) as an argument and returns its divergence

divergence[{xcomp_ ,  ycomp_}]  :=  
      Simplify[D[xcomp, x]  +  D[ycomp, y]]

divgradthreehole  =  
      divergence[gradfield[ThreeHolePotential]] // Simplify

\(Plotting\ the\ divergence\ of\ the\ gradient\ \((∇\(\(·\)\((∇\ f)\)\)\ is\ the\ ``Laplacian''\ ∇\^2\ f, \ 
        sometimes\ indicated\ with\ symbol\ Δf)\)\)

Plot3D[divgradthreehole, {x, -2, 2}, {y, -2, 2}, 
      PlotPoints -> 60]

It is no surprise that many of these differential operations already exist in Mathematica packages.

Using the Vector Analysis Package

<< Calculus`VectorAnalysis`

Coordinate Systems

Converting between coordinate systems

The spherical coordinates expressed in terms of the cartesian x,y,z

CoordinatesFromCartesian[{x, y, z}, 
      Spherical[r, theta, phi]]

The cartesian coordinates expressed in terms of the spherical r θ ϕ

CoordinatesToCartesian[{r, theta, phi}, 
      Spherical[r, theta, phi]]

The equation of a line  through the origin in spherical coodinates

Simplify[
      CoordinatesFromCartesian[{a t,  b t,  c t}, 
        Spherical[r, theta, phi]], t  >  0]

An example of calculating the positions of cities in cartesian and spherical coordinates.

<< Miscellaneous`CityData`

Boston is located at latitude 42° 21' 30" N  and longitude -71°,-3',-37" W

boston  =  CityData["Boston", CityPosition]

paris  =  CityData["Paris", CityPosition]

SphericalCoordinatesofCity[
        cityname_String]  :=  
{
6378.1 , 
((2 Pi)/360) 
          ToDegrees[CityData[cityname, 
                CityPosition][[1]]], 

((2 Pi)/360) 
          ToDegrees[ CityData[cityname, 
                CityPosition][[2]]]
}

SphericalCoordinatesofCity["Boston"]

CartesianCoordinatesofCity[cityname_String]  :=  
      CoordinatesToCartesian[SphericalCoordinatesofCity[cityname], 
        Spherical[r, theta, phi]]

CartesianCoordinatesofCity["Paris"]

MinimumTunnel[city1_String, city2_String]  := 
Norm[
        CartesianCoordinatesofCity[city1] - 
          CartesianCoordinatesofCity[city2]]

MinimumTunnel["Boston", "Paris"]

Minimum travel distance between Boston and Paris for a round earth

SphericalDistance[boston, paris] // N

SpheroidalDistance[boston, paris] // N

Using Vector Derivative Functions

SimplePot[x_ ,  y_ ,  z_,  n_]  := 
      (1/(((x^2  +  y^2  +  z^2))^((n/2))))

gradsp  = Grad[SimplePot[x, y, z, 1], Cartesian[x, y, z]]

\(The\ above\ is\ equal\ to\ \(r\& -> \)\/\((\(||\)\(r\& -> \)\(||\))\)\^3\)

\(Examples : \ 
      the\ same\ function\ and\ gradient\ in\ different\ coordinate\ systems\)

SimplePot[r_,  n_]  := (1/(r^n))

gradsphere = 
      Grad[SimplePot[r, 1],  
        Spherical[r, θ,  φ]]

Grad[SimplePot[r, 1],  Cylindrical[r, θ,  z]]

Grad[SimplePot[r, 1],  
      ProlateSpheroidal[r, θ, φ]]

GradSimplePot[x_,  y_,  z_,  n_]  := 
      Evaluate[Grad[SimplePot[x, y, z, n], Cartesian[x, y, z]]]

Div[GradSimplePot[x, y, z, n], Cartesian[x, y, z]] // 
      Simplify

\(The\ above\ is\ equal\ to\ \(n \((n - 1)\) \(r\& -> \)\)\/\((\(||\)\(r\& -> \)\(||\))\)\^\(2 + n\)\)

Div[GradSimplePot[x, y, z, 1], Cartesian[x, y, z]] // 
      Simplify

\(This\ makes\ ``sense''\ because\ the\ amount\ of\ stuff\ flowing\ into\ a\ sphere\ is\ like\ the\ gradient*4  π\ r\^2\ which\ is\ independent\ of\ the\ size\ of\ the\ sphere\ for\ n = 1\)

An example of Curl

There is a very useful free software tool for solving minimal surface (and many other) variational problems called Surface Evolver by Ken Brakke.  To use Surface Evolver to greatest possible advantage, a user should be adept at using results from vector analysis.  Mathematica's Vector Analysis package is very helpful aid for developing powerful Evolver codes.  The following example is extracted from the Surface Evolver manual.

LeavingKansas[x_,  y_,  z_ ,  n_] :=  
      ((z^n)/((x^2  +  y^2) (((x^2  +  y^2  +  z^2))^((n/2))))) {y, -x,  0}

LeavingKansas[x, y, z, 3]

<< Graphics`PlotField3D`

Visualize the vector field for n=3, note that the function will be singular near the z-axis

PlotVectorField3D[
      LeavingKansas[x, y, z, 3], {x, -1, 1}, {y, -1, 
        1}, {z, - .5,  .5}, VectorHeads -> True,  
      ColorFunction -> ((Hue[#* .66]) &), PlotPoints -> 15, 
      ScaleFactor -> 0.5]

We could make the function better behaved along the z-axis by brute force:

LeavingKansasNicely[x_,  y_,  z_ ,  n_]  :=  
      Module[{CindRadsq  =  
            x^2  +  y^2}, 
CindRadsq  =  
          If[CindRadsq  ≤  (10^(-4)), (10^(-4)),  
            CindRadsq,  CindRadsq]; 
 
        ((z^n)/(CindRadsq (((CindRadsq  +  z^2))^((n/2))))) {y, -x,  0}]

PlotVectorField3D[
      LeavingKansasNicely[x, y, z, 3], {x, -1, 1}, {y, -1, 
        1}, {z, - .5,  .5},  VectorHeads -> True, 
      ColorFunction -> ((Hue[#* .66]) &), PlotPoints -> 15, 
      ScaleFactor -> 0.5]

Or simply by avoiding the axis altogether and using the symmetry of the field

PlotVectorField3D[
      LeavingKansas[x, y, z, 3], {x,  .01, 1}, {y,  .01, 1}, {z,  .01,  .5}, 
      VectorHeads -> True,  
      ColorFunction -> ((Hue[#* .66]) &), PlotPoints -> 15, 
      ScaleFactor -> 0.5]

Calculate the curl of the function using the VectorAnalysis package--note that the coordinate system is specified as cartesian. 
For the particular case of n=3:

Curl[ LeavingKansas[x, y, z, 3], Cartesian[x, y, z]] // 
      Simplify

Define a new vector function for the curl for general n

Glenda[x_,  y_ ,  z_,  n_]  :=  
      Simplify[Curl[LeavingKansas[x, y, z, n], 
          Cartesian[x, y, z]]]

Demonstrate the assertion that the curl has a fairly simple form and is sphericaly symmetric for n=1

Glenda[x, y, z, n]

Glenda[x, y, z, 1]

\(The\ above\ is\ a\ vector\ field\ that\ points\ radially\ from\ the\ origin, \ with\ a\ magnitude\ that\ falls\ off\ like\ 1/
        r\^2\)

Visualize the curl for n=3

PlotVectorField3D[
      Evaluate[Glenda[x, y, z, 3]], {x, 0,  .5}, {y, 0,  .5}, {z, 0.1,  .5}, 
      VectorHeads -> True,  
      ColorFunction -> ((Hue[#* .66]) &), 
      PlotPoints -> 7]

Demonstrate that the divergence of the curl vanishes for the above function independent of n

DivCurl  =  Div[Glenda[x, y, z, n], Cartesian[x, y, z]]

Simplify[DivCurl]