RDYNLAB

This is a library for modeling Relativistic DYNamics in the LAB frame.

The General Problem of Relativistic Particle Dynamics

The general problem of modeling relativistic particle trajectories for a given force boils down to solving the following system of six differential equations:

[math]\frac{d\vec{x}}{dt}=\vec{v}[/math]

[math]\frac{d\vec{v}}{dt}=\vec{a}[/math]

where [math]\vec{x}[/math] is the position of the particle, [math]\vec{v}[/math] is the momentum of the particle, and [math]\vec{a}[/math] is the acceleration of the particle which can be derived from the force. Classically, we would just use [math]\vec{a}=\vec{F}/m[/math]. Such a description, however, is not consistent with special relativity. Indeed, starting from the relativistic relations

[math]\frac{\vec{v}}{c}=\frac{\vec{p}c}{E}[/math] [math]E^2-c^2\vec{p}\cdot\vec{p}=m_0^2c^4[/math],

and taking the time derivatives of both, we find

[math] \frac{1}{c}\frac{d\vec{v}}{dt}=\frac{c}{E}\frac{d\vec{p}}{dt}-c\vec{p}\frac{1}{E^2}\frac{dE}{dt}[/math]

[math] 2E\frac{dE}{dt}-2c^2\vec{p}\cdot\frac{d\vec{p}}{dt}=0[/math].

Solving this second equation for [math]dE/dt[/math] and plugging into the first, also noting that [math]d\vec{p}/dt=\vec{F}[/math], we obtain

[math] \frac{1}{c}\frac{d\vec{v}}{dt}=\frac{c}{E}\vec{F}-\frac{c}{E^3}\vec{p}\left(c^2\vec{p}\cdot\vec{F}\right)[/math].

Noting that [math]c\vec{p}=\frac{E}{c}\vec{v}[/math], we find

[math] \frac{E}{c^2}\frac{d\vec{v}}{dt}=\vec{F}-\frac{\vec{v}}{c}\left(\frac{\vec{v}}{c}\cdot\vec{F}\right)[/math].

We now decompose the force into pieces that are parallel and perpendicular to the velocity of the particle, i.e.

[math] \vec{F}=\vec{F}_{\parallel} + \vec{F}_{\perp} [/math]

where

[math] \vec{F}_{\parallel}=\left(\vec{F}\cdot\hat{v}\right)\hat{v} [/math].

This now gives us

[math] \frac{E}{c^2}\frac{d\vec{v}}{dt}=\left(\vec{F}_{\parallel}+\vec{F}_{\perp}\right)-\frac{\vec{v}}{c}\left[\frac{\vec{v}}{c}\cdot\left(\vec{F}_{\parallel}+\vec{F}_{\perp}\right)\right][/math].

The scalar product of the velocity and the perpendicular component of the force is zero. Furthermore the scalar product involving the velocity and the parallel component of the force, when multiplied again by the velocity, gives us [math]v^2/c^2\vec{F}_{\parallel}[/math], and thus

[math] \frac{E}{c^2}\frac{d\vec{v}}{dt}=\left(\vec{F}_{\parallel}+\vec{F}_{\perp}\right)-\frac{v^2}{c^2}\vec{F}_{\parallel}[/math].

Combining the terms with respect to the force vectors, we get

[math] \frac{E}{c^2}\frac{d\vec{v}}{dt}=\vec{F}_{\perp}+\left(1-\frac{v^2}{c^2}\right)\vec{F}_{\parallel}=\vec{F}_{\perp}+\frac{1}{\gamma^2}\vec{F}_{\parallel}[/math].

Last of all, noting that [math]E/c^2=\gamma m_0[/math]

where [math]m_0[/math] is the rest energy of the particle, and [math]\gamma[/math] is the usual relativistic factor, we obtain

[math] \frac{d\vec{v}}{dt}=\frac{1}{\gamma m_0}\vec{F}_{\perp}+\frac{1}{\gamma^3 m_0}\vec{F}_{\parallel}[/math].

This completes (almost) our set of differential equations.

Power Loss for Accelerated Charges

The code accounts for power radiated away as the charged particles accelerate. Such power loss is given by the Lienard formula,

[math] \frac{dE}{dt}=-\frac{\mu_0q^2\gamma^6}{6\pi c}\left(\vec{a}\cdot\vec{a}-\left|\frac{\vec{v}\times\vec{a}}{c}\right|^2\right)[/math]

However, noting that [math]E=\gamma m_0 c^2[/math], this can also be written as

[math] \frac{d\gamma}{dt}=\frac{1}{m_0 c^2}\frac{dE}{dt}=-\frac{\mu_0q^2\gamma^6}{6\pi m_0c^3}\left(\vec{a}\cdot\vec{a}-\left|\frac{\vec{v}\times\vec{a}}{c}\right|^2\right)[/math].

One can determine the acceleration of the particle using the last formula in the previous section to calculate this quantity.

In order to implement radiative power loss into the dynamical equations, we need to go back to the relationship between force and momentum, i.e.

[math]\vec{F}=\frac{d\vec{p}}{dt}=\frac{d\left(\gamma m_0 \vec{v}\right)}{dt}=m_0\frac{d\gamma}{dt}\vec{v}+\gamma m_0\frac{d\vec{v}}{dt}[/math].

The power losses come through the term involving [math]\frac{d\gamma}{dt}[/math]. Including this term into our previous differential equations yields

[math] \frac{d\vec{v}}{dt}=\frac{1}{\gamma m_0}\vec{F}_{\perp}+\frac{1}{\gamma^3 m_0}\vec{F}_{\parallel}-\frac{1}{\gamma}\frac{d\gamma}{dt}\vec{v}[/math].

Summary of Relevant Formulae

In summary, the problem of relativistic particle dynamics reduces to the solving of the following system of differential equations:

[math]\frac{d\vec{x}}{dt}=\vec{v}[/math]

[math] \frac{d\vec{v}}{dt}=\frac{1}{\gamma m_0}\vec{F}_{\perp}+\frac{1}{\gamma^3 m_0}\vec{F}_{\parallel}-\frac{1}{\gamma}\frac{d\gamma}{dt}\vec{v}[/math]

[math] \frac{d\gamma}{dt}=\frac{1}{m_0 c^2}\frac{dE}{dt}=-\frac{\mu_0q^2\gamma^6}{6\pi m_0c^3}\left(\vec{a}\cdot\vec{a}-\left|\frac{\vec{v}\times\vec{a}}{c}\right|^2\right)[/math],

where the latter equation is evaluated using the acceleration without the power loss term.

Numerical Method of Solving the Differential Equations

These equations are solved subject to an appropriate set of boundary conditions, usually the initial position and momentum of the particle.

In simple cases, such as when a charged particle travels through uniform electromagnetic fields, the set of dynamical equations can be solved analytically. In more complex cases, the system must be solved numerically.

The code solves the system of equations using an embedded 8th order Runge-Kutta Prince-Dormand method with 9th order error estimates with adaptive stepsize control, as implemented in the GNU Scientific Library [1].

I am unable to post code to the wiki, so I have copied and pasted it all here. You will need all of the files in order to compile the code. GSL must be installed before you can compile this library.

The code can be compiled into a shared object library using the autoconf, automake, and libtool tools. Optionally, you could just compile these source files into other programs.

Return to Simulating Particle Trajectories