A pedagogical review of the vacuum retarded dipole model of pulsar spin down

Abstract

Pulsars are rapidly spinning highly magnetised neutron stars. Their spin period is observed to decrease with time. An early analytical model for this process was the vacuum retarded dipole (VRD) by Deutsch (1955, AnAp, 18). This model assumes an idealised star and it finds that the rotational energy is radiated away by the electromagnetic fields. This model has been superseded by more realistic numerical simulations that account for the non-vacuum like surroundings of the neutron star. However, the VRD still provides a reasonable approximation and is a useful limiting case that can provide some qualitative understanding. We provide detailed derivations of the spin down and related electromagnetic field equations of the VRD solution. We also correct typographical errors in the general field equations and boundary conditions used by Deutsch (1955, AnAp, 18).

1. Introduction

Pulsars are highly magnetised rapidly rotating neutron stars that emit electromagnetic radiation across almost the entire spectrum, from radio to $\gamma$ -rays (Travelle 2011). Because of the extreme magnetic and electric fields, the environment around a pulsar proves to be interesting and is believed to be the main cause of a pulsar’s energy emission. This energy emission results in a decreasing spin-rate for the neutron star. An analytical approximation for this process is given by the Vacuum Retarded Dipole (VRD) model which was originally derived by Deutsch (1955).

The derivations of the results in Deutsch (1955) are very concise, with many details missing. The purpose of this review article is to rederive the main results of Deutsch (1955) and in doing so provide details on how this derivation takes place mathematically and physically. We also uncover three typographical errors in Deutsch (1955), the most important of which is in our Equation (31), the polar component of the electric field. It is worth noting the first papers we found that mentioned the error in the general field equations are Melatos (1997) and Michel & Li (1999). However, they do not show the derivations for the general field equations and how the typographical error is found. A reasonably detailed derivation of the VRD solution, without the typos, is shown in the textbook by Michel (1991). However, the method used by Michel is different from Deutsch (1955). While Deutsch (1955) uses methods outlined by Stratton (1941), Michel just considers multipole expansions at an inclination and does not derive the normalisations. Also the following working carried out by us goes into greater detail for the derivations. In this review we use methods outlined in Stratton (1941). In particular, we refer frequently to Chapters VII and VIII, which consider spherical waves and radiation respectively. We use Griffiths (2017) as a reference for standard electrodynamics theory.

In Section 2, we begin by outlying the definition of an idealised star, along with fundamental equations that apply due to this assumption, which include Maxwell’s equations. We also note down the Alfvén theorem and a brief description of a vector field. In Section 3, we define the magnetic field symmetry axis and its associated coordinate system. The magnetic field is assumed to take the form of a dipole about this axis. The transformation between the two coordinate systems is shown and used to derive the electromagnetic field boundary conditions in the rotation symmetry axis coordinates. In Section 4, we show the derivations for the time-dependent and static components of the electromagnetic fields using series expansions and matching coefficients to the boundary conditions. In Section 5, we simplify the form of the electromagnetic fields for a large radius, then the power and spin down equations are derived.

2. The internal fields of an idealised star

For the VRD solution, we consider the electric and magnetic fields of an idealised star. An idealised star is a star that is a sharply bounded sphere, is isotropic, is perfectly conducting, and has permittivity and permeability of a vacuum everywhere inside. The boundary is not continuous but rather an immediate change between the star’s medium and the surrounding vacuum. Which, in the case of a real star is not true. However, it greatly simplifies the situation.

2.1. Maxwell’s equation

The Maxwell equations in SI units are,

(1) \begin{gather} \nabla\cdot\boldsymbol{{E}}=\frac{\rho}{\epsilon_0}, \end{gather}

(2) \begin{align} \nonumber\\[-20pt] \nabla\cdot\boldsymbol{{H}}=0, \end{align}

(3) \begin{align} \nabla\times\boldsymbol{{E}}=-\mu_0\frac{\partial\boldsymbol{{H}}}{\partial t}, \end{align}

(4) \begin{align} \nabla\times\boldsymbol{{H}}=\boldsymbol{{J}}+\epsilon_0\frac{\partial\boldsymbol{{E}}}{\partial t}, \end{align}

where $\boldsymbol{{E}}$ is the electric field, $\boldsymbol{{H}}$ is the magnetic field (related to $\boldsymbol{{B}}$ by $\boldsymbol{{B}}=\mu_0\boldsymbol{{H}}$ ), $\epsilon_0$ is the permittivity of free space, and $\mu_0$ the permeability of free space. $\boldsymbol{{H}}$ will typically be referred to as the H field. In Deutsch (1955), there is a missing minus sign in his equivalent of Equation (3). This is clearly a typo as he takes into account the minus sign in later parts of his article.

2.2. The H field

Firstly, Deutsch (1955) describes the typical spherical polar coordinates around the star’s axis of rotation, $(r,\theta,\varphi)$ , with rotation axis $\boldsymbol{\omega}$ aligned with the z-axis and origin at the centre of the star. Here r is the radius from the origin, $\theta$ is the angle from the pole and $\varphi$ is the azimuth angle measured counterclockwise from the x-axis.

Let $\boldsymbol{{H}}(r,\theta,\varphi,t)$ be an arbitrary solenoidal function within the star. A solenoidal function is one for which the divergence is zero at all points in the field ( $\nabla\cdot\boldsymbol{{H}}=0$ ). Deutsch (1955) states that a necessary and sufficient condition that $\boldsymbol{{H}}$ be ‘frozen into’ the star is that

(5) \begin{equation}\boldsymbol{{H}}(r,\theta,\varphi,t)=r_0H_r(r,\theta,\lambda)+\theta_0H_\theta(r,\theta,\lambda)+\varphi_0H_\varphi(r,\theta,\lambda),\end{equation}

where $\lambda$ is an azimuthal coordinate measured from a meridian fixed in the star, and $r_0$ , $\theta_0$ and $\phi_0$ are the unit vectors in the respective coordinate directions. That is the time and azimuthal angle $\varphi$ can be replaced by $\lambda$ as the H field and the meridian move with the star. The variable $\lambda$ is defined more clearly in Equation (6) and Figure 2b. $H_r$ , $H_\theta$ and $H_\varphi$ give the magnitude of the field in each basis vector’s direction. The concept of the magnetic field being ‘frozen into’ a medium comes from Alfvén’s Theorem (Roberts 2007).

2.2.1. Vector fields

As $\boldsymbol{{H}}$ is a vector field, at a particular point $(r,\theta,\varphi)$ in space, the field can be represented as a vector with some magnitude in each unit vector’s direction. In Deutsch (1955), the directions are given by the subscript 0 terms, and the magnitudes are given by the functions proceeding the unit vector terms. In polar coordinates, $r_0$ is the radial component, pointing away from the centre of the sphere. The polar component is $\theta_0$ , its direction is tangent to the polar angle. Lastly, $\varphi_0$ is the azimuthal component which has a direction tangent to the azimuthal angle. This notation continues throughout D55. However, we will here on forth elect to use $\,{\hat{\textbf{r}}}$ , $\,\boldsymbol{\hat{\theta}}$ and $\,\boldsymbol{\hat{\varphi}}$ to represent these components of the vector field. A visual representation of the unit vectors of a field in polar coordinates is shown in Figure 1a. Hence, (5) in the notation used for the rest of this review becomes,

\begin{equation*} \boldsymbol{{H}}(r,\theta,\varphi,t)=H_r(r,\theta,\lambda)\,{\hat{\textbf{r}}}+H_\theta(r,\theta,\lambda)\,\boldsymbol{\hat{\theta}}+H_\varphi(r,\theta,\lambda)\,\boldsymbol{\hat{\varphi}}.\end{equation*}

Figure 1. (a) The unit vectors of a vector field at a point on a sphere in $(r, \theta, \varphi)$ . (b) The unit vectors of a vector field at a point on a sphere in $(r,\psi,\nu)$ .

3. The external field in vacuum

