Dipole of a Linear Charge Distribution

PROBLEM:

Consider the case of a linear charge distribution with the following charge density:

–a < x < a                   λ[x] = λo Sin[((π x)/a)]
x < –a    or  x > a        λ[x] = 0

Find the electric dipole moment for this charge distribution.

Note:  Adding a time-dependent factor yields a charge distribution found in dipole antennas for radio transmission.  Dipole antennas are perhaps the simplest possible antennas.  The dipole pattern of radiation is used to model more complex radiation phenomena.

(a)  Plot this charge density as a function of x for –2 a < x < 2 a.  You will have to choose numerical values for λo and a.  Comment on the result.

ClearAll["Global`*"]
λo=2; a = 2;
lambda[x_] = If[x >-a && x < a, λo Sin[((π x)/a)], 0]
Plot[lambda[x],{x,-2 a,2 a}, AxesLabel->{x, λ[x]}]

The plot conforms to the definition of λ[x]:  λ[x] = 0 when x < –a = –2 and x > a = 2.  When –2 = a < x < a = 2, λ[x] is a sign function with amplitude λo = 2 and period  2 a = 4.

(b)  Derive an expression (an integral) for the electric dipole of this charge distribution.  You may start by defining the dipole as an integral over the 3-D charge density ρ(x,y,z) or as the sum of the dipoles due to a large number of finite point charges.  Below, a prime is used to mark source variables.  For instance, x’ is the x-coordinate of the position of an element of charge.

(OverVector[p])= ∫(OverVector[r])' ρ[x',y',z']dτ'   or   (OverVector[p])=(∑_(i=1)^N)(OverVector[(Subscript[r, i])])'  (Subscript[Δq, i])

<Insert text here>
We view a linear charge distribution as the sum over very small charge elements with charges Δq_i at positions x’_i on the x-axis.  Each charge element contributes a vector-valued dipole to the total, equal to

Δ(Subscript[(OverVector[p]), i]) = (OverVector[r])' (Subscript[Δq, i])

where r’ is the position vector from the origin to the charge element.  (r’ is the source variable—the position of the charge element).  If our charge elements lie on the x-axis at positions (Subscript[x, i]^′), then

(OverVector[r])' = (Subscript[x, i]^′)  (OverHat[x])

The electric dipole has only an x-component.  In terms of λ[x’], the charge element Δq may be expressed as

(Subscript[Δq, i]) = λ[(Subscript[x, i]^('))] (Subscript[x, i]^('))  so that

			Δ(Subscript[(OverVector[p]), i]) = (Subscript[x, i]^(')) (OverHat[x]) λ[(Subscript[x, i]^('))] (Subscript[Δx, i]^(')) =λ[(Subscript[x, i]^('))] (OverHat[x]) (Subscript[x, i]^(')) (Subscript[Δx, i]^('))

The total dipole is computed by taking the limit as Δx → 0 (while the number of charge elements goes to infinity.  This is just the recipe for an integral over x.  One can avoid the If statement in λ[x] above by restricting the range of integration (normally -∞ to +∞) to the region(s) where λ[x] is nonzero.

(OverVector[p]) = λo (OverHat[x]) (Subscript[x, i]) (Subscript[Δx, i])  = (OverHat[x])(∫_(-∞)^∞)x' λ[x'] dx' = (OverHat[x]) (∫_(-a)^a)x' λo Sin[((π x')/a)] dx'

Note that Griffiths defines the dipole in terms of a 3-D charge density.

(OverVector[p]) = ∫r' ρ[x',y'z']dτ'=∫{x',y',z'}ρ[x',y'z'] dx' dy' dz'

where the integral is taken over the region of nonzero charge density, and {x’,y’,z’} is the position vector for the source coordinates.  

The 3-D definition is theoretically “more general”, but care is required in applying it to 1-D and 2-D charge distributions.  To use the 3-D definition in our 1-D case, one can write the three dimensional charge density ρ[x,y,z] in terms of delta functions.  For instance,

ρ[x',y',z'] = λo Sin[((π x')/a)]δ[y'] δ[z']  for -a < x < a.

Recall that the Dirac Delta function is defined so that

(∫_a^b)f[y'] δ[y'-c]dy'=({{{, ({{f[c], if a<c<b}, {0, if c < a  or c>b}})}})

so that

