Chapter 3: Emission from accelerated charges

Below is a set of python codes associated with Chapter 3 of Daniele Pelliccia and David M. Paganin, “Synchrotron Light: A Physics Journey from Laboratory to Cosmos” (Oxford University Press, 2025).

In order to run any of these python codes, you will need to include the following header file.

Plotting radial and transverse fields

To start off, let’s work out a way to plot the electric field generated by a charge +q located at the position \mathbf{r}_{0} = (x_{0}, y_{0}). This is given by the familiar Coulomb law. The field falls off as 1/r^{2} and is directed radially outwards from the charge. We can write this in vector form:

\mathbf{E}_{r} = \frac{1}{4 \pi \epsilon_{0}} \frac{q}{\vert \mathbf{r}-\mathbf{r}_{0} \vert^{2}}\mathbf{\hat{r}}.

To draw a vector plot, we need to find an expression for the field in Cartesian coordinates. The radial unit vector can be written as a superposition of the Cartesian unit vectors as follows:

\mathbf{\hat{r}} = (\cos \Theta) \mathbf{\hat{x}} + (\sin \Theta) \mathbf{\hat{y}}.

By inverting the usual transformations in polar coordinates:

x = r\cos \Theta,y = r \sin \Theta,

we get

\cos \Theta = x/r,\sin \Theta = y/r.

Finally, by substituting in the first expression, we get the field components in Cartesian coordinates:

\mathbf{E}_{r} = \frac{1}{4 \pi \epsilon_{0}} \frac{q}{\vert \mathbf{r}-\mathbf{r}_{0} \vert^{3}}\left ( (x-x_{0})\mathbf{\hat{x}} + (y-y_{0}) \mathbf{\hat{y}} \right).

One last thing to remember is that r = \sqrt{x^{2} + y^{2}}. We can now use Numpy and Matplotlib to define and plot these functions (see Fig. 3.4).

The parameters we need to pass to the streamplot function are the X and Y position, and the direction of the field E_{x} and E_{y}. The length of each arrow is proportional to the magnitude of the field at that position.

The amplitude of the field is encoded in the colour of the field lines. The large variation in magnitude of the electric field is not easily visualised using a linear scale. For this reason we define the colour scale to be the logarithm of the amplitude.

If we are only interested in a qualitative plot of the field lines, we can draw a streamplot without axes and colour bar.

The expression we derived for the transverse field is (we suppose for simplicity that the charge is centred at the origin)

\mathbf{E}_{t} = \frac{1}{4 \pi \epsilon_{0}} \frac{qa \sin\theta}{c^{2}r}\boldsymbol{\hat{\theta}}.

Note that the direction of the transverse field, \boldsymbol{\hat{\theta}}, is by construction perpendicular to the direction of the radial field. The corresponding transformation for the angular unit vector into Cartesian coordinates is

\boldsymbol{\hat{\theta}} = (\sin \theta) \mathbf{\hat{x}} - (\cos \theta) \mathbf{\hat{y}}.

By using the transformation for polar coordinates, the transverse field in Cartesian coordinates takes the form:

\mathbf{E}_{t} = \frac{1}{4 \pi \epsilon_{0}} \frac{qa y}{c^{2}r^{3}} \left( y\mathbf{\hat{x}} - x\mathbf{\hat{y}} \right).

Following exactly the same steps seen before, here we draw the streamplot of the transverse field. Note, as expected, the transverse field is zero when \sin \theta = 0, i.e. in the direction of the charge acceleration. The field shows the typical ‘double-lobe’ structure.

Spatial distribution of emitted power in a circular accelerator

As a warm-up, we plot the 2D contours of the emitted power in the frame of reference of the moving charge, using the proportionality relation

\frac{\mathrm{d}\mathcal{P}}{\mathrm{d}\Omega} \propto \frac{q^{2} a^{2}}{c^{3}} \sin^{2} \Theta.

If we neglect the various multiplicative constants, this contour will correspond to an equation, in polar coordinates, of the type: R = \sin^{2} \Theta. This is the trick to produce this kind of contour plot in Fig. 3.6, which we’ll extend to 3D later on.

In essence, we are going to define a grid in polar coordinates and then transform to Cartesian coordinates with the constraint R = \sin^{2} \Theta.

3D surface of constant power

Let’s now plot the 3D surface of constant received power (cf. Fig 3.13). The approach is exactly the same as we saw before.

Firstly we define spherical coordinates \theta and \phi. The we calculate the surface of constant power as in Eq. (3.66), which upon setting all the constants equal to one, reads:

\frac{\mathrm{d}\mathcal{P}}{\mathrm{d}\Omega} = \frac{1}{\left( 1 - \beta \cos \theta \right)^{4}} \left[ 1 - \frac{\sin^{2} \theta \cos^{2} \phi}{\gamma^{2} \left( 1 - \beta \cos \theta \right)^{2}} \right].

You can use \beta as a parameter and calculate the 3D contour as follows.

We can chart several curves at once (in 2D) by restricting ourselves to the plane \phi=0 and setting the calculation of the power distribution inside a function (see Fig. 3.15).

Synchrotron beam divergence in the ultra-relativistic case

The power per unit solid angle, in the limit \gamma \gg 1 , and setting all the multiplicative constants equal to one is

\frac{\mathrm{d}\mathcal{P}}{\mathrm{d}\Omega} = \frac{1}{(\theta^{2} + \gamma^{-2})^{4} } \left[ 1 - \frac{4\theta^{2}\gamma^{-2}}{(\theta^{2} + \gamma^{-2})^{2}} \right].

Here’s how to plot the graph in Fig. 3.14.

Spatial distribution of emitted power in a linear accelerator

In this case the relevant expression is eqn (3.87). As usual we set all the multiplicative constants equal to one to get

\frac{\mathrm{d}\mathcal{P}}{\mathrm{d}\Omega} = \frac{\sin^{2}\theta}{\left( 1 - \beta \cos \theta \right)^{6}}.

The above image coincides with Fig. 3.16(a).  Below is the corresponding plot for the 2D case, as shown in Fig. 3.17.