It is convenient to define a new set of polar coordinates $(r,\psi,\nu)$ with axis $\boldsymbol{{e}}$ such that the magnetic field is symmetrical about $\boldsymbol{{e}}$ (known as the magnetic field symmetry axis). The axis $\boldsymbol{{e}}$ will be in fixed rotation with the star due to Alfvén’s theorem. This new coordinate system is shown in Figure 2a. The magnetic field symmetry here refers to, for a fixed r and $\psi$ , the magnitude of the field does not change as $\nu$ is changed and the direction does not change relative to the local unit vectors. The angle between the rotation axis $\boldsymbol{\omega}$ and $\boldsymbol{{e}}$ is $\chi$ . This relationship is shown in Figure 2b. The axis $\boldsymbol{{e}}$ has polar angle $\chi$ and azimuthal angle $\omega t$ in $(r,\theta,\varphi)$ , $\boldsymbol{{e}}=(r,\chi,\omega t)$ , at some time t during the rotation of the star. The azimuthal angle $\lambda$ mentioned earlier is measured from $\boldsymbol{{e}}$ along $\varphi$ , so that

(6) \begin{equation} \lambda=\varphi-\omega t.\end{equation}

Figure 2. The relation between the two coordinate systems $(r,\theta,\varphi)$ and $(r,\psi,\nu)$ . x, y, and z are the Cartesian coordinates. (a) The coordinate system $(r,\psi,\nu)$ portrayed on a sphere. (b) The star’s rotation axis $\boldsymbol{\omega}$ and the magnetic field’s symmetry axis $\boldsymbol{{e}}$ . $\chi$ is the angle of inclination between $\boldsymbol{\omega}$ and $\boldsymbol{{e}}$ , and $\lambda$ is measured from the projection of $\boldsymbol{{e}}$ onto the xy-plane.

We shall now consider a magnetic field which within the star has the form,

(7) \begin{equation} \boldsymbol{{H}}=R_1(r)S_1(\psi)\,{\hat{\textbf{r}}}+R_2(r)S_2(\psi)\,\boldsymbol{\hat{\psi}}+R_3(r)S_3(\psi)\,\boldsymbol{\hat{\nu}},\end{equation}

where $R_i$ and $S_i$ are arbitrary functions of r and $\psi$ respectively, with $i=1,2$ , or 3. Note that $\boldsymbol{{H}}$ has no dependence on $\nu$ for the magnitude of the vectors and hence the magnitude is symmetrical about axis $\boldsymbol{{e}}$ as stated earlier. In this situation Deutsch (1955) assumes the interior magnetic field to be a dipole about $\boldsymbol{{e}}$ by taking $S_1(\psi)=\cos\psi$ .

3.1. Matrix transformation between vector field basis vectors

A vector transformation between the two coordinate systems mentioned is required whilst working out the boundary conditions for the fields in both coordinate systems. It is convenient to figure out the transformation now, then approach the boundary conditions. Define the transformation matrix $R:(\,{\hat{\textbf{r}}},\,\boldsymbol{\hat{\theta}},\,\boldsymbol{\hat{\varphi}})\to(\,{\hat{\textbf{r}}},\,\boldsymbol{\hat{\psi}},\,\boldsymbol{\hat{\nu}})$ which rotates between the two unit vector systems at a point in space. Figure 1a and b show the two sets of basis vectors at the same point. It is clear that $\,{\hat{\textbf{r}}}$ points radially out, normal to the sphere’s surface, in both figures. Hence, at the same point in space $\,{\hat{\textbf{r}}}$ is the same in both coordinate systems. However, the angular components undergo a rotation as $\,\boldsymbol{\hat{\theta}}, \,\boldsymbol{\hat{\varphi}}, \,\boldsymbol{\hat{\psi}}$ and $\,\boldsymbol{\hat{\nu}}$ all lie in a plane tangent to the sphere. The rotation between these unit vectors will be a two dimensional rotation about $\,{\hat{\textbf{r}}}$ . Let the angle of this rotation be $\xi$ , the transformation matrix R is then

(8) \begin{equation} R= \left(\begin{array}{c@{\quad}c@{\quad}c} 1 & 0 & 0\\ \\[-12pt] 0 & \cos\xi & -\sin\xi \\ \\[-12pt] 0 & \sin\xi & \cos\xi \end{array}\right).\end{equation}

This matrix transformation takes a basis of vectors at some point in space from $(\,{\hat{\textbf{r}}},\,\boldsymbol{\hat{\theta}},\,\boldsymbol{\hat{\varphi}})$ and rotates them anticlockwise by an angle $\xi$ with $\,{\hat{\textbf{r}}}$ as axis to the unit vectors $(\,{\hat{\textbf{r}}},\,\boldsymbol{\hat{\psi}},\,\boldsymbol{\hat{\nu}})$ . To reverse this transformation the angle of rotation is replaced by $-\xi$ .

The angle $\xi$ can be found in terms of already defined coordinate parameters and values. This is done by considering triangles on the surface of a sphere. The spherical triangle for this situation can be found in Figure F.1b in Appendix F. Then from solutions of spherical triangles (F2) and (F3) the following relationships can be found,

(9) \begin{equation} \begin{aligned} \sin\xi&=\frac{\sin\chi\sin\lambda}{\sin\psi},\\ \\[-10pt] \cos\xi&=\frac{\cos\chi-\cos\psi\cos\theta}{\sin\psi\sin\theta}. \end{aligned}\end{equation}

An interesting relationship that is obtained from the standard spherical triangle relation (F1) is

(10) \begin{equation} \cos\psi=\cos\chi\cos\theta+\sin\chi\sin\theta\cos\lambda,\end{equation}

which will be required momentarily.

3.2. Boundary conditions of the fields

As the surface of the star is a rigid boundary between the interior of the star and free space, the necessary components of the external electromagnetic fields must satisfy the boundary conditions at $r=a$ (see Griffiths 2017). These boundary conditions are set by the internal electromagnetic fields. From Equations (7.61) and (7.62) of Griffiths (2017), the magnetic fields orthogonal to the surface must be equal and the electric fields tangent to the surface must be equal. Hence, using (7), the first condition at the surface of the star that must be satisfied is,

\begin{equation*} H_r^{\text{ext}}=H_r^{\text{int}}=R_1(a)S_1(\psi).\end{equation*}

Next we look at the internal electric field. To do this we use Ohm’s law from Griffiths (2017), Equation (7.2), that relates the electric field and current density, $\boldsymbol{{J}}$ and is given by,

(11) \begin{equation}\boldsymbol{{E}}+\mu_0\boldsymbol{{V}}\times\boldsymbol{{H}}=\frac{\boldsymbol{{J}}}{\sigma}\end{equation}

where $\sigma$ is the conductivity and $\boldsymbol{{V}}$ is the local velocity of rotation of the star. Assuming $\sigma\rightarrow \infty$ gives

(12) \begin{equation} \boldsymbol{{E}}=-\mu_0\boldsymbol{{V}}\times\boldsymbol{{H}}. \end{equation}

The velocity in $(\,{\hat{\textbf{r}}},\,\boldsymbol{\hat{\theta}},\,\boldsymbol{\hat{\varphi}})$ is

(13) \begin{equation} \boldsymbol{{V}}=\boldsymbol{\omega}\times\boldsymbol{{r}}=r\omega\sin\theta\,\boldsymbol{\hat{\varphi}}.\end{equation}

However, the vector field (7) must be converted from $(\,{\hat{\textbf{r}}},\,\boldsymbol{\hat{\psi}},\,\boldsymbol{\hat{\nu}})$ to $(\,{\hat{\textbf{r}}},\,\boldsymbol{\hat{\theta}},\,\boldsymbol{\hat{\varphi}})$ so that $\boldsymbol{{V}}$ and $\boldsymbol{{H}}$ are defined on the same basis. This is done by rotating (7) using (8). However, as the reverse transformation is needed, the signs of $R_{2,3}$ and $R_{3,2}$ elements are switched.

\begin{align*} R(-\xi)\boldsymbol{{H}}&=R_1(r)S_1(\psi)\,{\hat{\textbf{r}}}\\ \\[-10pt] &\quad +\left[R_2(r)S_2(\psi) \cos\xi + R_3(r)S_3(\psi) \sin\xi \right]\,\boldsymbol{\hat{\theta}}\\ \\[-10pt] &\quad +\left[-R_2(r)S_2(\psi) \sin\xi + R_3(r)S_3(\psi) \cos\xi \right]\,\boldsymbol{\hat{\varphi}}. \end{align*}

