Interactive online version: .

[1]:
%matplotlib inline

Particle Tracker

An example of PlasmaPy’s particle tracker class.

Take a look at the docs to gyrofrequency() and ParticleTracker for more information

Initialize a :class:~plasmapy.plasma.grids.CartesianGrid object. This will be the source of electric and magnetic fields for our particles to move in.

[3]:
grid_length = 10
grid = CartesianGrid(-1 * u.m, 1 * u.m, num=grid_length)

Initialize the fields. We’ll take B in the x direction and E in the y direction, which gets us an E cross B drift in the negative z direction.

[4]:
Bx_fill = 4 * u.T
Bx = np.full(grid.shape, Bx_fill.value) * u.T

Ey_fill = 2 * u.V / u.m
Ey = np.full(grid.shape, Ey_fill.value) * u.V / u.m

ExB_drift(np.asarray([0, Ey_fill.value, 0]) * u.V / u.m, np.asarray([Bx_fill.value, 0, 0]) * u.T)
[4]:
$[0,~0,~-0.5] \; \mathrm{\frac{m}{s}}$

|ParticleTracker| takes arrays of particle positions and velocities of the shape [nparticles, 3], so these arrays represent one particle starting at the origin.

[5]:
x0 = [[0, 0, 0]] * u.m
v0 = [[1, 0, 0]] * u.m / u.s
particle = Particle("p+")

Initialize our stop condition and save routine. We can determine a relevant duration for the experiment by calculating the gyroperiod for the particle.

[6]:
particle_gyroperiod = 1 / gyrofrequency(Bx_fill, particle).to(u.Hz, equivalencies=u.dimensionless_angles())

simulation_duration = 100 * particle_gyroperiod
save_interval = particle_gyroperiod / 10

termination_condition = TimeElapsedTerminationCondition(simulation_duration)
save_routine = IntervalSaveRoutine(save_interval)

Initialize the trajectory calculation.

[7]:
simulation = ParticleTracker(grid, save_routine=save_routine, termination_condition=termination_condition, verbose=False)
/home/docs/checkouts/readthedocs.org/user_builds/plasmapy/envs/latest/lib/python3.12/site-packages/plasmapy/simulation/particle_tracker.py:633: RuntimeWarning: Quantities should go to zero at edges of grid to avoid non-physical effects, but a value of 4.00E+00 T was found on the edge of the B_x array. Consider applying a envelope function to force the quantities at the edge to go to zero.
warnings.warn(
/home/docs/checkouts/readthedocs.org/user_builds/plasmapy/envs/latest/lib/python3.12/site-packages/plasmapy/simulation/particle_tracker.py:633: RuntimeWarning: Quantities should go to zero at edges of grid to avoid non-physical effects, but a value of 2.00E+00 V / m was found on the edge of the E_y array. Consider applying a envelope function to force the quantities at the edge to go to zero.
warnings.warn(

We still have to initialize the particle’s velocity. We’ll limit ourselves to one in the x direction, parallel to the magnetic field B - that way, it won’t turn in the z direction.

Run the simulation.

[9]:

We can take a look at the trajectory in the z direction:

[10]:
t, x, v = save_routine.results()
particle_trajectory = x[:, 0]
particle_position_z = particle_trajectory[:, 2]

plt.scatter(t, particle_position_z)
[10]:
<matplotlib.collections.PathCollection at 0x7f29d43c7590>

or plot the shape of the trajectory in 3D:

[11]:
fig = plt.figure()

ax.plot(*particle_trajectory.T)
ax.set_xlabel("X")
ax.set_ylabel("Y")
ax.set_zlabel("Z")
[11]:
Text(0.5, 0, 'Z')

As a test, we calculate the mean velocity in the z direction from the velocity:

[12]:
v_mean = v[:, :, 2].mean()
print(
f"The calculated drift velocity is {v_mean:.4f} to compare with the "
f"expected E0/B0 = {-(Ey_fill/Bx_fill).value:.4f}"
)
The calculated drift velocity is -0.5076 m / s to compare with the expected E0/B0 = -0.5000