Generalizations of the Fundamental Theorem of Calculus (Part I)

Using Green's theorem in the plane to make the integration of the triangular charged surface patch simpler:

Suppose there is a uniformly charged surface  (σ≡charge/area=1) occupying an equilaterial triangle in the z=0 plane:

The electrical potential goes like  (1/r), therefore the potential of a unit charge located at (x,y,z) from a small surface patch at (ξ,η,0) is ((σ dξ dη)/r)= \(TraditionalForm\`\(dξ\ dη\)\/\@\(\((x - ξ)\)\^\(\(2\)\(\ \)\) + \((y - η)\)\^\(\(2\)\(\ \)\) + z\^\(\(2\)\(\ \)\)\)\)

r[x_ ,  y_ ,  z_ ]  :=  
        Sqrt[(x - ξ)^2  +  (y - η)^2  +  
            z^2];

Try to move the integral over the region to an integral over its boundary.  We may be able to find a closed form—and even if not the numerical integration should be much more efficient...

We will use Green's Theorem in the Plane:
∫∫(∂((Subscript[F, 2])/∂) x-∂((Subscript[F, 1])/∂) y)dxdy = ∫((Subscript[F, 1])dx + (Subscript[F, 2])dy)

Find a function, that when we take the derivative with respect to x, we get 1/r.
We will let (Subscript[F, 2]) be the function we find, and set (Subscript[F, 1])=0

\(The\ boundary\ defined\ by\ the\ upper\ right\ side\ of\ the\ triangle\ is\ \((ξ \((t)\)\  = \((1 - t)\)\/2, \ η \((t)\)\  = \(\(\@\(\(3\)\(\ \)\)\) t\)\/2)\)\ ; \ \((0\  < \ t\  < \ 1)\)\)
\(The\ other\ two\ sides\ \(are : \ \((ξ \((t)\)\  = \(-\(\(\(\ \)\(t\)\)\/2\)\)\ \ , \ η \((t)\)\  = \(\@\(\(3\)\(\ \)\)\) \((1 - t)\)\/2\ )\)\ \ \((0\  < \ t\  < \ 1)\)\ and\)\)
(ξ (t)  = ((-1)/2)  +  t  ,  η (t)  =  
          0)   (0  <  t  <  1)

To use Green's theorem in the plane  (∫_▲)∫\(\((∂F\_2\/∂x\  - ∂F\_1\/∂y)\)\ dx dy\) = (∫_△) (Subscript[F, 1]) dx + (Subscript[F, 2]) dy
Try to find a vector function ((Subscript[F, 1]), (Subscript[F, 2])) that matches the integrand over the domain.  Let (Subscript[F, 2]) = 0, then to find (Subscript[F, 1 ])which when differentiated by η gives 1/r

F1[x_, y_ ,  
        z_]  = -Integrate[
          (1/(Sqrt[(((x - ξ))^2) + (((y - η))^2) + (z^2)])), η]

The third (horizontal) boundary of the triangle patch looks like the easiest, let's see if an integral can be found over that patch:

Bottomside  =  
      F1[x, y, z] /. {ξ -> t - (1/2),  η -> 0} // 
        Simplify

BottomContribution  =  
      Integrate[Bottomside, {t, 0, 1}, 
        Assumptions -> 
          x ∈ Reals && y ∈  Reals && z > 0]

NEside  =  
      F1[x, y, z] /. {ξ -> ((1 - t)/2),  η -> (((Sqrt[3]) t)/2)} // Simplify

NWside = 
      F1[x, y, z] /. {ξ -> ((-t)/2),  η -> (((Sqrt[3]) (1 - t))/2)} // Simplify

(∫_△) (Subscript[F, 1]) dx =(∫_0^1) NEside ((-dt)/2)  +  (∫_0^1) 
            NWside ((-dt)/2)  +  (∫_0^1) Bottomside dt

integrand  = 
      Simplify[ ((-(NEside  +  NWside))/2)  + Bottomside]

PotXYZ  =  
      Integrate[integrand, {t, 0, 1}, 
        Assumptions -> 
          x ∈ Reals && y ∈  Reals && z > 0]

Pot[X_,  Y_,  Z_]  :=  
      NIntegrate[
        Evaluate[integrand /. {x -> X, y -> Y, z -> Z}], {t, 
          0, 1}]

ContourPlot[Pot[a, b, 1/10], {a, -1, 1}, {b, - .5, 1.5}, 
      Contours -> 8, 
      ColorFunction -> (Hue[1  - 0.66*#] &) ]

Surfaces

Surface of the form: (z=f(x,y))

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

Plot3D[
      Evaluate[GraphFunction[x, y]], {x, -1, 1}, {y, -1, 
        1}]

Surface of the form: (x(u,v), y(u,v), z(u,v))

SurfaceParametric[u_,  v_]  :=  {Cos[u] v,  u Cos[u + v],  
        Cos[u]}

ParametricPlot3D[
      Evaluate[SurfaceParametric[u, v]], {u, -2, 2}, {v, -2, 2}, 
      PlotRange -> {-1, 1}]

ParametricPlot3D[
      Evaluate[SurfaceParametric[u, v]], {u, -2, 2}, {v, -2, 2}, 
      ViewPoint -> {1.632,  0.625,  -0.972}]

Table[
      ParametricPlot3D[
        Evaluate[SurfaceParametric[u, v]], {u, -ep, ep}, {v, -ep, 
          ep}, PlotRange -> {{-4, 4}, {-4, 4}, {-1, 1}}, 
        PlotPoints -> {1 + Round[ep/ .125], 
            1 + Round[ep/ .125]}], {ep,  .125, 4.25,  .125}]

Surface of the form: (F(x,y,z) = constant)

<< Graphics`ContourPlot3D`

ConstFunction[x_,  y_,  z_,  p_]  :=  
      (x^2)  - p x y  +  (y^2)  +  (z^2)

ContourPlot3D[
      ConstFunction[x, y, z, 2], {x, -1, 1}, {y, -1, 1}, {z, -1, 
        1}, Contours -> {6}]

cpa  = 
      ContourPlot3D[
        (x^2)  - 4 x y  +  (y^2)  +  y x (z^2), {x, -3, 
          3}, {y, -3, 3}, {z, -3, 3}, Contours -> {0, 2, 8}, 
        PlotPoints -> {5, 7}, 
        DisplayFunction -> Identity]

cpb  =  
      ContourPlot3D[
        (x^2)  - 4 x y  +  (y^2)  +  y x (z^2), {x, -3, 
          3}, {y, -3, 3}, {z, -3, 3}, Contours -> {0, 2, 8}, 
        PlotPoints  ->  {5, 7}, 
        ContourStyle ->  {{Hue[0]}, {Hue[ .25]}, {Hue[ .5]}}, 
        Lighting -> False, DisplayFunction -> Identity]

Show[GraphicsArray[{cpa, cpb}]]

Animate[
      ContourPlot3D[
        (x^2)  - 4 x y  +  (y^2)  +  y x (z^2), {x, -3, 
          3}, {y, -3, 3}, {z, -3, 3}, Contours -> {i}, 
        PlotPoints -> {5, 7}, 
        DisplayFunction -> Identity], {i, -2, 8,  .5}]

Surface Integration

Example: Integrating an Orientation-Dependent Surface Tension over the surface of a cube and a sphere

Vector function for surface of sphere of radius R, u is longitude (u ∈ (0, 2π)) and v is latitude (v ∈ (-(π/2), (π/2)))

spheresurf[u_ ,  v_ ]  :=  
      R {Cos[v] Cos[u] ,  Cos[v] Sin[u] ,  Sin[v]}

The coordinate vectors at (u,v)

Ru[u_ ,  v_]  =  
      D[spheresurf[u, v], u] // Simplify
Rv[u_ ,  v_]  =  
      D[spheresurf[u, v], v] // Simplify

<< Calculus`VectorAnalysis`

The normal vector:

NormalVector[u_ ,  v_ ]  =  
      CrossProduct[Ru[u, v], Rv[u, v]] // Simplify

NormalMag =  
      Sqrt[DotProduct[NormalVector[u, v], NormalVector[u, v]] // Simplify] // 
        PowerExpand

The Unit Normal: No surprise here, it is just the vector for a unit sphere...

UnitNormal[u_,  v_]  =  NormalVector[u, v]/NormalMag

Now, suppose the surface energy ( (Joules/(meter^2))) depends on the local orientation of the surface. This is case for many single crystals. If the single crystal is cubic and the (100) planes have lower energy than the (111) planes,  the following model for orientation-dependent surface tension can be employed.

SurfaceTension[nvec_]  := 
      1  +  (Subscript[gamma, 111])*nvec[[1]]^2 
          nvec[[2]]^2 nvec[[3]]^2

<< Graphics`ParametricPlot3D`

Here is a plot of our model aniostropic surface energy-- the lobes at [111] indicate that those orientations have higher energy:

SphericalPlot3D[
      SurfaceTension[UnitNormal[u, v]] /. (Subscript[gamma, 111])  ->  12, {u, 0, 
        2  Pi}, {v, -Pi/2, Pi/2}]

The total energy of a sphere

SphereEnergy  = 
      Integrate[
        Integrate[
          SurfaceTension[UnitNormal[u, v]] Cos[v], {u, 0, 
            2  π}], {v, -(π/2), (π/2)}]

Compare to Cube of same volume (CubeSide^3) = ((4 π)/3)

CubeSide = ((((4 π)/3))^(1/3))

CubeEnergy  =  
      6 ( (CubeSide^2) SurfaceTension[{1, 0, 0}])

\(Find\ the\ condition\ on\ gamma\_111\ such\ that\ the\ two\ surfaces\ have\ the\ same\ energy\)

EqualEnergies  =  
      Solve[CubeEnergy == SphereEnergy,  (Subscript[gamma, 111])] // 
        Flatten

If the anisotropy factor is greater than about 25, then a cube is much lower energy

N[(Subscript[gamma, 111]) /. EqualEnergies]

This computation compared two fixed shapes of the same volume.  To find which shape of all the possible shapes, one would use the Wulff Theorem