Applying the previous result to (12) with (13) and using Figure F.1b of Appendix F we can find the E field inside the star,

\begin{align*} \boldsymbol{{E}}&=-\mu_0\boldsymbol{{V}}\times(R(-\xi)\boldsymbol{{H}})\\ \\[-8pt] &=\mu_0 r\omega\sin\theta\left(R_2(r)S_2(\psi) \cos\xi + R_3(r)S_3(\psi) \sin\xi \right)\,{\hat{\textbf{r}}}\\ \\[-8pt] &\quad -\omega\mu_0rR_1(r)S_1(\psi)\sin\theta\,\boldsymbol{\hat{\theta}}\\ \\[-8pt] &\quad +0\,\boldsymbol{\hat{\varphi}}, \end{align*}

and evaluating this at $r=a$ gives the tangential components as

\begin{equation*} \begin{aligned} E_\theta&=-\omega\mu_0aR_1(a)S_1(\psi)\sin\theta,\\ \\[-9pt] E_\varphi&=0. \end{aligned}\end{equation*}

All but the third boundary condition match what Deutsch (1955) obtained. There is a typographical error in Deutsch (1955)’s boundary conditions (Deutsch’s Equation (10)), the third condition should be written as $E_\varphi^{\text{ext}}=E_\varphi^{\text{int}}=0$ rather than $E_r^{\text{ext}}=E_\varphi^{\text{int}}=0$ .

A more useful form of the boundary conditions can be found. In particular, the field equations are solved for when $S_1=\cos\psi$ . So by recalling (10), the boundary conditions are then

(14) \begin{align} H_r^{\text{ext}}=H_r^{\text{int}}&=R_1(a)(\cos\chi\cos\theta+\sin\chi\sin\theta\cos\lambda), \nonumber\\ E_\theta^{\text{ext}}=E_\theta^{\text{int}}&=-\omega\mu_0aR_1(a)(\cos\chi\cos\theta+\sin\chi\sin\theta\cos\lambda)\sin\theta \nonumber\\ &=-\frac{1}{2}\omega\mu_0aR_1(a) (\sin{2\theta}\cos\chi+(1-\cos{2\theta})\sin\chi\cos\lambda), \nonumber\\ E_\varphi^{\text{ext}}=E_\varphi^{\text{int}}&=0. \end{align}

From the boundary conditions in the $(r,\theta,\phi)$ coordinate system, there is a static (time-independent) component and a time-dependent component for the fields. The terms in (14) which relate to the time-dependent component have a dependency on $\lambda$ , which depends on time from (6).

4. Derivation of general field equations

The general field equations for the VRD solution external to the star are provided in the Appendix of Deutsch (1955). However, with minimal explanation given in Deutsch (1955) about the derivation of the field equations, we wish to show how they are found. We use results found in Chapters VII and VIII of Stratton (1941) for this section. Chapter VII proves beneficial in supplying general solutions to the vector wave equation, which must be satisfied by the time evolving electromagnetic fields around the neutron star. We apply these solutions to the neutron star situation by firstly considering the propagation factor of the waves outside the star. In Stratton (1941), on p. 392, the propagation factor (complex wavenumber) is given as,

(15) \begin{equation} k^2=\epsilon\mu\omega^2+i\sigma\mu\omega,\end{equation}

and as we are now considering the fields external of the star, in a vacuum, $\epsilon=\epsilon_0$ , $\mu=\mu_0$ , $\sigma=0$ and $c=1/\sqrt{\mu_0\epsilon_0}$ . Hence,

(16) \begin{equation} k = \omega/c.\end{equation}

Then for convenience define,

(17) \begin{equation} \begin{gathered} \rho=kr=(\omega/c)r,\\ \\[-8pt] \alpha=ka=(\omega/c)a, \end{gathered}\end{equation}

where $\omega$ is the rotational velocity of the star, c the speed of light, and a is the radius of the star. Here $\rho$ and $\alpha$ are unitless numbers. A trivial relationship that will be required is

\begin{equation*} \frac{\partial}{\partial\rho}=\frac{\omega}{c}\frac{\partial}{\partial r}.\end{equation*}

Next consider the EM fields surrounding the star in two parts, a static time-independent component and a time-dependent component, where the complete field will be the linear combination of the two components (Griffiths 2017). The time-independent component will be axially symmetric about the axis $\boldsymbol{\omega}$ and hence have no $\varphi$ dependence ( $m=0$ , m will be defined clearer shortly). Whilst the time-dependent component will be non-axially symmetric and hence have a dependence on $\varphi$ ( $m\neq0$ ). This can be expressed as

\begin{equation*}\begin{gathered} \boldsymbol{{H}} = \boldsymbol{{H}}_{\text{static}}(r,\theta)+\boldsymbol{{H}}_{\text{time-dependent}}(r,\theta, \lambda),\\ \\[-8pt] \boldsymbol{{E}} = \boldsymbol{{E}}_{\text{static}}(r,\theta)+\boldsymbol{{E}}_{\text{time-dependent}}(r,\theta, \lambda).\end{gathered}\end{equation*}

We first consider the time-dependent component in the following section on the wave solution. Then use a multipole expansion to provide the static component.

4.1. Wave solution

Due to the rotating nature of the magnetic field symmetry axis around the star’s axis of rotation, the vector fields $\boldsymbol{{E}}$ and $\boldsymbol{{H}}$ satisfy the same vector differential wave equation from Stratton (1941), Section 7.1, Equation (1), which is given as

(18) \begin{equation} \nabla^2\boldsymbol{{C}}-\mu\epsilon\frac{\partial^2\boldsymbol{{C}}}{\partial t^2}-\mu\sigma\frac{\partial \boldsymbol{{C}}}{\partial t}=0,\end{equation}

where $\boldsymbol{{C}}$ represents either vector $\boldsymbol{{E}}$ or $\boldsymbol{{H}}$ . We can express $\boldsymbol{{E}}$ and $\boldsymbol{{H}}$ as a series of solutions $\textbf{M}_n$ and $\textbf{N}_n$ which satisfy the vector differential wave Equation (18) for $\boldsymbol{{C}}$ . From Stratton (1941), Section 7.1, Equation (12) this representation is,

(19) \begin{equation} \begin{gathered} \boldsymbol{{E}}=-\sum_n\left(a_n\textbf{M}_n+b_n\textbf{N}_n\right),\\ \\[-8pt] \boldsymbol{{H}}=-\frac{k}{i\omega\mu_0}\sum_n\left(a_n\textbf{N}_n+b_n\textbf{M}_n\right) \end{gathered}\end{equation}

where $\textbf{M}_n$ and $\textbf{N}_n$ have more recently been referred to as the Magnetic and Electric harmonics respectfully. In particular, the functions are expressed as $\textbf{M}_n=\nabla\times\boldsymbol{{a}}\psi$ and $\textbf{N}_n=(1/k)\nabla\times\textbf{M}_n$ , where $\boldsymbol{{a}}$ is a constant vector and $\psi$ is a scalar function satisfying $\nabla^2\psi+k^2\psi=0$ , which is known as the Helmholtz equation. The derivation for the general solutions of $\textbf{M}_n$ and $\textbf{N}_n$ is provided in Stratton (1941) and the results are quoted below. These solutions will also contain a time dependence. Due to the linearity of the wave equation, the time dependency can be split off such that $\textbf{M}=\textbf{m}e^{-i\omega t}$ and $\textbf{N}=\textbf{n}e^{-i\omega t}$ without loss of generality. Then from Stratton (1941), Section 7.11, Equations (11) and (12), and rewriting them in a slightly more useful form with $\rho$ instead of kr,

(20) \begin{equation} \begin{aligned} \textbf{m}_{\genfrac{}{}{0pt}{}{e}{o} {m\,n}}=&\mp \frac{m}{\sin \theta}h^{(1)}_{n}(\rho) P_{n}^{m}(\cos \theta) \genfrac{}{}{0pt}{}{\sin}{\cos}(m \varphi)\,\boldsymbol{\hat{\theta}}\\ \\[-8pt] &-h^{(1)}_{n}(\rho)\frac{\partial P_{n}^{m}}{\partial \theta}\genfrac{}{}{0pt}{}{\cos}{\sin}(m\varphi)\,\boldsymbol{\hat{\varphi}}, \end{aligned}\end{equation}

