In a simulation, individual particles can be tagged and tracked. This would allow, for example, gathering and plotting particle trajectories.
Warp provides a facility for tracing particles. Add the following to the input file, just after importing warp:
from warp.particles.singleparticle import TraceParticle
The class TraceParticle can be then be used. An instance of the class needs to be created, passing in the starting coordinates and velocities of the particles to be traced. It is recommended that a separate species be created to be used for tracking. This will be more efficient than using one of the simulation species, which will likely have a large number of particles already.
electrons_tracked = Species(type=Electrons)
ntrack = 100
track_radius = 1.*mm
track_vthermal = sqrt(5.*jperev/electrons_tracked.mass)
electron_tracker = TraceParticle(js=electrons_tracked.jslist[0],
x=track_radius*random.random(ntrack),
y=track_radius*random.random(ntrack),
z=0.,
vx=track_vthermal*random.normal(scale=track_vthermal,size=ntrack),
vy=track_vthermal*random.normal(scale=track_vthermal,size=ntrack),
vz=track_vthermal*random.normal(scale=track_vthermal,size=ntrack))
This will randomly fill a square sheet of particles at z=0., with a thermal energy of 5 eV. Note that the weight of the electrons_tracked was not set and will be zero, so they will not contribute to the charge density and self-field calculation. However, they will contribute to any built in moments calculation that includes all species.
There are several additional options that can be set. For example, the input parameter savedata controls how often trajectories are saved, and lsavefields flags sets whether the E and B fields that the particles experience are saved.
During that calculation, the positions and velocities of the tracked particles will be saved. This can be accessed and plotted at any time.
electron_tracker.getz(i=27)
will return the history of the z position of the 27th particle.
electron_tracker.pzx(i=18,color=red)
will plot the z versus x trajectory of the 18th particle, coloring it red.
If tracer particles are lost, for example by hitting a wall, the trajectory up to that point will be saved.
More information can be seen by using "doc(TraceParticle)" or by viewing the source, singleparticle.py.
Automatically generated ID numbers
Particles can be tagged directly, by either giving them a identification number or by giving them a user generated tag. Tagging is necessary to keep track of individual particles since their position in the particle arrays will change during the course of the simulation. ID numbers, sometimes also called "social security numbers", are setup by including the following line in the input file, before generate and before any particles are created.
top.ssnpid = nextpid()
After this, all particles will automatically be given a unique ID. So, for example, early in the calculation, you can do:
idsave = electrons.getssn(zl=0.49,zu=0.51)
xsave = electrons.getx(zl=0.49,zu=0.51)
This gets the ID's and x positions of all of the electrons that have z position within the range 0.49 and 0.51 meters. Then, later in the calculation, you can do:
xnewsave = zeros(len(idsave))
for i in range(len(idsave)):
xnewsave[i] = electrons.getx(ssn=idsave[i])
This saves the new x position of the same particles whose ID was saved earlier. Note that the above code should check if the particles were lost. If the particle was lost, then electrons.getx(ssn=idsave[i]) will return an empty array and the assignment to xnewsave[i] would raise an error.
User generated tags
Note that the above code can be very slow if there are a large number of particles and a large number of particles whose ID's are saved. A faster, though somewhat more complicated method, is to setup a user generated tag. In this case, the particles of interest are given tags. Add the following line to the input file before generate and before any particles are created:
electrons.addpid('flag1')
This creates a flag for each particle. This flag will stay with the particles as they are moved around in memory and among processors in a parallel simulation. The flag is a real number and can have any value (though note that untagged particles will have a value of zero for the flag). The above example, saving the original x and the new x position, can be done using the following code.
ii = electrons.selectparticles(zl=0.49,zu=0.51)
electrons.flag1[ii] = arange(1,len(ii)+1)
xsave = electrons.xp[ii]
The first line gets the current indices of the particles that have their z position in the specified range. The next sets the flag for those particles - each is given a unique number that can be used for sorting later. Note that the unique numbers must be non-zero since that is the default value of the flag for the rest of the particles. The third line saves the x position of those particles. The example uses fancy numpy indexing - please see the numpy documentation for more information.
Then, after the simulation is run, the following will collect the new x positions.
ii = nonzero(electrons.flag1)[0]
id = nint(electrons.flag1)
xnewsave = zeros_like(xsave)
xnewsave[id-1] = electrons.xp[ii]
The first line gets the indices of the particles that had been flagged. nonzero is a numpy function which returns the indices of the elements of the input array which are nonzero. The next line gets the unique number originally given to each particle. In the last line, this unique number is used as to index the array xnewsave, putting the x's in the same order as they are in xsave. Note that the "-1" is needed in "id-1" since the numpy arrays are zero-index based. Also, nint is used to convert flag1, which is floating point, to a usable index, which must be integer.
The most important thing to note here is that the particles at the end of the simulation will almost certainly be in a different order than they were at the beginning. The id, which holds the unique numbers, obtained here afterwards will be in a different order than they were at the beginning. Care is needed to ensure that the particles are obtained in the correct order. This example uses fancy numpy indexing to do this, for example the [id-1] indexing of xnewsave.