(∫_a^b)δ[y']dy'=({{{, ({{1, if a<0<b}, {0, if 0 < a  or  0 > b}})}})

The Dirac Delta functions in ρ[x’, y’, z’] ensure that the integral over ρ is zero unless y’ = z’ = 0; that is, everywhere except along the x-axis.  They also ensure that when you integrate over dy’ and dz’, you get a charge density of λ[x’] along the x’-axis.  

This new charge density be substituted directly into the integral expression for the dipole.

(OverVector[p])= ∫(OverVector[r])' ρ[x',y',z']dτ' = (∫_(x'=-a)^(x'=a))(x' (OverHat[x]) + y' (OverHat[y]) + z' (OverHat[z])) λo Sin[x']δ[y'] δ[z']dx ' dy' dz'

(c)  Compute the electric dipole of this charge distribution.  Since M interprets the prime in x’ as a derivative, you should either remove the prime or substitute another variable.

Below I have used the 1-D integral from (b) to calculate the dipole.  Print[] statements are used to add labels to the output.  This makes the results easier to read and interpret.

ClearAll["Global`*"]
$Assumptions-> λo ∈ Reals && a ∈ Reals; 
lambda[x_] = λo Sin[((π x)/a)];
Print["λ[x] = ", lambda[x]]
Print["p = ", Integrate[x lambda[x], {x,-a, a}]]

So the dipole of the λoSin[π x/a] charge distribution is  (((2 (a^2) λo)/π)).  The direction of the dipole is  given by the unit vector pointing from the negative charge to the positive charge, that is, in the positive x-direction.  In M, this is the vector (OverVector[p]) = {((2 (a^2) λo)/π), 0, 0}.  Obviously, one could also write it as  (OverVector[p]) =  ((2 (a^2) λo)/π)(OverHat[x]).

To compute the dipole using the 3-D charge density and delta functions, you can use M’s delta function command.  For instance, δ[y’] = DiracDelta[y’].   With this substitution, M can evaluate the dipole directly using the 3-D charge density from Griffiths.  Griffiths writes the position vector as {x’,y’,z’}, but the primes must be omitted to avoid confusing M.

ClearAll["Global`*"]
$Assumptions-> λo ∈ Reals && a ∈ Reals; 
rho[x_,y_,z_] = λo Sin[((π x)/a)] DiracDelta[y] DiracDelta[z];
dipole = Integrate[{x,y,z} rho[x,y,z], {x,-a,a}, {y,-∞,∞}, {z,-∞,∞}];
Print["(OverVector[p]) = ", dipole]

The result is a vector, as expected, with only one nonzero component (the x-component) and agrees with what we wrote above.   Also expected, the magnitude of this x-component agrees with the dipole calculation above.

(d)  To get a feel for how the dipole works, it may help to construct a point-charge model of this charge distribution, where all the positive charge (all the charge for x > 0) is located at one point and all the negative charge (all the charge for x < 0) is located at another point.  

Begin by calculating the amount of positive charge in the region where x is positive?  (By symmetry, this will equal the amount of negative charge in the regions where x is negative.)

ClearAll["Global`*"]
$Assumptions-> λo ∈ Reals && a ∈ Reals; 
lambda[x_] = λo Sin[((π x)/a)];
Qrighthand = Integrate[lambda[x], {x,0,a}];
Print["Total positive charge for x > 0 = ",Qrighthand]

The total charge for 0 > x > a is (2a/π)λo ~ (0.64 a λo).   To check the reasonableness of this result, one can estimate the area under the triangle formed by the two ends and the peak of the charge distribution; that is, the points (0,0), (a/2,λo), and (a,0).  This area is (0.5 a λo).  The true total charge will be a little higher, because the λoSin[πx/a] function equals or exceeds the value assumed by the triangle at every point.

(e)  We could approximate the potential due to the above charge distribution by assuming that a charge –(2a/π)λo is located at x = –d and a charge (2a/π)λo is located at x = +d.  (This is a physical dipole.)  What would d have to be to produce the same total dipole as our λoSin[πx/a] charge distribution?  Comment on your result.

<Insert text here>
The dipole due to a pair of point charges with charge ±q at points ±d is p = 2qd.  (I assume the positive charge is located at x = d, the negative charge is located at point x = –d, and q is the magnitude of the two point charges.)  For the Sin[πx/a] charge distribution, we have shown that  p = 2 a^2 λo/π.  We want 2qd = (2a/π)a λo, where q = (2a/π)λo.   Therefore:

d = (p/(2q))= ((2 (a^2) λo / π)/(4 a λo /π)) = (a/2)

The charge must be located at x = ±a/2.  Interestingly, the λoSin[πx/a] function has it’s negative peak at x = –a/2 and it’s  positive peak at x = a/2.  In terms of the resulting dipole, it is as if the entire charge in the positive side of the charge distribution (2aλo/π) were concentrated at x = a/2, and the entire charge in the negative side of the charge distribution (–2aλo/π) were concentrated at x = –a/2.

Note that we do not expect the electric field of the λoSin[πx/a] charge distribution to be identical to the field of an electrical dipole, ideal or not.  

With this information, we have what we need to calculate the potential of the λoSin[πx/a] at three levels of approximation:

(i)	The “exact” potential found by integrating the λoSin[πx/a] charge distribution.
(ii)	A better approximation, assuming a physical dipole with the entire positive charge located at 
	x = a/2 and the entire negative charge located at x = –a/2.
(iii)  An ideal dipole located at the origin.

(g) Compare contour plots of the “exact” potential due to the λoSin[π x/a] potential found by direct integration with contour plots of the corresponding physical and ideal dipoles.  You must choose suitable values for λo and a.  I chose λo = 1 C/m and a = 1 m.  

When computing the “exact” potential, you will probably find it necessary to perform the integration numerically.  This approach often fails near the origin, where r’ ~ 0 and the value of the integral diverges.  This problem can be sidestepped by performing the numerical integration (NIntegrate[]) inside the ContourPlot[] function.  The ContourPlot function is designed to ignore “off-scale” function values.  It will (usually) ignore the problematic results. 

Try it.

<Input text here>  
Far away from the dipole, we expect the potential distributions to be similar.  This is the principle motivation for the multipole expansion.  To observe differences, we should employ distance scales comparable to the “size” of the charge distribution, which extends from –a to a along the x-axis.  Below I plot the potential in a 2 by 2 square centered on the charge distribution.  The “exact” potential distribution is shown first.  

Note that the potential exists in all three dimensions.  For plotting purposes, I show the potential in the x-y plane only.  The z-dependence can be determined by rotating the plot around the x-axis or  adding a z-term to the expression for r, so that r = Sqrt[x^2 + y^2 + z^2] .  I have replaced x’ by ξ below to avoid confusing M with primes.  For this charge distribution, all the charge is located at points for which y’ = z’ = 0.   

The differential for computing the “exact” potential is:

dV = (dq/((r-r'))) = ((λ[x] dx)/(Sqrt[(((x-ξ))^2) + (y^2)]))

I numerically integrate this over the non-zero portion of the charge distribution with M.  It is not a trivial integral.  It may take a few minutes, depending on your computer.

ClearAll["Global`*"]
λo = 1;  a = 1;
ContourPlot[ NIntegrate[((λo Sin[((π ξ)/a)])/(Sqrt[(((x-ξ))^2) + (y^2)])), {ξ,-a, a}], {x,-2,2}, {y,-2,2}, PlotLabel->"Potential Near Origin for Sin[πx/a] Charge Distribution in x-y Plane"]

The second level of approximation, that of the physical dipole, for the x-y plane only, is given by:

V2q[x,y] = ((+ 2 (a^2) λo / π)/(Sqrt[(((x - a/2))^2)+ (y^2)])) + ((- 2 (a^2) λo / π)/(Sqrt[(((x - a/2))^2)+ (y^2)]))

with its contour potential plot:

ContourPlot[((+ 2 (a^2) λo)/(π Sqrt[(((x - a/2))^2)+ (y^2)])) + ((- 2 (a^2) λo )/(π Sqrt[(((x + a/2))^2)+ (y^2)])), {x,-2,2}, {y,-2,2}, PlotLabel->"Potential Near Origin for Equivalent Physical Dipole in x-y Plane"]

In the x-y plane, the potential due to an ideal dipole located at the origin is given by

V[r] = ((p cosθ)/(r^2))=((p x)/(r^3))=((p x)/((((x^2)+(y^2)))^3))

The corresponding contour plot is:

ContourPlot[((2 (a^2)λo x/π)/((((x^2) + (y^2)))^(3/2))), {x,-2,2}, {y,-2,2}, PlotLabel->"Potential Near Origin for Equivalent Ideal Dipole in x-y Plane"]

The most significant differences between the three plots are near the origin.  Although the potential along the y-axis (x = 0) is all three cases, the potential at nearby points (x ≠ 0, both x, y ≪ a) increases in the following order:

(ii)	Physical dipole.  The nearest charge is relatively far away, at x = ± a/2.
(i)	Sine distribution.  Some charge is near origin, but not much.
(iii)	Ideal dipole:  All the charge of the dipole is located very close to the origin.

Higher potentials are found at points closer to points with non-zero charge.  As you move away from the origin along the y-direction, with x ≠ 0 but x ≪ a/2, the magnitudes of potentials of the three charge distributions rank the same way.  

However, as you move far away from the origin along the ±x axes, (x > a/2 and y ≪ a/2), the potentials tend to be lowest for the ideal dipole and highest for the Sin[x] distribution.  The ranking is (from lowest to highest):  

(iii)  Ideal dipole.  Charge concentrated at the origin. 
(ii)	Physical dipole.  Charge concentrated at x = ± a/2.
(i)	Sine distribution.  Some charge is found all the way out to x = ±a.

In general, one might expect that the best approximation to the “true” potential near the origin is found using approximations where the charge density at nearby points is most similar to the true charge distribution.  This depends on exactly where you look.  Any approximations you make in this case should be checked carefully.  

At large distances from the origin (at large distances from the charge), the true potentials and the two approximations should be quite similar.  The same three potentials are compared at greater distances from the origin in the 20 by 20 square plots below.

ClearAll["Global`*"]
λo = 1;  a = 1;
ContourPlot[NIntegrate[((λo Sin[((π ξ)/a)])/(Sqrt[(((x-ξ))^2) + (y^2)])) , {ξ,-a, a}], {x,-20,20}, {y,-20,20}, PlotLabel->"Potential Far from Origin for Sin[πx/a] Charge Distribution"]

ContourPlot[((+ 2 (a^2) λo)/(π Sqrt[(((x - a/2))^2)+ (y^2)])) + ((- 2 (a^2) λo )/(π Sqrt[(((x + a/2))^2)+ (y^2)])), {x,-20,20}, {y,-20,20}, PlotLabel->"Potential Far from Origin for Equivalent Physical Dipole"]

ContourPlot[((2 (a^2)λo x/π)/((((x^2) + (y^2)))^(3/2))), {x,-20,20}, {y,-20,20},PlotLabel->"Potential Far from Origin for Equivalent Ideal Dipole in x-y Plane"]

Outside the high-potential white area of the 20 x 20 contour plots, the three sets of potentials are nearly identical.  You don’t have to get far from the charged wire before its potential is well approximated by the potential of an ideal dipole.  Multipole expansions can be quite precise for points well outside the region occupied by charge, and are significantly easier (faster) to calculate.

Again, these plots show the potential in the x-y plane only.  One can use symmetry to determine the potential as a function of all three variables {x, y, z}.

For the hardcore:

Can there be higher order moments ?    Let’s test for an octopole moment
In Jackson’s notation, the m = 0 term (no ϕ-dependence) of the octopole moment  is given by:

qlm = ∫Ylm(θ,ϕ)(((r'))^3)ρ((OverVector[x])')(d^3)x'
where
Y30(θ,ϕ) = (Sqrt[(7/(4 π))])((5/2)(cos^5)θ -(3/2)(cos^3)θ)

The Sin[πx/a] distribution is confined to the case of θ = 0 and θ = π, where 
	r’ = x’ for  θ = 0
	r’ = –x’ for θ = π.

With some cancelling minus signs, this yields
qlm = (Sqrt[(7/(4 π))])∫((5/2) -(3/2)(cos^3)θ)(((x'))^3)λ((OverVector[x])')dx'  for θ'=1
	= (Sqrt[(7/(4 π))])∫(-(5/2) +(3/2))(((x'))^3)λ((OverVector[x])')dx'  for θ'= π
q30 = (∫_(-a)^a)(((x'))^3)λ((OverVector[x])')dx'

For λ[x’] = λo Sin[π x / a]

ClearAll["Global`*"]
$Assumptions-> λo ∈ Reals && a ∈ Reals; 
lambda[x_] = λo Sin[((π x)/a)];
Print["λ[x] = ", lambda[x]]
Print["q30 (magnitude) = ", Sqrt[(7/(4 π))] Integrate[(x^3) lambda[x], {x,-a, a}]]

So q30 ≠ 0. 

 I expect some non-zero, higher order terms.  These are needed to make the expansion converge on the “true” potential.