(21) \begin{equation} \begin{aligned} \textbf{n}_{\genfrac{}{}{0pt}{}{e}{o} {m\,n}}=&\frac{n(n+1)}{\rho} h^{(1)}_{n}(\rho) P_{n}^{m}(\cos \theta) \genfrac{}{}{0pt}{}{\cos }{\sin }(m \varphi)\,{\hat{\textbf{r}}}\\ \\[-10pt] &+\frac{1}{\rho} \frac{\partial}{\partial \rho}\left[\rho\,h^{(1)}_{n}(\rho)\right] \frac{\partial}{\partial \theta} P_{n}^{m}(\cos \theta) \genfrac{}{}{0pt}{}{\cos}{\sin}(m \varphi)\,\boldsymbol{\hat{\theta}}\\ \\[-10pt] &\mp \frac{m}{\rho \sin \theta} \frac{\partial}{\partial \rho}\left[\rho\, h^{(1)}_{n}(\rho)\right] P_{n}^{m}(\cos \theta)\genfrac{}{}{0pt}{} {\sin}{\cos}(m \varphi)\,\boldsymbol{\hat{\varphi}}, \end{aligned}\end{equation}

where $h^{(1)}_n$ is the spherical Bessel function of the third kind (refer to Appendix C), $P_{n}^{m}(\cos \theta)$ are the associated Legendre polynomials (refer to Appendix D), and m and n are associated to a particular solution to the Helmholtz equation. In particular, m defines the azimuthal dependence of the E and H fields. Whilst n controls the polar dependence. These functions (20) and (21) are now known as vector spherical harmonics. The subscript e or o define the even and odd solution, this chooses the sign of some terms and whether the trigonometric function for $\varphi$ is $\sin$ or $\cos$ . From here on $h^{(1)}_n$ will be referred to as $h_n$ . However, Deutsch (1955) does not use the even-odd notation from Stratton (1941) but rather Deutsch (1955) multiplies the function by a complex exponential then takes the real part of the result to give the wave solutions. This can be expressed as an operator that first multiplies a function f(x) by $e^{im\varphi}$ , then takes the real part of the result,

(22) \begin{equation}\begin{aligned} f(x)\genfrac{}{}{0pt}{}{\cos}{\sin}(m \varphi)&=\Re(f(x)e^{im\varphi}),\\ \mp f(x)\genfrac{}{}{0pt}{}{\sin}{\cos}(m \varphi)&=-\Re(f(x)ie^{im\varphi}),\end{aligned}\end{equation}

where f(x) may be complex-valued function and $\Re$ denotes the operator which returns the real part of an expression. The operation (22) essentially selects $\cos(m\varphi)$ or $\sin(m\varphi)$ depending on the complex-valued nature of f(x).

We wish to solve for the field equations in polar coordinates around the rotation axis $\boldsymbol{\omega}$ . Due to the dependence of $\varphi$ in our electromagnetic boundary conditions (14) only being $\cos(\varphi)$ and $\sin(\varphi)$ , we are limited to $m=1$ for our Equations (20) and (21). Otherwise terms which contain $\cos(m\varphi)$ and $\sin(m\varphi)$ , where $m\neq1$ , will appear in the solutions of (19) and cause disagreement of the external and internal electromagnetic fields at the surface of the star. By also comparing the series expression (19) and the boundary conditions (14) we see that if a trigonometric function of $\theta$ does not appear in the boundary conditions, its weighting coefficient will be zero (as will be noticed in the derivations shortly). Hence, it is apparent from inspection of the boundary conditions and the Associated Legendre Polynomials that only $n=1$ and $n=2$ are required. So evaluating (20) and (21) for $m=1$ , $n=1$ and $n=2$ and including the new notation from (22) gives,

(23) \begin{align} \textbf{m}_{1,1} &= i h_1(\rho)e^{i\varphi}\,\boldsymbol{\hat{\theta}} - h_1(\rho)\cos\theta e^{i\varphi}\,\boldsymbol{\hat{\varphi}}, \nonumber\\[3pt] \textbf{m}_{1,2} &= 3 i h_2(\rho)\cos\theta e^{i\varphi}\,\boldsymbol{\hat{\theta}} - 3h_2(\rho)\cos{2\theta}e^{i\varphi}\,\boldsymbol{\hat{\varphi}},\nonumber\\[3pt] \textbf{n}_{1,1} &= 2\frac{h_1(\rho)}{\rho}\sin\theta e^{i\varphi}\,{\hat{\textbf{r}}} + \frac{1}{\rho}(\rho h_1^{\prime}(\rho)+h_1(\rho))\cos\theta e^{i\varphi}\,\boldsymbol{\hat{\theta}}\nonumber\\[3pt] &\quad + \frac{i}{\rho}(\rho h_1^{\prime}(\rho)+h_1(\rho))e^{i\varphi}\,\boldsymbol{\hat{\varphi}},\nonumber\\[3pt] \textbf{n}_{1,2}&= 9\frac{h_2(\rho)}{\rho}\sin{2\theta} e^{i\varphi}\,{\hat{\textbf{r}}} + \frac{3}{\rho}(\rho h_2^{\prime}(\rho)+h_2(\rho))\cos{2\theta} e^{i\varphi}\,\boldsymbol{\hat{\theta}}\nonumber\\[3pt] &\quad + \frac{3i}{\rho}(\rho h_2^{\prime}(\rho)+h_2(\rho))\cos\theta e^{i\varphi}\,\boldsymbol{\hat{\varphi}}, \end{align}

where we have neglected the operator that extracts the real part of the result. This will be added during the calculation for the coefficients in (19).

4.2. Multipole expansion

A static field symmetric about the rotation axis can be expressed using a multipole expansion. Upon inspection of the static terms in the boundary conditions (14) and the multipole expansions given in Griffiths (2017) it can be seen, with some working, that only a magnetic dipole and electric quadrupole are present in this situation based on the conditions that multipole terms are linearly independent of each other and that the boundary conditions will uniquely determine the fields (electromagnetism uniqueness theorem, see Griffiths 2017).

4.2.1. Magnetic dipole

We can begin by taking a magnetic dipole and check whether it satisfies the H field boundary condition from (14). From Equation (5.88) of Griffiths (2017), a magnetic dipole is given as,

(24) \begin{equation} \boldsymbol{{H}}_{\text{dip}}(\boldsymbol{{r}})=\frac{m}{4\pi r^3}\left(2\cos\theta\,{\hat{\textbf{r}}}+\sin\theta\,\boldsymbol{\hat{\theta}}\right),\end{equation}

where m is the magnetic dipole moment of the field. Now compare the radial component of (24) to the static term of the radial H field boundary condition from (14).

\begin{equation*} \frac{m}{4\pi a^3}2\cos\theta = R_1(a)\cos\chi\cos\theta.\end{equation*}

We see that the function depending on $\theta$ in the boundary condition is satisfied fully by the dipole term. Hence, a magnetic dipole is the only static term present due to the boundary conditions. This implies the ratio $m/(4\pi)$ will satisfy,

\begin{equation*} \frac{m}{4\pi} = \frac{1}{2}a^3R_1(a)\cos\chi.\end{equation*}

Therefore, the static component of the magnetic field will be

(25) \begin{equation} \boldsymbol{{H}}_{\text{static}}(\boldsymbol{{r}})=\frac{1}{2}\frac{a^3}{ r^3}R_1(a)\cos\chi\left(2\cos\theta\,{\hat{\textbf{r}}}+\sin\theta\,\boldsymbol{\hat{\theta}}\right).\end{equation}

4.2.2. Electric quadrupole

We begin by assuming the solution for the static electric field to be a quadrupole, after which we show that it is the only component required to satisfy the boundary conditions due to the dependency on $\theta$ . From Griffiths (2017), Equation (3.65), the quadrupole potential term is

\begin{equation*} V_{\text{quad}}(r)=\frac{B_2}{r^3}\frac{1}{2}\left(3\cos^2\theta-1\right).\end{equation*}

