Contour and surface integrals

Contour integrals of the scalar type

These are integrals of the type

I=(∫_C)f[r]dl

where C is a contour and dl is an infinitesimal element of the contour length. If the contour can be described by {y[x], z[x]), one can use

dl=(Sqrt[1+(((((dy)/(dx))))^2)+(((((dz)/(dx))))^2)])dx.

For more complicated contours one can rewrite the contour integral as an integral over a parameter t and use the expression of the contour in the parametric form {x[t], y[t], z[t]}. In this case the length element becomes

dl=(Sqrt[(((((dx)/(dt))))^2)+(((((dy)/(dt))))^2)+(((((dz)/(dt))))^2)])dt.

Perimeter of an ellipse

As an example, let us calculate the perimeter of the ellipse

((x^2)/(a^2))+((y^2)/(b^2))=1,

Here, for the upper part of the ellipse, one has

y=b(Sqrt[1-((x^2)/(a^2))]),   ((dy)/(dx))=b((-x/(a^2))/(Sqrt[1-((x^2)/(a^2))]))

Thus the perimeter is given by

P=2(∫_(-a)^a)(Sqrt[1+(((((dy)/(dx))))^2)])dx=2(∫_(-a)^a)(Sqrt[1+(((b^2)(x^2)/(a^4))/(1-((x^2)/(a^2))))])dx=2(∫_(-a)^a)(Sqrt[(((a^2)-(1-((b^2)/(a^2)))(x^2))/((a^2)-(x^2)))])dx.

In the case of a circle, a = b = R, one easily obtains P = 2πR. However, in the general case the integral cannot be expressed via elementary functions.

2Integrate[(Sqrt[(((a^2)+(-1+((b^2)/(a^2)))(x^2))/((a^2)-(x^2)))]),{x,-a,a},Assumptions->{0<a<b}]
2Integrate[(Sqrt[(((a^2)+(-1+((b^2)/(a^2)))(x^2))/((a^2)-(x^2)))]),{x,-a,a},Assumptions->{0<b<a}]

The solution is one of elliptic integrals that received their name from the ellipse. For the circle one has

P=4 a EllipticE[0]

In the case of b -> 0 one has

P=4 a EllipticE[1]

that is an expected result as well. A similar result arises in the case a -> 0

P=4 Limit[a EllipticE[1-((b^2)/(a^2))],a->0]

Let us check that the solution is symetric in a and b

FullSimplify[a EllipticE[1-((b^2)/(a^2))]==b EllipticE[1-((a^2)/(b^2))],Assumptions->{0<b<a}]

This is too difficult for Mathematica that seemingly does not know relations between elliptic integrals. Let us define the function

fab[a_,b_]:=a EllipticE[1-((b^2)/(a^2))]-b EllipticE[1-((a^2)/(b^2))]

This function is zero for any set of arguments:

fab[1,2]//FullSimplify
fab[1,2.]

that is zero up to the Mathematica's standard precision. Let us now plot P/(4a) that has particular cases

