I am working on solving the Biot-Savart Law equation for the magnetic field around a charged ring of uniform current density. The expression that Mathematica gives is rather nasty, as expected. Notice that the limits of integration are not included in the solution. Mathematica rejects the bounds of integration by silently refusing to even consider a solution (literally, it sits there for several minutes "Running..." and then just stops like nothing ever happened). My idea was that since Mathematica didn't mind handing over the indefinite integral, that I could turn it into a function of theta and manually apply the first fundamental theorem of calculus. For some reason, I can't seem to convince Mathematica that I want this to be a function that I want evaluated at 2Pi and at 0 so that they can be subtracted. Here is the code up to the integration. Thanks in advance for any feedback.

l[Theta_]:={s Cos[Theta],s Sin[Theta],0};
B=(m/4*Pi) Integrate[(Integrand),Theta,GenerateConditions ->False]
  • 1
    $\begingroup$ It looks like a ghastly integral. I am not so sure the results, symbolically, would be of any use. Is there a reason not to make this into a numerical problem wherein, for given (x,y,z), you use NIntegrate to get a value? $\endgroup$ – Daniel Lichtblau Sep 21 '14 at 0:33
  • 1
    $\begingroup$ I don't recall off the top of my head whether there exists a closed-form solution for the magnetic field of a current loop, although there are closed-form solutions for the on-axis field strength, and they are typically derived in intro E&M classes. For the off-axis case, if a closed-form solution does not exist, there exist multipole-expansions which are available in closed-form, but again, they are series expansions. $\endgroup$ – DumpsterDoofus Sep 21 '14 at 2:19
  • 1
    $\begingroup$ See my answer at physics.stackexchange.com/questions/107632/…, I don't know if this helps, but it's at least a usable series approximation. $\endgroup$ – DumpsterDoofus Sep 21 '14 at 2:24

The integrals can in fact be done exactly, but only if you make some use of the symmetries of the problem first.

  • The circular ring geometry implies that the magnetic field will look the same in any vertical plane going through the rotation axis (which we call the z axis). Therefore, we don't need to specify three independent variables x, y and z to do the calculation. Instead, set y = 0 and look only at the xz plane.
  • The component of the magnetic field perpendicular to the xz plane must be zero, so we don't have to calculate the integral for it.
  • Use spherical polar coordinates to specify the point of interest, instead of Cartesian coordinates
  • In the non-zero integrals over the ring, the integration variable Theta appears only in the form Cos[Theta] which is symmetric around Theta = π. Therefore, all integrals can be reduced to the domain {Theta, 0, π} if we multiply the result by 2 at the end.

With these points in mind, the calculation goes like this (using exactly your setup): replace the Integrand by int2 which uses spherical coordinates with the azimuthal angle set to zero (xz plane), then observe that we only need to integrate int[[1]] and int[[3]] to get Bx and Bz:

s = 1;
m = 4*Pi*10^7;
l[Theta_] := {s Cos[Theta], s Sin[Theta], 0};
R = {x, y, z} - l[Theta];
r = R/Sqrt[R.R];
Deel = D[l[Theta], Theta];
Cross[Deel, r];
Integrand = %/R.R;
int2 = Simplify[
  Integrand /. 
   Thread[{x, y, 
       z} -> ρ {Cos[ϕ] Sin[θ], 
        Sin[ϕ] Sin[θ], Cos[θ]} /. ϕ -> 
      0], ρ > 0 && θ > 0 && ϕ > 0 && Theta > 0]

$$\left\{\frac{\rho \cos (\theta ) \cos (\text{Theta})}{\left(\rho ^2-2 \rho \sin (\theta ) \cos (\text{Theta})+1\right)^{3/2}},\\ \frac{\rho \cos (\theta ) \sin (\text{Theta})}{\left(\rho ^2-2 \rho \sin (\theta ) \cos (\text{Theta})+1\right)^{3/2}},\\ \frac{1-\rho \sin (\theta ) \cos (\text{Theta})}{\left(\rho ^2-2 \rho \sin (\theta ) \cos (\text{Theta})+1\right)^{3/2}}\right\}$$