where $B_2$ is the weighting coefficient from the series expansion in Griffiths (2017).

From Equation (2.23) of Griffiths (2017), the electric field is the negative gradient of the potential. Hence, the polar component of the E field is

\begin{equation*} E_{\text{quad},\,\theta}=-\frac{1}{r}\frac{\partial V}{\partial\theta} = \frac{B_2}{r^4}3\cos\theta\sin\theta.\end{equation*}

The functions depending on the variable $\theta$ match those of the boundary condition (14). Therefore, the boundary condition only requires the static component of the electric field to be a quadrupole. Comparing the above result to (14) we can uniquely determine $B_2$ for the situation,

\begin{equation*}\begin{gathered} \frac{B_2}{a^4}3\cos\theta\sin\theta = -\omega\mu_0 a R_1(a)\cos\chi\cos\theta\sin\theta,\\ \\[-9pt] \implies B_2 = -\frac{1}{3}\omega\mu_0 a^5 R_1(a) \cos\chi.\end{gathered}\end{equation*}

and on applying the gradient function to the potential in spherical coordinates (see Appendix E),

\begin{equation*} \begin{aligned} E_{\text{quad},\,r}&=-\frac{\partial V}{\partial r} =-\frac{1}{2}\omega\mu_0 a\frac{a^4}{r^4}R_1(a)\cos\chi\left(3\cos^2\theta-1\right),\\ \\[-9pt] E_{\text{quad},\,\theta}&=-\frac{1}{r}\frac{\partial V}{\partial\theta} = - \omega \mu_0 a R_1(a)\frac{a^4}{r^4}\cos\chi\cos\theta\sin\theta,\\ \\[-9pt] E_{\text{quad},\,\varphi}&=-\frac{1}{r\sin\theta}\frac{\partial V}{\partial\varphi}=0. \end{aligned}\end{equation*}

The trigonometric identities $3\cos^2\theta-1=\frac{1}{2}\left(3\cos{2\theta}+1\right)$ and $\cos\theta\sin\theta=\frac{1}{2}\sin 2\theta,$ can be used to write the quadrupole terms in a similar form as Deutsch (1955),

(26) \begin{equation}\begin{aligned} E_{\text{quad},\,r}&=-\frac{1}{4}\omega\mu_0 a\frac{a^4}{r^4}R_1(a)\cos\chi\left(3\cos2\theta+1\right),\\ \\[-9pt] E_{\text{quad},\,\theta}&= - \frac{1}{2}\omega \mu_0 a R_1(a)\frac{a^4}{r^4}\cos\chi\sin2\theta,\\ E_{\text{quad},\,\varphi}&=0.\end{aligned}\end{equation}

4.3. Derivation of general field equations

We are now in a position to derive the equations in the Appendix of Deutsch (1955). The general field equations will be the sum of the multipoles derivations (25) and (26), and the series expressions (19),

(27) \begin{equation} \begin{gathered} \boldsymbol{{H}}=\boldsymbol{{H}}_{\text{dip}}-\frac{1}{ic\mu_0}(a_1\textbf{n}_{1,1}+a_2\textbf{n}_{1,2}+b_1\textbf{m}_{1,1}+b_2\textbf{m}_{1,2})e^{-i\omega t},\\ \\[-9pt] \boldsymbol{{E}}=\boldsymbol{{E}}_{\text{quad}}-(a_1\textbf{m}_{1,1}+a_2\textbf{m}_{1,2}+b_1\textbf{n}_{1,1}+b_2\textbf{n}_{1,2})e^{-i\omega t},\\ \end{gathered} \end{equation}

where $a_1$ , $a_2$ , $b_1$ and $b_2$ are constants to be determined by the boundary conditions (14). These constants may not necessarily be real valued.

4.3.1. H field radial component

By substituting the radial components from (23) and (25) into (27) we have,

\begin{equation*}\begin{aligned} H_r=&\frac{a^3}{ r^3}R_1(a)\cos\chi\cos\theta\\ &- \frac{1}{ic\mu_0}\left(2a_1\frac{h_1(\rho)}{\rho}\sin\theta e^{i\varphi}+9a_2\frac{h_2(\rho)}{\rho}\sin{2\theta} e^{i\varphi}\right)e^{-i\omega t},\end{aligned}\end{equation*}

which, when $r=a$ ( $\rho=\alpha$ ), will be equal to the boundary condition on $H_r$ from (14). We can use $e^{i\varphi}e^{-i\omega t}=e^{i(\varphi-\omega t)}=e^{i\lambda}$ to write the exponent in terms of $\lambda$ and that $e^{i\lambda}=\cos\lambda+i\sin\lambda$ . The constants can then be determined by matching trigonometric functions. For $\sin\theta$ ,

\begin{equation*} \Re\left(-\frac{1}{ic\mu_0}a_1\frac{h_1(\alpha)}{\alpha}e^{i\lambda}\right)=R_1(a)\sin\chi\cos\lambda,\end{equation*}

where the $\Re$ operator has come from (22). Hence, $a_1$ must be purely imaginary to cancel out i in the bottom of the fraction so that $\cos\lambda$ is obtained instead of $\sin\lambda$ .

(28) \begin{gather} -\frac{2}{ic\mu_0}a_1\frac{h_1(\alpha)}{\alpha}\cos\lambda=R_1(a)\sin\chi\cos\lambda, \nonumber\\ \implies a_1=-\frac{R_1(a)ic\mu_0\sin\chi}{2}\frac{\alpha}{h_1(\alpha)} .\end{gather}

As there is no $\sin{2\theta}$ term in the boundary condition, we will have $a_2=0$ . The gives the radial component of the H field as the real part of

\begin{equation*} H_r=R_1(a)\left(\frac{a^3}{r^3}\cos\chi\cos\theta+\frac{\alpha}{h_1(\alpha)}\frac{h_1(\rho)}{\rho}\sin\chi\sin\theta e^{i\lambda}\right).\end{equation*}

4.3.2. E field polar component

Now considering the $E_\theta$ component which can be found by combing the $\,\boldsymbol{\hat{\theta}}$ components of results (23) and (26) into (27),

\begin{equation*} \begin{aligned} E_\theta=&- \frac{1}{2}\omega \mu_0 a R_1(a)\frac{a^4}{r^4}\cos\chi\sin2\theta - a_1 i h_1(\rho) e^{i\varphi} e^{-i\omega t}\\ \\[-9pt] &-b_1\frac{1}{\rho}(\rho h_1^{\prime}(\rho)+h_1(\rho))\cos\theta e^{i\varphi}e^{-i\omega t}\\ &-b_2\frac{3}{\rho}(\rho h_2^{\prime}(\rho)+h_2(\rho))\cos{2\theta} e^{i\varphi}e^{-i\omega t}, \end{aligned}\end{equation*}

and when $r=a$ this is equal to the $E_\theta$ component from (14). Once again, matching functions of $\theta$ and using $e^{i\varphi}e^{-i\omega t}=e^{i\lambda}$ , we find that for the term excluding any function of $\theta$ ,

\begin{equation*}\begin{gathered} \Re\left(-\frac{c\mu_0R_1(a)}{2}\frac{\alpha}{h_1(\alpha)}h_1(\alpha)\sin\chi e^{i\lambda}\right)\\ =-\frac{1}{2}\omega\mu_0aR_1(a)\sin\chi\cos\lambda,\end{gathered}\end{equation*}

then it is convenient to change $\alpha$ to $(\omega/c)a$ ,

\begin{equation*} -\frac{\omega\mu_0aR_1(a)}{2}\sin\chi\cos\lambda=-\frac{1}{2}\omega\mu_0aR_1(a)\sin\chi\cos\lambda,\end{equation*}

the LHS and RHS match, which validates the result for $a_1$ found in (28). As the boundary condition does not contain any $\cos\theta$ terms,

(29) \begin{equation}\begin{gathered} -b_1\frac{1}{\rho}(\rho h_1^{\prime}(\rho)+h_1(\rho))\cos\theta e^{i\lambda}=0,\\ \implies b_1=0.\end{gathered}\end{equation}

Lastly, the $\cos{2\theta}$ terms,