(P/(4a))={({{1,, b<<a}, {π/2,, b=a}, {b/a,, b>>a}})

as a function of rba ≡ b/a.

Plot[{rba,EllipticE[1-(rba^2)],π/2},{rba,0,4},PlotStyle->{Dashed,Thick,Dotted}]

The ratio of the area of the ellipse S = πab to the a square of its perimeter is

(S/(P^2))=((π )/16)(b/a)(1/((EllipticE[1-((b^2)/(a^2))])^2))

This means that we have expressed the area through the perimeter and the ratio rba ≡ b/a, as two independent variables:

S=((π )/16)(b/a)((P^2)/((EllipticE[1-((b^2)/(a^2))])^2))

Let us plot this function vs b/a for a fixed P = 1 to see that its maximum corresponds to b/a = 1

Plot[((π )/16)rba((P^2)/((EllipticE[1-(rba^2)])^2))/.P->1,{rba,0,5}]

Center of mass of a parabola segment

Let us calculate the position of the center of mass of the parabola segment

y=((k (x^2))/2),  -a<x<a

The vertical coordinate of the center of mass is given by the contour integral

(Subscript[y, CM])=(1/L)(∫_C)ydl=(1/L)(∫_C)y(Sqrt[1+(((((dy)/(dx))))^2)])dx=(k/L)(∫_(-a)^a)((x^2)/2)(Sqrt[1+(k^2)(x^2)])dx

where

L=(∫_C)dl=(∫_(-a)^a)(Sqrt[1+(k^2)(x^2)])dx

is the length of the parabola segment. Calculation of the integrals yields

L=Integrate[(Sqrt[1+(k^2)(x^2)]),{x,-a,a},Assumptions->{a>0,k>0}]

and

FullSimplify[yCM=(k/L)Integrate[((x^2)/2)(Sqrt[1+(k^2)(x^2)]),{x,-a,a},Assumptions->{a>0,k>0}],
Assumptions->k>0]

Particular value:

yCM=yCM/.{a->1,k->1.}

Show[
Plot[((x^2)/2),{x,-1,1},AspectRatio->0.25],
Graphics[Disk[{0,yCM},0.015]]

]

Length of a cycloide

Find the length of a cycloid

x[ϕ]=Sin[ϕ]+a ϕ;                          y[ϕ]=Cos[ϕ];

in the interval a∈{0,2π}. Here is a plot of the cycloid.

x[ϕ_]:=Sin[ϕ]+a ϕ;                          y[ϕ_]:=Cos[ϕ];
ParametricPlot[{x[ϕ]/.a->1/2,y[ϕ]},{ϕ,0,4π}]  (* Cycloide *)

One can see that for a < 1 the dependence y[x] is non-monotonic. In this case parametric representation is the only convenient way to calculate the contour integral.

Simplify[((((∂_ϕ) x[ϕ]))^2)+((((∂_ϕ) y[ϕ]))^2)]

For a < 1 the length of the cycloid is given by

2Integrate[(Sqrt[((((∂_ϕ) x[ϕ]))^2)+((((∂_ϕ) y[ϕ]))^2)]),{ϕ,0,π},Assumptions->{a∈Reals,a≠1,a≠0}]

L1[a_]=2Integrate[(Sqrt[((((∂_ϕ) x[ϕ]))^2)+((((∂_ϕ) y[ϕ]))^2)]),{ϕ,0,π},Assumptions->0<a<1]

and for a > 1 it is given by

L2[a_]=2Integrate[(Sqrt[((((∂_ϕ) x[ϕ]))^2)+((((∂_ϕ) y[ϕ]))^2)]),{ϕ,0,π},Assumptions->1<a]

The full dependence for a > 0 can be represented by

L[a_]:=If[0<a<1,L1[a],L2[a]]
Plot[((L[a])/(2π)),{a,0,2},PlotRange->{0,2}]

L[a] has a singularity at a = 1. that becomes better visible if we plot the derivative of L[a]

L1Der[a_]=(∂_a)L1[a];  (* Immediate set is mandatory! *)
L2Der[a_]=(∂_a)L2[a];
LDer[a_]:=If[0<a<1,L1Der[a],L2Der[a]]
Plot[LDer[a],{a,0,2}]

Note that LDer cannot be determined as

LDer[a_]:=If[0<a<1,(∂_a)L1[a],(∂_a)L2[a]]
LDer[2]   (* with delayed set, Mathematica tries to calculate derivatives when they are needed, i.e., when a is set to 2. But it cannot differentiate with respect to 2 *)

or as

LDer[a_]=If[0<a<1,(∂_a)L1[a],(∂_a)L2[a]] (* This remains a delayed set because it depends on a condition, and delayed set does not work *)
LDer[2]

Contour integrals of the work / circulation type

An example of a contour integral is work along the trajectory C that is a line in 3d space

W=(∫_C)F·dr

Such integrals can be introduced in a usual way by splitting of the trajectory into a sum of many small finite displacements (Subscript[Δr, n]) and taking the limit of the sum

W=(lim_(N->∞))(∑_(n=1)^N)(Subscript[F, n])·(Subscript[Δr, n])

Practically contour integrals can be calculated in the parametric form as ordinary integrals over a parameter. The natural parameter for the mechanical integral above is time t. Changing the variable to , one obtains

W=(∫_(Subscript[t, 1])^(Subscript[t, 2]))F·((dr)/(dt))dt=(∫_(Subscript[t, 1])^(Subscript[t, 2]))F·vdt

that is, the work is the integral of the power P = F·v over time. As an example, consider the work of the gravity force

F=mg=-(Subscript[mge, z])

along a contour being a part of a circle

z=R Cos[ϕ],     x=R Sin[ϕ],      (Subscript[ϕ, 1])≤ϕ≤(Subscript[ϕ, 2])

to compute the work, one could consider motion of the particle along the circular trajectory with ϕ = ϕ[t] and write

W=(∫_(Subscript[t, 1])^(Subscript[t, 2]))((Subscript[F, x])(Subscript[v, x])+(Subscript[F, y])(Subscript[v, y])+(Subscript[F, z])(Subscript[v, z]))dt=(∫_(Subscript[t, 1])^(Subscript[t, 2]))-m g ((dz)/(dt))dt=(∫_(Subscript[t, 1])^(Subscript[t, 2]))-m g ((dz)/(dϕ))((dϕ)/(dt))dt=(∫_(Subscript[t, 1])^(Subscript[t, 2]))m g R Sin[ϕ]((dϕ)/(dt))dt.

The result does not depend on the angular velocity dϕ/dt with which the particle moves along the trjectory, because one can change the integration variable to ϕ

W=(∫_(Subscript[ϕ, 1])^(Subscript[ϕ, 2]))m g R Sin[ϕ]dϕ.

In fact, ϕ could be used as the parameter for the contour integration instead of t from the very beginning. Moreover, one can firther change the integration variable to z and finally calculate the integral

W=-(∫_(Subscript[ϕ, 1])^(Subscript[ϕ, 2]))m g dz=-m g((Subscript[z, 2])-(Subscript[z, 1]))=-m g R(Cos[(Subscript[ϕ, 2])]-Cos[(Subscript[ϕ, 1])]).

Of course, this was a toy example because the gravity force is a potential force and its work can be related to the difference of potential energies,

W=-((Subscript[U, 2])-(Subscript[U, 1]))

This relation can be proven in a general form. If the force is potential, F = -∇U(r), then the infinitesimal work is simply the negative differential of the potential energy:

dW=F·dr= -∇U(r)·dr=-(((∂U)/(∂x))(Subscript[e, x])+((∂U)/(∂y))(Subscript[e, y])+((∂U)/(∂z))(Subscript[e, z]))·dr=-(((∂U)/(∂x))dx+((∂U)/(∂y))dy+((∂U)/(∂z))dz)=-d U

Thus

W=(∫_C)-∇U(r)·dr=-(∫_C)d U=-((Subscript[U, 2])-(Subscript[U, 1])).

Work of any potential force over a closed contour is zero. If instead of the force one has a vector function that is not a gradient of a scalar, the integral of this type over a closed contour can be non-zero. In this case it can be called circulation.

Circulation of the magnetic field

As an illustration of the circulation, consider the magnetic field produced by an infinite thin straight wire directed along the z axis and carrying a current I in the positive z direction. The magnetic field is directed in the x-y plane and in the polar coordinate system is given by

B=a((Subscript[e, ϕ])/r),      a=(((Subscript[μ, 0])I)/(2π)),

where r is the distance from the wire that is put into the origin of the coordinate system, and (Subscript[e, ϕ]) is the unit vector perpendicular to r (in the x-y plane) and directed in the direction of increasing of ϕ, in accordance with the screw rule

(Subscript[e, ϕ])=-Sin[ϕ](Subscript[e, x])+Cos[ϕ](Subscript[e, y]).

Let us first calculate the contour integral

(∫_C)B·dr

with the contour C being a circle of radius R centered at the wire and going in the positive ϕ direction. Since for this contour dr = (Subscript[e, ϕ])dr =  (Subscript[e, ϕ])Rdϕ, this integral simplifies to

(∫_C)B·dr=(∫_0^(2π))(a/R)((Subscript[e, ϕ])·(Subscript[e, ϕ]))Rdϕ=(∫_0^(2π))(a/R)Rdϕ=2πa=(Subscript[μ, 0])I,

independently of R.

In fact, the contour integral (circulation of the magnetic field) depends only on the current that crosses the surface bounded by the contour and not on the shape of the contour (Ampere's law). To illustrate this, let us use the contour that consists of an infinite horizontal line y = -(Subscript[y, 1]) going to the right and another horizontal line y = (Subscript[y, 2]) going to the left. With

(Subscript[B, x])=B·(Subscript[e, x])=-a((Sin[ϕ])/r)=-a (y/((x^2)+(y^2)))

the contour integral is given by

(∫_C)B·dr=(∫_(-∞)^∞)(Subscript[B, x])(Subscript[|, y=-(Subscript[y, 1])])dx+(∫_∞^(-∞))(Subscript[B, x])(Subscript[|, y=(Subscript[y, 2])])dx=(∫_(-∞)^∞)((a (Subscript[y, 1]))/((x^2)+((Subscript[y, 1])^2)))dx+(∫_(-∞)^∞)((a (Subscript[y, 2]))/((x^2)+((Subscript[y, 2])^2)))=(∫_(-∞)^∞)((a )/((u^2)+1))du+(∫_(-∞)^∞)((a )/((u^2)+1))du=πa+πa=2πa=(Subscript[μ, 0])I,

the same value as that along the circular contour.

Surface integrals of the scalar type

These are the integrals of the type

(∫_S)f[r]dS

where S is a surface in a 3d space. Surface can be parametrized by two parameters, r = r(u,v). The surface element dS created by changing u and v by du and dv is given by the absolute value of the vector product

dS=|((∂r)/(∂u))×((∂r)/(∂v))|dudv

For instance, the surface of a sphere can be parametrized by the two angles ϑ and ϕ, and the surface element is elementarily given by

dS=(R^2) Sin[ϑ]dϑdϕ ,

where R is the radius of the sphere. With

r={R Sin[ϑ]Cos[ϕ],R Sin[ϑ]Sin[ϕ],R Cos[ϑ]}

and

((∂r)/(∂ϑ))={R Cos[ϑ]Cos[ϕ],R Cos[ϑ]Sin[ϕ],-R Sin[ϑ]}
((∂r)/(∂ϕ))={-R Sin[ϑ]Sin[ϕ],R Sin[ϑ]Cos[ϕ],0}

one obtains

((∂r)/(∂ϑ))×((∂r)/(∂ϕ))={-(R^2) Cos[ϑ] Cos[ϕ] Sin[ϑ],-(R^2) Cos[ϑ] Sin[ϑ] Sin[ϕ],(R^2) ((Sin[ϑ])^2)}

that can be computed as

Cross[{R Sin[ϑ]Cos[ϕ],R Sin[ϑ]Sin[ϕ],R Cos[ϑ]},{-R Sin[ϑ]Sin[ϕ],R Sin[ϑ]Cos[ϕ],0}]//Simplify

The absolute value of this vector is

Simplify[
(Sqrt[%.%])
, Assumptions->{R>0, 0<ϑ<π/2}]

the same result. Let us now compute |((∂r)/(∂ϑ))×((∂r)/(∂ϕ))|in one step:

r={R Sin[ϑ]Cos[ϕ],R Sin[ϑ]Sin[ϕ],R Cos[ϑ]};
Simplify[
(Sqrt[Cross[(∂_ϑ) r,(∂_ϕ) r].Cross[(∂_ϑ) r,(∂_ϕ) r]])
, Assumptions->{R>0, 0<ϑ<π/2}]

In most cases, like the sphere, one can write down the expression for dS without using the general formula. For instant, for the paraboloid

z=((k ((x^2)+(y^2)))/2)=((k (ρ^2))/2)

ParametricPlot3D[{ρ Cos[ϕ],ρ Sin[ϕ],((k (ρ^2))/2)/.k->1},{ϕ,0,2π},{ρ,0,3}]

in the cylindrical coordinate system one has

dS=(Sqrt[1+(((((dz)/(dρ))))^2)])dρ ρdϕ=(Sqrt[1+(((k ρ))^2)])ρdρ dϕ

Excersice: derive this from the general formula. The area of the paraboloid in the region ρ < a is

S=2π Integrate[(Sqrt[1+(((k ρ))^2)])ρ,{ρ,0,a},Assumptions->{a>0,k>0}]

This result can be simplified by the Simplify command applied locally

((2 (-1+Simplify[(Sqrt[1+(a^2) (k^2)])+(a^2) (k^2) (Sqrt[1+(a^2) (k^2)])]) π)/(3 (k^2)))

whereas Simplify applied to the whole expression does not work. In the limit ak <<1 the area of the paraboloid reduces to the area of the circle ((πa)^2),as it should be. In fact, this integral is very simple because ρdρ = (1/2)d(ρ^2). Similar integral describing the length of a parabola above is much more complicated.

Surface integrals of the flux type

These are integrals of the type

(∫_S)F[r]·dS,     dS=ndS

where n is a unit vector perpendicular to the surface at the point r. If the surface is parametrized by two parameters, r = r(u,v), then dS generated by changing u and v by du and dv is given by

dS=(((∂r)/(∂u))×((∂r)/(∂v)))dudv

and the flux integral becomes

(∫_S)F[r]·(((∂r)/(∂u))×((∂r)/(∂v)))dudv,     dS=ndS

As in the case of surface integrals of the scalar type, in most cases the expression for dS can be written without calculation. For instance, for the sphere of radius R one has

dS=(r/r)(R^2) Sin[ϑ]dϑdϕ

Flux of the electric field

As an illustration we calculate the flux of the Coulomb field produced by a point charge Q

E=a(r/(r^3)),       a=(Q/(4(Subscript[πϵ, 0])))

through the sphere of radius R:

Φ=(∫_S)E[r]·dS=(∫_0^(2π))(∫_0^π)a(r/(R^3))·(r/R)(R^2) Sin[ϑ]dϑdϕ=2π a(∫_0^π)Sin[ϑ]dϑ=4π a=(Q/(Subscript[ϵ, 0])).

Gauss theorem states that the flux of the electric field depends only on the total charge inside the surface but not on the position of the charge (charges) and the surface shape. To illustrate this, let us calculate the flux of the same Coulomb field through the paraboloid

z=((k ((x^2)+(y^2)))/2)-(Subscript[z, 0])=((k (ρ^2))/2)-(Subscript[z, 0])

In the integral

Φ=(∫_S)E[r]·dS

dS is perpendicular to the paraboloid at a given point of its surface. It is convenient to make a general calculation of dS in the cylindrical coordinate system

r={ρ Cos[ϕ],ρ Sin[ϕ],z}={ρ Cos[ϕ],ρ Sin[ϕ],((k (ρ^2))/2)-(Subscript[z, 0])}.

The vector product (((∂r)/(∂ρ))×((∂r)/(∂ϕ)))can be computed as follows

r={ρ Cos[ϕ],ρ Sin[ϕ],((k (ρ^2))/2)-(Subscript[z, 0])};
Simplify[
Cross[(∂_ϕ) r,(∂_ρ) r]
]

This vector is directed out of the paraboloid. Then, the electric field on the paraboloid becomes

E=a(r/(r^3))=a(({ρ Cos[ϕ],ρ Sin[ϕ],((k (ρ^2))/2)-(Subscript[z, 0])})/((((ρ^2)+(((((k (ρ^2))/2)-(Subscript[z, 0])))^2)))^(3/2)))

Thus  E·(((∂r)/(∂ρ))×((∂r)/(∂ϕ))) is given by

a(({ρ Cos[ϕ],ρ Sin[ϕ],((k (ρ^2))/2)-(Subscript[z, 0])})/((((ρ^2)+(((((k (ρ^2))/2)-(Subscript[z, 0])))^2)))^(3/2))) .{k (ρ^2) Cos[ϕ],k (ρ^2) Sin[ϕ],-ρ}//Simplify

Now the flux through the paraboloid is given by

Φ=2π(∫_0^∞)((4 a ρ (k (ρ^2)+2 (Subscript[z, 0])))/(((4 (ρ^2)+(k^2) (ρ^4)-4 k (ρ^2) (Subscript[z, 0])+4 (Subscript[z, 0]^2)))^(3/2)))dρ=2π(∫_0^∞)((2 a  (k u+2 (Subscript[z, 0])))/(((4 u+(k^2) (u^2)-4 k u (Subscript[z, 0])+4 (Subscript[z, 0]^2)))^(3/2)))du

that is

Φ= Integrate[((2π 2 a  (k u+2 (Subscript[z, 0])))/(((4 u+(k^2) (u^2)-4 k u (Subscript[z, 0])+4 (Subscript[z, 0]^2)))^(3/2))),{u,0,∞},Assumptions->{(Subscript[z, 0])>0,k>0}]

This result is correct and in accordance with the Gauss theorem saying that the flux must be the same as the flux through the sphere. In the case (Subscript[z, 0])<0 the charge is outside the paraboloid and the flux is zero:

Φ= Integrate[((2π 2 a  (k u+2 (Subscript[z, 0])))/(((4 u+(k^2) (u^2)-4 k u (Subscript[z, 0])+4 (Subscript[z, 0]^2)))^(3/2))),{u,0,∞},Assumptions->{(Subscript[z, 0])<0,k>0}]

In the case (Subscript[z, 0])<0 the charge is exactly on the surface of the paraboloid and the flux is

Φ= Integrate[((2π 2 a  (k u+2 (Subscript[z, 0])))/(((4 u+(k^2) (u^2)-4 k u (Subscript[z, 0])+4 (Subscript[z, 0]^2)))^(3/2))),{u,0,∞},Assumptions->{(Subscript[z, 0])==0,k>0}]

The same result arises if the charge is placed on the surface of a sphere.

Vector integral relations: Stokes  and Gauss theorems

There are important vector integral relations between integrals over a region and over the border of this region. The simplest example of these relations is the formula

(∫_a^b)(∂_x) F[x]dx=F[b]-F[a]

and a slightly more complicated formula for the integral of a gradient over a contour. In both cases, the region is one-dimensional and its border is zero-dimensional, thus there is no integral over the border, just values of the functions at the border points. Next of these relations is the Stokes' theorem that relates the flux through an open surface (a 2d integral) with a work-type integral over a contour bounding this surface. It has the form

(∫_S)▽×F·dS=(∫_C)F·dr .

There is a one-dimensional integral on the right but a two-dimensional integral on the left, however, of a differentiated function, the curl. The physical implication of the Stokes theorem is the Ampere's law, where F is the magnetic field B and ▽×F is proportional to the density of current flowing through the surface.

Next there is the Gauss theorem relating a volume integral with a flux through the closed surface bounding this volume

(∫_V)▽·FdV=(∫_S)F·dS .

Here there is a two-dimensional integral on the right but a three-dimensional integral on the left, however, of a differentiated function, the divergence. The physical implication of the Gauss theorem is the relation between the electric-field flux out of the volume V through its surface S, F being the electric field E, and the total electric charge in the volume, the charge density being proportional to ▽·E.

There is a vector variant of the Gauss theorem

(∫_V)▽×FdV=-(∫_S)F×dS

arising in the theory of magnetism.