Bx = 
 2 Assuming[ρ > 0 && Pi/2 > θ > 0, 
   Integrate[int2[[1]], {Theta, 0, Pi}]]

$$\frac{2 \cot (\theta ) \left(\left(\rho ^2+1\right) E\left(\frac{4 \rho \sin (\theta )}{\rho ^2+2 \sin (\theta ) \rho +1}\right)-\left(-2 \rho \sin (\theta )+\rho ^2+1\right) K\left(\frac{4 \rho \sin (\theta )}{\rho ^2+2 \sin (\theta ) \rho +1}\right)\right)}{\left(-2 \rho \sin (\theta )+\rho ^2+1\right) \sqrt{2 \rho \sin (\theta )+\rho ^2+1}}$$

Bz = 
 2 Assuming[ρ > 0 && Pi/2 > θ > 0, 
   Integrate[int2[[3]], {Theta, 0, Pi}]]

$$\frac{2 \left(\left(-2 \rho \sin (\theta )+\rho ^2+1\right) K\left(\frac{4 \rho \sin (\theta )}{\rho ^2+2 \sin (\theta ) \rho +1}\right)-\left(\rho ^2-1\right) E\left(\frac{4 \rho \sin (\theta )}{\rho ^2+2 \sin (\theta ) \rho +1}\right)\right)}{\left(-2 \rho \sin (\theta )+\rho ^2+1\right) \sqrt{2 \rho \sin (\theta )+\rho ^2+1}}$$

Here, $E$ and $K$ are the elliptic integrals EllipticE and EllipticK.

   Bz} /. {θ -> ArcTan[z, x], ρ -> Sqrt[
    x^2 + z^2]}, {x, -2, 2}, {z, -2, 2}]


This is a plot of the exact result.


(I know this is an old question, but it just got nudged and I saw another way to do it.)

You can also take advantage of the vector potential to do this: $$ \vec{A}(\vec{r}) = \frac{\mu_0 I}{4\pi} \oint \frac{d\vec{r}'}{|\vec{r} - \vec{r}'|} $$ Using much the same symmetry arguments as given in the answer above, we can first look at the vector potential in the $xy$-plane. Its components will be $$ A_x = - \frac{\mu_0 I}{4\pi} \int_0^{2\pi} \frac{\sin \phi}{\sqrt{R^2 + s^2 -2 R s \cos \phi + z^2}} d\phi $$ $$ A_y = \frac{\mu_0 I}{4\pi} \int_0^{2\pi} \frac{\cos \phi}{\sqrt{R^2 + s^2 -2 R s \cos \phi + z^2}} d\phi, $$ where $R$ is the radius of the ring, and $s, z$ are the locations of the field point in cylindrical coordinates. (We have used symmetry to set the field point at $\theta = 0$. Mathematica knows how to do both of these integrals, though it does take a few moments for it to solve them:

{Ax, Ay} = Integrate[{-Sin[ph], Cos[ph]}/Sqrt[s^2 + R^2 - 2 R s Cos[ph] + z^2], 
  {ph, 0, 2 Pi}, Assumptions -> {s > 0, R > 0, R != s, Element[z, Reals]}]

enter image description here

To find the magnetic field, we then take $\vec{B} = \nabla \times \vec{A}$. Working in cylindrical coordinates (i.e., substituting $A_x \to A_s$ and $A_y \to A_\phi$,) we have $$ \vec{B} = - \frac{\partial A_\phi}{\partial z} \hat{s} + \frac{1}{s} \frac{\partial (s A_\phi)}{\partial s} $$ Mathematica can take these derivatives as well:

B = Simplify[{-D[Ay, z], 1/s D[s Ay, s]}];
StreamPlot[B /. R -> 1, {s, -2, 2}, {z, -2, 2}, AspectRatio -> Automatic]

enter image description here


You might also be interesting in looking at web-Mathematica (version 5) solution below.


  • $\begingroup$ FYI the link is broken. $\endgroup$ – chris Mar 25 '20 at 16:07
  • 1
    $\begingroup$ @chris Corrected, it is however, old interactive project written on webMathematica 5, and after system upgrade not fully working. $\endgroup$ – Acus Mar 25 '20 at 16:21

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.