\begin{gather*} \Re\left(-b_2\frac{3}{\alpha}(\alpha h_2^{\prime}(\alpha)+h_2(\alpha))\cos{2\theta} e^{i\lambda}\right)\\ =-\frac{1}{2}\omega\mu_0aR_1(a)\cos{2\theta}\sin\chi\cos\lambda,\end{gather*}

again we want $\cos\lambda$ instead of $\sin\lambda$ so $b_2$ must be purely real. Hence,

(30) \begin{gather}-b_2\frac{3}{\alpha}(\alpha h_2^{\prime}(\alpha)+h_2(\alpha))\cos{2\theta}\cos\lambda \nonumber \\=\frac{1}{2}\omega\mu_0aR_1(a)\cos{2\theta}\sin\chi\cos\lambda, \nonumber \\\implies b_2=-\frac{1}{2}\omega\mu_0aR_1(a)\frac{1}{3}\frac{\alpha}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)}\sin\chi .\end{gather}

All the results for the constants thus far lead to,

(31) \begin{equation}\begin{aligned} E_\theta=&\frac{1}{2}\omega\mu_0 a R_1(a) \left[-\frac{a^4}{r^4}\cos\chi\sin{2\theta}\right.\\ &+ \left.\left(\frac{\alpha}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)}\frac{\rho h_2^{\prime}(\rho) + h_2(\rho)}{\rho}\cos{2\theta} - \frac{h_1(\rho)}{h_1(\alpha)}\right)\sin\chi e^{i\lambda}\right].\end{aligned}\end{equation}

Deutsch (1955) has a typographical error in his Appendix for Equation (31). The $\frac{\alpha}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)}\frac{\rho h_2^{\prime}(\rho) + h_2(\rho)}{\rho}$ term is instead written as $\frac{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)}{\alpha}\frac{\rho}{\rho h_2^{\prime}(\rho)+h_2(\rho)}$ .

4.3.3. E field azimuthal component

The final boundary condition is on the $E_\varphi$ component. We combine the results from above with azimuthal component of (27), which yields

\begin{equation*} \begin{aligned} E_\varphi=&-\left(-a_1 h_1(\rho)\cos\theta+b_2\frac{3i}{\rho}(\rho h_2^{\prime}(\rho)+h_2(\rho))\cos\theta\right)e^{i\lambda}\\ =&\frac{1}{2}\omega\mu_0aR_1(a) \left(\frac{\alpha}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)}\frac{\rho h_2^{\prime}(\rho)+h_2(\rho)}{\rho} - \frac{h_1(\rho)}{h_1(\alpha)}\right)\\&\quad i\sin\chi\cos\theta e^{i\lambda}, \end{aligned} \end{equation*}

and clearly at $r=a$ ( $\rho=\alpha$ ) this reduces to zero and matches the required boundary condition $E_\varphi(a,\theta,\varphi)=0$ , implying that the values for the constants are consist with all the boundary conditions in this situation.

4.3.4. Remaining vector field components

Noting that the expressions for the constants $a_1$ and $b_2$ can be rewritten using $\alpha=(\omega/c)a$ , the field components are given by the real part of the following,

\begin{equation*} \begin{aligned} H_\theta&=\frac{1}{2}\frac{a^3}{ r^3}R_1(a)\cos\chi\sin\theta - \frac{1}{i c\mu_0}\left(a_1\frac{\rho h_1^{\prime}(\rho)+h_1(\rho)}{\rho}\cos\theta + b_2(3ih_2(\rho)\cos\theta)\right)e^{i\lambda}\\ &=\frac{1}{2}R_1(a)\left[\frac{a^3}{r^3}\cos\chi\sin\theta + \left(\frac{\alpha}{h_1(\alpha)}\frac{\rho h_1^{\prime}(\rho)+h_1(\rho)}{\rho}\right.\right.\\ &\quad\left.\left. + \frac{1}{c}\omega a\frac{\alpha}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)} h_2(\rho)\right)\sin\chi\cos\theta e^{i\lambda}\right]\\ &=\frac{1}{2}R_1(a)\left[\frac{a^3}{r^3}\cos\chi\sin\theta+\left(\frac{\alpha}{h_1(\alpha)}\frac{\rho h_1^{\prime}(\rho)+h_1(\rho)}{\rho} +\frac{\alpha^2h_2(\rho)}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)} \right)\sin\chi\cos\theta e^{i\lambda}\right],\\ \\[-9pt] H_\varphi&=-\frac{1}{i c\mu_0}\left(a_1\frac{i}{\rho}(\rho h_1^{\prime}(\rho)+h_1(\rho)) - 3 b_2 h_2(\rho)\cos{2\theta}\right)e^{i\lambda}\\ &=\frac{1}{2}R_1(a)\left(i\frac{\alpha}{h_1(\alpha)}\frac{\rho h_1^{\prime}(\rho)+h_1(\rho)}{\rho}\sin\chi-\frac{1}{i}\frac{\alpha^2}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)}h_2(\rho)\sin\chi\cos{2\theta}\right) e^{i\lambda}\\ &=\frac{1}{2}R_1(a)\left(\frac{\alpha}{h_1(\alpha)}\frac{\rho h_1^{\prime}(\rho)+h_1(\rho)}{\rho}+\frac{\alpha^2}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)}h_2(\rho)\cos{2\theta}\right)i\sin\chi e^{i\lambda},\\ \\ E_r &= -\frac{1}{4}\omega\mu_0 a\frac{a^4}{r^4}R_1(a)\cos\chi\left(3\cos2\theta+1\right)-9b_2\frac{h_2(\rho)}{\rho}\sin{2\theta} e^{i\lambda}\\ &=\frac{1}{2}\omega\mu_0aR_1(a)\left(-\frac{1}{2}\frac{a^4}{r^4}\cos\chi(3\cos{2\theta}+1)+3\frac{\alpha}{\alpha h_2^{\prime}(\alpha)+h_2(\alpha)}\frac{h_2(\rho)}{\rho}\sin\chi\sin{2\theta} e^{i\lambda}\right). \end{aligned} \end{equation*}

The components of the general field equations are gathered in Appendix A for easier reference. A plot^a of the magnetic field is shown in Figure 3.

Figure 3. The magnetic field from (A1) using typical star parameters is plotted as streamlines showing the direction and the log magnitude of the magnetic field strength represented by the colour. Directly vertical is the star’s rotation axis, and the magnetic field axis is set as $\chi=\pi/4$ . Note the dashed red circle is the light cylinder, where $\rho=1$ or $r=c/\omega$ .

5. The radiation fields of the star

Far from a star, the electromagnetic fields take on their radiation form, allowing us to find the fields which ‘leak’ out of the magnetosphere and propagate away from the electromagnetic source. To find the radiation form of the electromagnetic fields we consider the limit in which $r\gg c/\omega$ (Griffiths 2017). Hence, from (17) if $r\gg c/\omega$ then $\rho\gg1$ . Then the limit of the spherical Bessel functions when $\rho$ is large (C4) can be substituted into (A1). However, we find that $\alpha=(\omega/c)a\ll1$ , which means $h_1(\alpha)$ can be approximated using (C3). Using $c\simeq3\times10^{8}\textrm{ms}^{-1}$ and the typical star values given by Deutsch (1955),

\begin{align*} a\simeq10^{9}\textrm{m},\quad \omega\simeq10^{-5}\textrm{rad}\,\textrm{s}^{-1}, \\[-25pt] \end{align*}

we get $\alpha\simeq 3\times 10^{-5}\ll 1$ . While a typical neutron star has

\begin{equation*} a\simeq10^{4}\textrm{m},\quad \omega\simeq2\pi\,\textrm{rad}\,\textrm{s}^{{-1}}, \end{equation*}

which gives $\alpha\simeq 2\times 10^{-4}\ll 1$ . We can also express this as a scaled equation,

\begin{equation*} \alpha = 2\times 10^{-4}\left(\frac{a}{10^{4}\textrm{m}}\right)\left(\frac{\omega}{2\pi\,\textrm{rad}\,\textrm{s}^{{-1}}}\right)\end{equation*}

The fastest recorded neutron star has $a\simeq1.6\times10^{4}\textrm{m}$ and $\omega\simeq{4.5}\times10^{3}\textrm{rad}\,\textrm{s}^{{-1}}$ (Hessels et al. 2006). We see from these values that $\alpha\simeq 0.2 < 1$ . Though this limit is not as strong as slower neutron stars, it shows an appreciable upper limit that is still less than one. Suggesting that taking the approximation $\alpha\ll1$ in the following derivations is valid for slower stars but can still apply to faster stars, but with some error.

5.1. Radiation form derivation

To derive the radiation form we take the limit in which $r\gg c/\omega$ . In taking this limit, we only consider terms that will be of the largest magnitude for each electromagnetic field component, neglecting terms containing high powers of r in their denominator. Recall $\Re$ denotes the real part of an expression from (22) and $\lambda=\varphi-\omega t$ from (6). Working firstly on the H field’s radial component, after neglecting smaller magnitude terms and substituting the limit forms of $h_1(\rho)$ and $h_1(\alpha)$ in, $H_r$ becomes

\begin{equation*}\begin{aligned} H_r&\simeq\Re\left({R_1(a)\frac{e^{i\rho}}{\rho^2}}\frac{\alpha^3}{i}\sin\chi\sin\theta e^{i\lambda}\right)\\ &=\Re\left(R_1(a)\frac{\omega}{c}\frac{a^3}{r^2}\frac{1}{i}\sin\chi\sin\theta e^{i(\lambda+\rho)}\right)\\ &=\Re\left(R_1(a)\frac{\omega}{c}\frac{a^3}{r^2}\frac{1}{i}\sin\chi\sin\theta\left[\cos(\lambda+\rho)+i\sin(\lambda+\rho)\right]\right)\\ &=\frac{\omega}{c}a^3R_1(a)\frac{1}{r^2}\sin\chi\sin\theta\sin\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right],\end{aligned}\end{equation*}

where $\lambda+\rho=\varphi-\omega t + \omega r/c=\omega(r/c-t)+\varphi$ . Repeating the same limit and substituting the approximations in for the polar component yields,

\begin{equation*}\begin{aligned} H_\theta&\simeq\Re\left(\frac{1}{2}R_1(a)\left[\left(\frac{\alpha^2}{\alpha\frac{9i}{\alpha^4}-\frac{3i}{\alpha^3}}\right)\frac{ie^{i\rho}}{\rho}+\frac{\alpha^3}{-i}\left(\frac{-ie^{i\rho}}{\rho}-\frac{e^{i\rho}}{\rho^2}\right)\right]\sin\chi\cos\theta e^{i\lambda}\right)\\ &=\Re\left(\frac{1}{2}R_1(a)\left[\frac{\alpha^5}{6\rho}+\frac{\alpha^3}{\rho}+\frac{\alpha^3}{i\rho^2}\right]\sin\chi\cos\theta e^{i(\lambda+\rho)}\right). \end{aligned}\end{equation*}

The term with the largest magnitude in the line above is the $\alpha^3/\rho$ term as $\alpha\ll1$ and $\rho\gg1$ . So,

\begin{equation*} \begin{aligned} H_\theta&\simeq\Re\left(\frac{1}{2}R_1(a)\frac{\alpha^3}{\rho}\sin\chi\cos\theta e^{i(\lambda+\rho)}\right)\\ &=\Re\left(\frac{1}{2}R_1(a)\frac{\alpha^3}{\rho}\sin\chi\cos\theta \left[\cos(\lambda+\rho)+i\sin(\lambda+\rho)\right]\right)\\ &=\frac{1}{2}\frac{\omega^2}{c^2}a^3R_1(a)\frac{1}{r}\sin\chi\cos\theta\cos\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right]. \end{aligned}\end{equation*}

The azimuthal component has similar intermediate steps as above. Hence skipping intermediate working gives,

\begin{equation*}\begin{aligned} H_\varphi&\simeq\Re\left(\frac{1}{2}R_1(a)\left[\left(\frac{\alpha^2}{\frac{9i}{\alpha^4}-\frac{3i}{\alpha^3}}\right)\frac{ie^{i\rho}}{\rho}\cos{2\theta}-\frac{\alpha^3}{i}\left(\frac{-ie^{i\rho}}{\rho}+\frac{-e^{i\rho}}{\rho^2}\right)\right]i\sin\chi e^{i\lambda}\right)\\ &=\Re\left(\frac{1}{2}R_1(a)\left[\left(\frac{\alpha^6}{9-3\alpha}\right)\frac{1}{\rho}\cos{2\theta}+\frac{\alpha^3}{\rho}+\frac{\alpha^3}{i\rho^2}\right]i\sin\chi e^{i(\lambda+\rho)}\right). \end{aligned}\end{equation*}

The largest term in the line above is the $\alpha^3/\rho$ term. Removing all smaller magnitude terms yields,

\begin{equation*} \begin{aligned} H_\varphi&\simeq\Re\left(\frac{1}{2}R_1(a)\frac{\alpha^3}{\rho}i\sin\chi e^{i(\lambda+\rho)}\right)\\ &=\frac{1}{2}\frac{\omega^2}{c^2}a^3R_1(a)\frac{1}{r}\sin\chi\sin\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right]. \end{aligned}\end{equation*}

Now consider the E field. Using similar intermediate steps and the same principle as before where the largest order term is the only term considered and the limit forms of $h_1(\rho)$ and $h_1(\alpha)$ are used. The E field’s radial component will be,

\begin{equation*} \begin{aligned} E_r&\simeq\Re\left(\frac{1}{2}\omega\mu_0aR_1(a)3\left(\frac{\alpha}{\alpha\frac{9i}{\alpha^4}-\frac{3i}{\alpha^3}}\right)\frac{ie^{i\rho}}{\rho^2}\sin\chi\sin{2\theta}e^{i\lambda}\right)\\ &=\Re\left(\frac{3}{2}\omega\mu_0aR_1(a)\frac{\alpha^4}{6}\frac{1}{\rho^2}\sin\chi\sin{2\theta}e^{i(\lambda+\rho)}\right) \end{aligned}\end{equation*}

however, due to the $\alpha^4/\rho^2$ term, the radial component is equivalently equal to zero when compared to the other E field components, this implies that $E_r\approx0$ . The polar component is,

\begin{equation*}\begin{aligned} E_\theta&\simeq\Re\left(\frac{1}{2}\omega\mu_0aR_1(a)\left[\frac{\alpha^4}{6i}\left(\frac{-e^{i\rho}+\frac{ie^{i\rho}}{\rho}}{\rho}\right)\cos{2\theta}-\frac{\alpha^2}{i\rho}e^{i\rho}\right]\sin\chi e^{i\lambda}\right)\\ &=\Re\left(\frac{1}{2}\omega\mu_0aR_1(a)\left[\frac{\alpha^4}{6i}\left(\frac{-1}{\rho}+\frac{i}{\rho^2}\right)\cos{2\theta}-\frac{\alpha^2}{i\rho}\right]\sin\chi e^{i(\lambda+\rho)}\right)\\ &\simeq-\frac{1}{2}\frac{\omega^2\mu_0}{c}a^3R_1(a)\frac{1}{r}\sin\chi\sin\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right]. \end{aligned}\end{equation*}

By using some of the previous intermediate steps from above, the azimuthal component is,

\begin{equation*}\begin{aligned} E_\varphi&\simeq\Re\left(\frac{1}{2}\omega\mu_0aR_1(a)\left[\frac{\alpha^4}{6i}\left(\frac{-1}{\rho}+\frac{i}{\rho^2}\right)-\frac{\alpha^2}{i\rho}\right]i\sin\chi\cos\theta e^{i(\lambda+\rho)}\right)\\ &\simeq-\frac{1}{2}\frac{\omega^2\mu_0}{c}a^3R_1(a)\frac{1}{r}\sin\chi\cos\theta\cos\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right]. \end{aligned}\end{equation*}

The results of the radiation field equations are gathered in Appendix B as Equations (B1) and (B2). The results obtained here match the results obtained by Deutsch (1955).

These results are interesting as we can see the retarded time affect, $(\frac{r}{c}-t)$ , on the magnetic field in the radiation zone. That is we would not expect an instantaneous change to the electromagnetic fields when it is far from the surface from the star as electromagnetic ‘news’ travels at the speed of light (Griffiths 2017). Instead we have a delayed affect which depends on the distance from the star.

5.2. Power equation—Deutsch (1955) Equation (15)

The power P passing through a surface is the integral of the Poynting vector ( $\boldsymbol{{S}}=\boldsymbol{{E}}\times \boldsymbol{{H}}$ ) over said surface. From Griffiths (2017), Equation (11.1),

(32) \begin{equation} P(r,t)=\oint{\boldsymbol{{S}}\cdot d\boldsymbol{{a}}}.\end{equation}

This is the power passing through the surface, so a negative must be introduced to give the rate at which energy is radiated away from the system,

(33) \begin{equation} -\frac{dW}{dt}=\oint{\boldsymbol{{S}}\cdot d\boldsymbol{{a}}}.\end{equation}

For the E and H fields Griffiths (2017) advocates taking the fields when r is very large, when they are in radiation form. Therefore, by taking the field Equations (B1) and (B2) (which apply when $r\gg c/\omega$ ) leads to the Poynting vector,

\begin{equation*} \boldsymbol{{S}}=\begin{pmatrix} 0\\ E_{\theta}\\ E_{\varphi} \end{pmatrix} \times \begin{pmatrix} H_r\\ H_{\theta}\\ H_{\varphi} \end{pmatrix} =\begin{pmatrix} E_{\theta}H_{\varphi}-E_{\varphi}H_{\theta}\\ E_{\varphi}H_r\\ E_{\theta}H_r \end{pmatrix}\end{equation*}

and by also taking the surface to be a sphere of radius r, this implies that $S_r\parallel\boldsymbol{{a}}$ , $S_\theta\perp\boldsymbol{{a}}$ and $S_\varphi\perp\boldsymbol{{a}}$ , as $\boldsymbol{{a}}$ is the normal of the surface (see Figure 1a). So (33) becomes,

\begin{equation*} -\frac{dW}{dt}=\int_0^{2\pi}\int_0^\pi(E_{\theta}H_{\varphi}-E_{\varphi}H_{\theta})r^2\sin(\theta)\,d\theta\,d\varphi.\end{equation*}

Upon applying expressions for the H and E fields from (B1) and (B2) we get,

(34) \begin{equation}\begin{aligned} -\frac{dW}{dt}&=\int_0^{2\pi}\int_0^\pi \frac{1}{4}\frac{\omega^4\mu_0}{c^3}a^6 R_1^2 (a)\sin^2(\chi)\left\{\sin^2\left[\omega\left(\frac{r}{c} - t\right) + \varphi\right]\right.\\ &\quad \left. + \cos^2(\theta)\cos^2\left[\omega\left(\frac{r}{c} - t\right) + \varphi\right]\right\}\sin(\theta)\,d\theta\,d\varphi.\end{aligned}\end{equation}

Now we evaluate the integrals over the functions that depend on $\theta$ and $\varphi$ . The double integral can firstly be expressed as,

(35) \begin{gather} \begin{aligned} \int_0^{2\pi}\int_0^\pi &\sin(\theta)\sin^2\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right]\\ &+\sin(\theta)\cos^2(\theta)\cos^2\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right]\,d\theta\,d\varphi, \end{aligned}\end{gather}

and using the integration identities,

\begin{gather} \int_0^{2\pi}\sin(\theta)\,d\theta =2, \nonumber\\ \int_0^{2\pi}\sin(\theta)\cos^2(\theta)\,d\theta=\frac{2}{3}, \nonumber\\ \begin{gathered} \int_0^\pi\sin^2\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right]\,d\varphi\\ =\int_0^\pi\cos^2\left[\omega\left(\frac{r}{c}-t\right)+\varphi\right]\,d\varphi=\pi, \end{gathered} \nonumber\end{gather}

which leads to the integral (35) yielding $8\pi/3$ . On combining this result with (34) we get

(36) \begin{equation} -\frac{dW}{dt}=\frac{2\pi\omega^4\mu_0}{3c^3}a^6R_1^2(a)\sin^2(\chi).\end{equation}

This is the same result obtained by Deutsch (1955) in his Equation (15). The RHS is always positive. Hence, $dW/dt$ must always be negative. That is, energy is always being radiated away from the star. The energy that is radiated comes from the rotational kinetic energy causing the star’s spin down.

5.2.1. Spin down rate—Deutsch (1955) Equation (16)

To derive the spin down rate of a star in the VRD solution, we require two standard relationships (see for example Moebs, Ling, & Sanny 2016, Equations (10.31) and (11.8)),

\begin{equation*} \boldsymbol{\tau}=\frac{d\boldsymbol{{L}}}{dt}, \qquad \frac{dW}{dt}=P=\boldsymbol{\tau}\cdot\boldsymbol{\omega}, \quad \implies \quad \frac{dW}{dt}=\frac{d\boldsymbol{{L}}}{dt}\cdot\boldsymbol{\omega},\end{equation*}

where $\boldsymbol{{L}}$ is the angular momentum and $\boldsymbol{\tau}$ is the torque. As the angular velocity and the angular momentum are in the same direction, we can rearrange the equation above to become,

\begin{equation*} \frac{d\boldsymbol{{L}}}{dt}=\frac{\boldsymbol{\omega}}{\left\lVert{\boldsymbol{\omega}}\right\rVert^2}\frac{dW}{dt},\end{equation*}

which is the result obtained by Deutsch (1955).

However, another typo is present in Deutsch (1955). It is not obviously shown that his $\omega$ is actually a vector and that the $\omega$ in the bottom of the fraction is the norm of $\omega$ . Explaining why he has $\omega/\omega^2$ rather than just $1/\omega$ . Lastly, Deutsch (1955) displays it with a negative so it can be compared easier with (36),

(37) \begin{equation} -\frac{d\boldsymbol{{L}}}{dt}=\frac{\boldsymbol{\omega}}{\left\lVert{\boldsymbol{\omega}}\right\rVert^2}\left(-\frac{dW}{dt}\right).\end{equation}

It can be seen from this equation that the angular momentum also decreases with time. The star will be brought to rest by its own radiation as Deutsch (1955) states.

6. Conclusion

By considering an idealised star and assuming a dipolar form for the internal magnetic field, the general field equations of the electromagnetic fields external of the star were found analytically. These general field equations were then able to be reduced to their radiation limit ( $r\gg c/\omega$ ). From the radiation form, the relationship for the power radiated by the rotating magnetic field was found.

It was found that the power equation matches the equation found in Deutsch (1955) and the derivation for which is shown in this review article for scrutiny. However, a typo appears in Deutsch (1955)’s general field equations. In particular, the $E_\theta$ (31) component appears to evaluate the wrong fraction at $\alpha$ . This error in Deutsch (1955) causes trouble when reducing the general field equations to their radiation form. When using the $E_\theta$ component found in this review article the reduction to radiation form matches Deutsch (1955).

A plot of the magnetic field using the general field equations was also produced, showing the summation of the magnetic dipole and the wave solution. This generalises a similar plot given in Deutsch (1955) which was zoomed in around the $r\ll c/\omega$ region.

A 1969 paper by Goldreich and Julian argues that a pulsar cannot be surrounded by a vacuum but rather a magnetosphere containing plasma (Goldreich & Julian 1969). A related approach can be taken when estimating the power equation. The Force-Free (FF) magnetosphere, which is the solution to Maxwell’s Equations for a rotating star with a dipole field. It ignores the plasma pressure and requires that the electric field parallel to the magnetic field is zero everywhere (Contopoulos, Kazanas, & Fendt 1999). This solution requires approximations just like the VRD solution. The most realistic way is to use numerical simulations. The results from Kalapotharakos et al. (2018) find their power equations ranging from the near FF solution to near the VRD solution given in Equation (36). The VRD solution proves to be an interesting limiting case of numerical simulations.

The VRD solution discussed in this review article has been extended in a number of ways. For example, general relativistic corrections to the Deustch solution can be added (Rezzolla et al. 2001; Pétri 2013). Also, considering an electromagnetic field that consists of more multipole terms than just a dipole adds extra correction terms that are potentially required to match reality (Pétri 2015; Bonazzola, Mottez, & Heyvaerts 2020).