Particle Trajectory Calculations

Common or useful tips and techniques for controlling particles or recording data on particle trajectories, often via user programming or other techniques.

Overview

SIMION calculates trajectories of particles during the “Fly’m”. Normally this is done according to the Lorentz Force Law within the electromagnetic fields previously calculated or defined, but you can apply user programming code to alter the equations of motion. See Chapter 8 “Flying Particles” and the “Calculation Methods” appendix of the SIMION User Manual for details.

Here are some common or useful tips and techniques for controlling particles or recording data on particle trajectories, often via user programming or other techniques.

Background

How is time defined

See Time.

Questions

Why are my particles stopping?

The reason why a particle terminates can be determined from Data Recording. In Data Recording, enable recording of “Events” when “Ion’s Splat” and set “Format” to “Verbose”. Splat reasons are listed in the SIMION User Manual Chapter 8 “Flying Particles” “8.6.6 Selecting What to Record”. The ion_splat reserved variable in user programs also contains this information (“L.2.4 SIMION Reserved Variables”).

2D views of a 2D cylidrical symmetry PA only show the Z=0 cross sectional plane of the electrodes. But particles may splat at other Z values. 3D views may show more clearly where particles splat.

User programs can also cause premature splatting by setting the ion_splat reserved variable to a non-zero value.

Count/detection

Counting number of particles

A user program may want to know the number of particles in the simulation, such as for computing a percent ion transmission:

simion.workbench_program()
local nhit = 0
function segment.terminate()
  if ion_px_mm > -5 and ion_py_mm < 5 then  -- acceptance region
    nhit = nhit + 1
  end
end
function segment.terminate_run()
  print('% transmission:', nhit / sim_ions_count * 100)
end

The sim_ions_count reserved variable in SIMION 8.2 obtains the number of particles being flown.

In SIMION 8.0/8.1 (a similarly in 7.0), a workaround is this:

local nions = 0
function segment.initialize()
  nions = ion_number
end

Explanation: The initialize segment gets called for each particle that is created inside a PA instance. Assuming all particles are created inside a PA instance, after all initialize segments are run, nions will contain the last ion number, which also is the number of particles. You will need to ensure particles are created inside PA instances to make this work (or add an empty PA instances there).

If using grouped flying, to count the number of particles still flying (in the previous time step), you might do this:

local n = 0
local last_time    -- be sure to set to nil prior to any rerun
function segment.other_actions()
  last_time = last_time or ion_time_of_flight
  if ion_time_of_flight ~= last_time then
    print("particles flying =", n)
    n = 1
  else
    n = n + 1
  end
  last_time = ion_time_of_flight
end

See also SIMION Example: count_transmission.

Counting number of times a particle crosses a plane

To count the number of times a particle crosses the plane x = 100 mm, you can of course enable Data Recording to record on X crossing and post-process the output file in like Excel.

Alernately, you could do this in a user program with something like this:

simion.workbench_program()
local xcross = 100
local xlast = {}
local count = 0
function segment.other_actions()
  if ((xlast[ion_number] or ion_px_mm) < xcross) == (ion_px_mm >= xcross) then
    count = count + 1
  end
  xlast[ion_number] = ion_px_mm
end

Explanation:

  • xlast stores the x position of the particle in the previous time-step, which can be compared to the x position of the particle in the current time-step (ion_px_mm). You may use a table (keyed by ion_number) for this variable, to store different x positions for each particle, so the code will work even if you have multiple particles flying simultaneously via Group flying.

  • count counts the total number of crossing for all particles, but if you have multiple particles, you might want to make it a table as well.

  • In the first time-step, xlast[ion_number] will be undefined (nil). To prevent a run-time error (nil compared to a number), we use (xlast[ion_number] or ion_px_mm), which means use xlast[ion_number] if it is defined; otherwise, use ion_px_mm. In the first time step, the conditional code will never be executed because ion_px_mm < xcross and ion_px >= xcross will never be the same.

  • The construct (xlast[ion_number] < xcross) == (ion_px_mm >= xcross) means that the two comparisons either both return a true result or both return a false result. In other words, the particle was to the left the X plane in the previous time-step and now is to the right of the X plane; otherwise, the particle was to the right of the X plane in the previous time-step and now is to the left of the X plane. This allows the code to count regardless from which side of the plane the particle crosses from. This is more general than writing (xlast[ion_number] < xcross) and (ion_px_mm >= xcross), which would only would count particles crossing from the left.

Be careful if your particle has small oscillations superimposed on the large circular motion. If so, it may be possible for the particle to cross the plane multiple times per circular cycle. One way perhaps to overcome that is to decrement the counter if the particle crosses the plane in the reverse direction:

local was_left = xlast[ion_number] or ion_px_mm) < xcross
local is_right = ion_px_mm >= xcross
if was_left and is_right then  -- cross from left
  count = count + 1
elseif not was_left and not is_right then  -- cross from right
  count = count - 1
end
-- alternately, test sign of velocity ion_vx_mm to determine direction.

You can also use the test plane library SIMION Example: test_plane, which makes it easier to do harder things, as well as planes not orthogonal to axes or multiple planes:

simion.workbench_program()
local TP = simion.import 'testplanelib.lua'
local count = 0
local test1 = TP(100,0,0, 1,0,0, function()  -- (point and normal vector)
  count = count + 1
end)
function segment.tstep_adjust()
  test1.tstep_adjust()
end
function segment.other_actions()
  test1.other_actions()
end

Counting splats within region or percent transmission

See SIMION Example: count_transmission.

Recording which electrode number a particle hit

You can generally tell which electrode a particle hit by the splat position. The splat position is available via Data Recording or in user programs via the (ion_px_mm, ion_py_mm, ion_pz_mm reserved variables in mm workbench coordinates or (ion_px_gu, ion_py_gu, ion_pz_gu in PA grid units).

For example:

simion.workbench_program()
function segment.other_actions()
  if ion_splat ~= 0 and ion_px_mm > 10.5 and ion_px_mm < 11.5 then
    print 'hit'
  end
end

If the electrode have different voltages, you might also to be able to tell from the voltage at the current particle position (“Potential: V” in Data Recording or ion_volts in user programs).

Below is a more general solution that automatically infers the electrode number by consulting the PA# definition file. It requires SIMION 8.1. Add the PA# file to the workbench in the same position as the PA0 file but with a lower PA priority number (e.g. #1 using L+/L- buttons) so that the PA0 file masks the PA# from view of the particles. Then create a user program like this:

simion.workbench_program()
local padef = simion.wb.instances[1].pa
local function round(x) return math.floor(x + 0.5) end
function segment.terminate()
  local xg = round(ion_px_abs_gu)
  local yg = round(ion_py_abs_gu)
  local zg = round(ion_pz_abs_gu)
  local num = (padef:potential(xg,yg,zg))
  print(ion_number, num, 'at', xg,yg,zg, ion_volts)
end

That uses the current particle position in PA array units. The position is rounded to the nearest grid point so that you are not just outside the electrode in non-electrode space. It looks up the potential at that same point in the .PA# array (here assumed to be PA instance #1 as seen in simion.wb.instances[1].pa). The .PA# array actually contains electrode numbers rather than potentials.

We use padef:potential which returns a non-interpolated potential at a PA grid point. padef:potential_vc is not used since it interpolates, so it may for example return 5.5 V on a 6 V electrode. (post)

Storing data on particles for statistical computation

Say you want to compute statistics on particles such as min, max, or mean splat points. In simple cases like this you can update statistics incrementally as particles splat:

simion.workbench_program()
local ymin = math.huge  -- positive infinity
local ymax = -math.huge -- negative infinity
local ysum = 0
local nparticles
function segment.initialize()
  nparticles = ion_number  -- gets last particle number
end
function segment.terminate()
  ymin = math.min(ymin, ion_py_mm)
  ymax = math.max(ymax, ion_py_mm)
  ysum = ysum + ion_py_mm
end
function segment.terminate_run()  -- SIMION 8.1
  local ymean = ysum / nparticles
  print('y min,max,mean=', ymin, ymax, ymean)
end

Each call to an other_actions or terminate segment only has data for a single particle at a single point in time.

At times it’s useful for those segments to store that data away in a Lua table (perhaps keyed by ion_number) so that a terminate_run segment can subsequently compute statistics on ALL the data previously collected into that table. The above example can be equivalently written:

simion.workbench_program()
local ys = {}  -- table keyed by ion_number
function segment.terminate() -- store data
  ys[ion_number] = ion_py_mm
end
function segment.terminate_run()  -- SIMION 8.1
  -- Compute stats from stored data.
  local ymin = math.huge
  local ymax = -math.huge
  local ysum = 0
  local nparticles = 0
  for i, y in pairs(ys) do
    ymin = math.min(ymin, ion_py_mm)
    ymax = math.max(ymax, ion_py_mm)
    ysum = ysum + ion_py_mm
    nparticles = nparticles + 1
  end
  local ymean = ysum / nparticles
  print('y min,max,mean=', ymin, ymax, ymean)
end

SIMION Example: emittance and SIMION Example: count_transmission do things like that.

Associating data with particles (or storing per-ion data in tables)

For each particle, SIMION maintains data like mass, charge, velocity, etc. You might want to associate your own data on individual particles for Data Recording or statistical computation purposes. For example, you might want to maintain the number of times each particle crosses a plane. The ion color (ion_color) or (if not using charge repulsion) the charge weighting factor CWF (ion_cwf) are normally convenient variables for this since they often don’t affect the trajectory calculation, so you can use them for your own purposes. The CWF is also available in Data Recording.

Often a more flexible option is use your own Lua tables keyed by ion_number. The following example prints the starting and ending X position for each particle:

simion.workbench_program()
local xs = {}
function segment.initialize()
  xs[ion_number] = ion_px_mm
end
function segment.terminate()
  print(ion_number, ion_px_mm, xs[ion_number])
end

In contrast, the following would not work as desired because all the ions overwrite the single x variable:

-- does not work correctly
simion.workbench_program()
local x
function segment.initialize()
  x = ion_px_mm
end
function segment.terminate()
  print(ion_number, ion_px_mm, x)
end

Here is another example:

-- Records particles whose y values ever exceed 10 mm.
simion.workbench_program()
local outside = {}
function segment.other_actions()
  if ion_py_mm > 10 then
    outside[ion_number] = true
  end
end
function segment.terminate_run()
  for i, value in pairs(outside) do
    print(i, value)
  end
end

Another application of these tables is ensure that if particles cross a test plane multiple times that only the first crossing is recorded:

simion.workbench_program()
local TP = simion.import 'testplanelib.lua'
local seen1 = {}  -- table keyed by ion_number
local test1 = TP(31,0,0, 1,0,0, function()
  if not seen1[ion_number] then
    seen1[ion_number] = true
    mark()
    print('In test plane X: n=', ion_number, 'x=', ion_px_mm)
  end
end)
function segment.tstep_adjust()
  test1.tstep_adjust()
end
function segment.other_actions()
  test1.other_actions()
end

Using the table, these also all work with Grouped particle flying too.

Record the initial and final particle positions on the same line of an output file

One approach is to attached a user program to your workbench as follows:

-- SIMION 8.0
simion.workbench_program()
local xs = {}
function segment.initialize()
  xs[ion_number] = ion_px_mm
end
function segment.terminate()
  print(xs[ion_number], ion_time_of_flight)
end

Here, xs is a table that stores data keyed by particle number. It is consulted when outputting data upon particle splats.

An alternate approach is to post-process your output file. I’ll assume that you have Data Recording options set to record to the file “out.txt” and that “ion number” is checked under “What data to record”. Now, in your project folder containing your IOB file, create a file named “postprocess.lua” containing the following code:

-- postprocess.lua, SIMION 8.0

-- open files
collectgarbage()  -- ensure unused files are closed
local fin = assert(io.open("out.txt"))
local fout = assert(io.open("out2.txt", "w"))

-- read and merge lines based on first integer
-- on line (presumably ion number)
local saved = {}
for line in fin:lines() do
  local n = line:match("%d+")
  if n then
    saved[n] = (saved[n] and saved[n].." , " or "") .. line
  end
end

-- output file
for _,v in pairs(saved) do
  fout:write(v, "\n")
end
fin:close()
fout:close()

print "done"

Do a Fly’m to generate the “out.txt” file. Then to generate a text file “out2.txt” in the format you want, enter this in the command bar at the bottom of the SIMION window:

dofile "postpocess.lua"

Measure beam

Measure length of particle trajectory

To calculate the length of a particle trajectory in mm, see SIMION Example: lens_properties - particle_distance_traveled.iob (SIMION 8.1, 2016-12-09).

Measure beam diameter

Beam diameter can of course be measured with the mouse, hovering the mouse over the edges of the beam and reading coordinates in the bottom status bar.

To calculate the diameter of a beam programmatically, there are various approaches. You can calculate min/max values (like examples above) or statistics like standard deviation on splat points or test planes (see SIMION Example: test_plane), either from Data Recording results or user programs.

If you don’t know on which plane the min/max diameter would be but want to also compute that location too, then that is more involved. Using the beam method of charge repulsion (even with repulsion amount of zero) to force space-coherent integration can simplify this calculation since the particles will all be on the same plane on every time-step. Otherwise more generally you’ll need to store the full trajectory path (x,y,z on every time-step for every particle) for later analysis. See SIMION Example: lens_properties beam_minmax_diameter.iob (SIMION 8.1, 2016-12-09) and also confusion_disc.iob.

Measure aberration coefficients

In theory, aberration coefficients can be calculated from field data or from ray tracing. The particle data is available via Data Record but is more conveniently available via user programming, especially if plots of aberration coefficients v.s. system parameters (e.g. lens voltages) are to be created. See SIMION Example: lens_properties (in SIMION 8.0 and simplified significantly in SIMION 8.1) for various ray tracing examples, such as

  • pq.iob - Generates a P-Q curve for a two-cylinder lens. Image position (P) and object position (Q) relative to the reference plane (x=0) are computed and plotted for combinations of lens voltage ratio (V2/V1) and lens linear magnification (M).

  • zoom_lens_curve.iob - Generates a lens focal curve for a three-cylinder “zoom” lens. plots curves of focus locus, as well as linear magnification and spherical and chromatic aberration coefficients. This is an extension of surface_plot.iob.

  • zoom_lens_curve_75.iob - similar to lens_curves.iob but with different lens parameters.

  • zoom_afocal.iob - tunes lens voltages for a five-cylinder afocal zoom lens.

  • confusion_disc.iob - calculates disc of least confusion and plots beam size along axial distance.

Lens also has articles on this.

Measure emittance

See Emittance. Data recording can record data necessary for emittance calculations, such as position and velocity at particle start, splat, or X/Y/Z crossing planes. This can also be computed more seamlessly by a user program. See SIMION Example: emittance for an example calculating beam emittance programmatically.

Measure current

Q: How do I calculate current in the beam?

It depends what you intend the particles you simulate to represent.

Perhaps each particle trajectory you trace shall represent a continuous beam of 2 microamps of current (this current needn’t be actually entered into SIMION). Simply sum the currents of the trajectories to get the total current of the overall beam. That’s most likely what you want. And if you want current density then that’s current divided by area (A/mm2).

Alternately, perhaps each particle you trace shall represent a point charge of some number of coulombs of charge that moves through the system over time (TOF). Current is defined as charge per time. It’s even possible in this simulation that the current could vary over time. One way you might simulate this is to emit some number N of particles each of charge C at different times over some period T. Current is N*C/T. The T here does not necessarily need to be TOF here, though sometimes it’s nice to make T >> TOF if the beam is continuous rather than a short pulse.

Consider the following scenario. 5000 electrons start at the cathode at x=0 at time t=0 with the exact same velocity and they all reach the target on x=10 mm all at the exact same time t=20 microseconds. In this scenario, the current at the target will always be zero, except at time t=20 microseconds, at which time the current (delta charge per delta time where deta time = 0) is effectively infinite. Alternately, say the particles have a slight variance in velocity and they reach the target with a spread of about 0.5 microseconds, in which case the current at the target is about (5000 e) / (0.5 microseconds). The TOF (20 microseconds) is not a factor in the calculation of the magnitude of current here.

Related discussion: Monte Carlo Method. and More Notes on the Meaning of Particle Charge.

Phase plots

Phase plots are typically plots of beam position v.s. velocity at a certain cross-sectional plane. This date is available via Data Recording. See Emittance. SIMION Example: excel and SIMION Example: plot have examples of doing this programatically.

Particle definition

Defining particle starting conditions

Particles can be defined in a particle definition file (ION File, FLY File, FLY2 File). See Particle Initial Conditions. They can also be redefined in an initialize segment like SIMION Example: random (e.g. change variables like mass, charge, and velocity). New particles can even be added during a run: simion.experimental.add_particles().

SIMION Example: field_emission has an example of starting particles automatically along a curved potential contour line.

Reproducible particle trajectories when using randomization

When randomization is utilized in the SIMION simulation–such as randomized particle initial conditions (statistical distributions in the FLY2 File), randomized collisions or various other Monte Carlo Method techniques–simulation runs might not be reproducible. That is, two runs under identical conditions might give slightly different results because each run selects a different sample of the particles to trace from the population of all particles physically in the system.

At times it can be useful to disable randomization to make runs reproducible. This can make it easier to compare the effect of changing some parameter.

You have a few options to disable randomization in order to make runs reproducible:

  • (1) Convert the .FLY2 definition to an .ION definition. To do so, just change “How are particles defined?” from “Grouped (.FLY2)” to “Individually (.ION)” on the Particles Define screen. Any statistical distributions in the FLY2 are converted to a simple list of particles, which will remain identical between Fly’ms.

  • (2) Alternately, in the .FLY2 definitions, use “sequences” rather than “distributions”. Sequences are deterministic (e.g. equally spaced positions) rather than randomly sampled.

  • (3) If your user programming code (or FLY2) uses the random number generator, you may want to “seed” the random number generator at the beginning of each run (e.g. segment.initialize_run() or segment.flym()) to ensure the same random numbers are generated on each run. Beware that there are two random number generators: SIMION’s random number generator (simion.rand() and simion.seed()) and Lua’s random number generator (math.random() and math.randomseed()).

Adjusting for KE for starting potential

The following workbench user program subtracts the voltage at each particle starting point from its initial KE. This type of correction can be useful if particles originate from a cathode surface but you start the particle trajectories a short distance away from the cathode surface, perhaps to avoid field calculation distortions near curved surfaces. It can be helpful to correct the KE for the potential drop over which the particle starting point was displaced.

simion.workbench_program()
local Vref = 0  -- reference voltage
function segment.initialize()
  local volts = simion.wb:epotential(ion_px_mm, ion_py_mm, ion_pz_mm)
         -- note: ion_volts not accessible from initialize segment; why?
  local eV = Vref - volts
  assert(eV >= 0)
  local new_speed = ke_to_speed(eV, ion_mass)
  local old_speed = rect3d_to_polar3d(ion_vx_mm, ion_vy_mm, ion_vz_mm)
  assert(old_speed ~= 0)
  ion_vx_mm = ion_vx_mm * new_speed / old_speed
  ion_vy_mm = ion_vy_mm * new_speed / old_speed
  ion_vz_mm = ion_vz_mm * new_speed / old_speed
end
-- note: could use simion.wb:efield to obtain the initial direction
--  vector instead

Particle control

Stop particles after given time

Attach this user program to your workbench to splats particles after a given length in time:

simion.workbench_program()
adjustable max_time = 10   -- microseconds
function segment.other_actions()
  if ion_time_of_flight >= max_time then
    ion_splat = -1
  end
end

The other_actions segment is called on every time step, and setting ion_splat to a non-zero value causes it to splat.

If you need to terminate precisely at a given time, you can use a tstep_adjust segment too. This is often not essential but does improve the preciseness of the termination time without needing to increase the trajectory quality (T.Qual). Use data recording to record TOF on ion splat to verify this works.

simion.workbench_program()
adjustable max_time = 2    -- max time of flight in microseconds

-- Ensure a time-step boundary exactly coincides with max_time (precise)
function segment.tstep_adjust()
  if ion_time_of_flight < max_time then
    ion_time_step = min(ion_time_step, max_time - ion_time_of_flight)
  end
end

local splat_next = {}  -- flag splat on next time-step
function segment.other_actions()
  -- Splat particle if previously scheduled.
  if splat_next[ion_number] then ion_splat = -1 end

  -- If end of current time-step exceeds max_time,
  -- then finish current time-step but schedule splat on next step.
  if ion_time_of_flight >= max_time then
    splat_next[ion_number] = true
  end
end

Stop particles at a certain location

You can attach a user program like this to your workbench (IOB):

simion.workbench_program()
function segment.other_actions()
  if ion_px_mm >= 50 then  -- or "ion_time_of_flight" for time in microseconds
    ion_splat = 1
  end
end

Note however that it’s likely that a time-step straddles the position x = 50 mm, in which case your particle will stop slightly after 50 mm (which might not be ok with you). A simple way to reduce that small skip is to just increase your trajectory quality (TQual) factor. Or you could write a tstep_adjust segment to ensure that time steps end exactly at the desired end positions. More simply, you can use the test plane library (see SIMION Example: test_plane) available in 8.1 (or 8.0 Check for Updates), which handles this all for you.

Yet another solution, which involves no user programming and avoids the time-step inaccuracies, is to add a “stopper electrode” to your workbench. To do so, create a small PA file that contains only electrode points (no empty space) and add it to your workbench at the location where you want to splat ions. The PA must have the highest priority number in the PA Instances list on the PAs tab on the View screen; otherwise your particles might never see it. This stopper electrode, thankfully, will not affect your fields in the rest of your system because fields between PA instances do not interact in SIMION, so your particles will (correctly) not observe any effects from the stopper electrode until they actually hit it.

Another possibility is to refine your PA and then use the Modify Find/Replace functions to convert some non-electrode points to electrode points. However, be careful never to refine the array again, even if SIMION prompts you, because this electrode will affect the fields around it if you re-refine. Therefore, this approach may be error prone. Setting simion.pas pa.refinable to false can suppress the prompts to refine though.

Get/set kinetic energy of current particle

In SIMION 8.2 you can obtain or change the current particle’s kinetic energy in a user program via the ion_ke reserved variables. Internally these read/write to the velocity (ion_vx_mm, ion_vy_mm, ion_vz_mm).

In SIMION 8.0/8.1 (a similarly in 7.0), you can use these functions:

local function get_ke()
  local speed = math.sqrt(ion_vx_mm^2 + ion_vy_mm^2 + ion_vz_mm^2)
  return speed_to_ke(speed, ion_mass)
end

local function set_ke(ke)
  local speed = math.sqrt(ion_vx_mm^2 + ion_vy_mm^2 + ion_vz_mm^2)
  assert(speed ~= 0)
  local new_speed = ke_to_speed(ke, ion_mass)
  local scale = new_speed/speed
  ion_vx_mm = ion_vx_mm * scale
  ion_vy_mm = ion_vy_mm * scale
  ion_vz_mm = ion_vz_mm * scale
end

Flying particles in the reverse direction

As mentioned in Chapter 5 of the SIMION User Manual: “Try to back-fly your ions from their termination points (using reversed directions and termination energies). The ability to duplicate ion trajectories bi-directionally (at least approximately) means that the modeling is reasonably self-consistent. However, neither test guarantees that the results are really correct.”

This backflying can be done with Lua code. The ion color can be used to give these different colors and avoid backflying the backflown particles.

function segment.other_actions()
  if ion_splat ~= 0 and ion_color == 1 then
    ion_splat = 0
    ion_vx_mm = - ion_vx_mm
    ion_vy_mm = - ion_vy_mm
    ion_vz_mm = - ion_vz_mm
    ion_color = 3
  end
end

Time step control

See Time Step Size on adjustment of time step size (trajectory quality, T.Qual) in the trajectory calculation.

Data Record periodically, with constant size time steps, possibly for FFT

Sometimes you need particle data at periodic intervals. This is especially critical for Fast Fourier Transform (FFT) analysis.

The easiest way is to enable Time Markers (on the View screen Particles tab) for a given period in microseconds and then in Data Recording set When to Record to All Markers. Time markers force time steps to end precisely on multiples of the given period and can also be used as a Data Recording event.

A similar thing can alternately be done with user programming by executing print statements at the appropriate times (ion_time_of_flight) in an other_actions segment or invoking mark() within other_actions to trigger Data Recording when the Data Recording All Markers event is selected.

If using the programming method, care is required to ensure that the data output occurs on the precise times. other_actions is only executed at the end of time-steps, so you want to ensure that time-steps end at the precise times in which you want to print data, such as by using a tstep_adjust segment (or very large T.Qual values). Note, however, that numerical error (e.g. 0.9999999999999 does not equal 1) may require introducing a small tolerance factor in time calculation in a tstep_adjust segment.

The following code ensures that time step boundaries fall on multiples of time step_size. For example, 0.001, 0.002, 0.0003, etc. microseconds. SIMION may still insert additional time-steps unless you set the trajectory quality (TQual) to zero:

-- SIMION 8.0/8.1
simion.workbench_program()
adjustable step_size = 0.001  -- microseconds
function segment.tstep_adjust()
   ion_time_step = math.min(ion_time_step, step_size)
end

Custom physics

Custom forces

Normally SIMION applies forces on particles according to the Lorentz Force Law, but you can customize this as described on Lorentz Force Law.

Randomly convert ion to electron

To convert one particle to another particle (e.g. ion is replaced with an electron), just change the charge and mass (ion_charge and ion_mass variables). You probably also want to change the velocity (ion_vx_mm, ion_vy_mm, ion_vz_mm) or KE (Get/set kinetic energy of current particle) too.

simion.workbench_program()

adjustable _mean_free_path_mm = 0.1

function segment.other_actions()
  if ion_px_mm > 10 and ion_px_mm < 20 then
    local speed_ion = sqrt(ion_vx_mm^2 + ion_vy_mm^2 + ion_vz_mm^2)
    local collision_prob = 1 - exp(- speed_ion * ion_time_step / effective_mean_free_path_mm)
    -- Test for collision
    if rand() < collision_prob then
      ion_charge = -1
      ion_mass = 0.00054857990946
      -- TODO: also change ion_vx_mm, ion_vy_mm, ion_vz_mm appropriately.
    end
  end
end

The above example performs conversions at random times with a given region. Here, collision probability is derived from a mean-free-path, similar to that done in SIMION Example: collsion_hs1 (also Collision Model HS1).

You can kind-of simulate this type of thing without user programs too, if you just do a second Fly’m with the child particles. Note that Data Recording has a “.ION” option which saves particles to its output file in ION format that can be read back (optionally after manual editing) into the Particles screen for a subsequent simulation.

Fragment particles

What if one particle fragments into multiple particles or a particle hitting a surface causes multiple secondary emission particles? Techniques other than the above are required.

In SIMION 8.0, the simplest method is to fly two identical AB++ particles and then at a certain point in time (other_actions segment in a user program), change one particle to A+ and the other to B+. You should change not just the mass and charge but also the velocities (momentum) of the particles. For example:

simion.workbench_program()
function segment.other_actions()
  if ion_px_mm > 10 and ion_mass == 100 then
    if ion_number == 1 then
      ion_mass = 99
      ion_charge = 1
      ion_vx_mm = -0.1
      ion_vy_mm = 0
      ion_vz_mm = 0
    else
      ion_mass = 1
      ion_charge = 1
      ion_vx_mm = 10
      ion_vy_mm = 0
      ion_vz_mm = 0
    end
  end
end

SIMION also has an SIMION Example: child_particles though is a little more complicated.

More recently, SIMION 8.2 also has a newer API to create particles during the Fly’m, which is a more flexible way of doing these things: simion.experimental.add_particles().

Secondary emission

Seconary emissions effects is more involved topic with multiple approaches possible for simulating it. Specifically, the surface normal might not be obtained from the jagged PA surface but rather via the nearby (smooth) potential contour, a GEM file defintion, STL file definition, perhaps even surface enhancement info, or just a simple analytical equation provided. There are also various statistics that could applied to secondary energies and angles. See SIMION Example: secondary and SIMION Example: child_particles for some starting points.

Collisional effects

For things like ion-gas collision or other randomized collisions, see Ion-Gas Collisions.

Multiple runs or flym segment

Performing multiple runs is often used for things like Optimization. Examples with multiple runs include SIMION Example: tune (e.g. tune_series.iob/lua). In SIMION 8.0, this was done via the sim_rerun_flym variable. In SIMION 8.1, there is a new segment.flym() segment that makes this nicer.

Here’s a few additional notes on the flym segment. When click the Fly'm button the segment.flym executes. If you didn’t define a segment.flym, the default segment.flym is essentially as follows, which does a does a run:

function segment.flym() run() end

But you can have any code you want inside the segment.flym. This prints “hello” without doing any runs:

function segment.flym() print "hello" end

This prints “hello” and does two runs:

function segment.flym()
  print "hello"
  run() run()
end

Trajectory Display

Clearing trajectories

The View and Retain checkboxes on the Particles tab control whether calculated particle trajectories are displayed on the screen and/or retained on disk. On redrawing the screen, particles will be cleared on the screen if they are not retained on disk. View and retain options can also be controlled programmatically via the sim_trajectory_image_control reserved variable. Redrawing occurs on window resizing or by issuing a simion.redraw_screen() command. Having Rerun enabled (or sim_rerun_flym reserved variable) also effectively disables retaining, although doing multiple runs instead via the segment.flym() segment does not cause that.

Displaying only some trajectories

This will cause only every 10th particle trajectory to be displayed.

simion.workbench_program()
function segment.other_actions()
  if ion_number % 10 == 0 then
    sim_trajectory_image_control = 0
  else
    sim_trajectory_image_control = 3
  end
end
function segment.terminate_run()
  sim_trajectory_image_control = 1  -- keep viewing
end

It uses the sim_trajectory_image_control reserved variable to suppress viewing and retaining of most trajectories. This technique may be useful if you are computing statistics on a lot of trajectories but only care to display some of them. It may be equivalently written using the modulus % operator:

simion.workbench_program()
function segment.other_actions()
  sim_trajectory_image_control = (ion_number % 10 == 0) and 0 or 3
end
function segment.terminate_run()
  sim_trajectory_image_control = 1
end

This example preserves only every 10th time-step (for all particles):

simion.workbench_program()
local n = 0
function segment.other_actions()
  n = (n + 1) % 10
  sim_trajectory_image_control = (n == 0) and 0 or 3
end
function segment.terminate_run()
  sim_trajectory_image_control = 1
end

Another option, without programs, is to just do another Fly’m with a different and smaller set of particle definitions intended only for visualization, not statistics.

Marker Color

If you set the marker color to 15 (transparent/white), SIMION will use the particle color for the marker color.

Particle trails

This interactively shows particle traces only for the previous 10 time-steps. The effect is similar to Dots most flying but the particle has a short tail (like a tadpole).

simion.workbench_program()
function segment.initialize_run()
  sim_trajectory_image_control = 1
end
local n = 0
function segment.other_actions()
  n = n + 1
  if n % 10 == 0 then
    simion.sleep(1.0)
    simion.redraw_screen()
  end
end

Explanation: The sim_trajectory_image_control = 1 turns off trajectory retaining so that simion.redraw_screen() results in previous particle trajectories disappearing. This occurs every 10 time-steps (using the modulus operator % to compute the remainder of the counter divided by 10). simion.sleep(0.5) adds a 0.5 second delay, which can be adjusted.

Fields seen by particles

Overlapping PA instances

A particle in SIMION might not see the fields from all potential array instances in the workbench. A particle only sees the electric fields from the electric potential array instance (if any) with the highest priority number among the electric PA instances containing the current particle position. Likewise, a particle only sees the magnetic fields from the magnetic potential array instance (if any) with the highest priority number among the magnetic PA instances containing the current particle position. The priority numbers are displayed in the PA Instances list on the PAs tab on the View screen. They can be changed with the L- and L+ buttons.

Sometimes, you might not want the fields in the highest priority PA instance to be effective within the entire volume of that PA instance. Perhaps you broke a system into to PA instances that overlap and the fields are not valid on the ends of both PA instances (see SIMION Example: bender_cut). One solution is to crop that PA instance (using the Modify screen Crop function after refining the PA, as shown in SIMION Example: bender_cut). An alternative and simpler solution is to just utilize a user program to apply fields from a different PA instance, like in the example below:

-- Apply E field from different PA instance (requires SIMION 8.2EA-20170127)
simion.workbench_program()
local pa1 = simion.wb.instances[1]
local pa2 = simion.wb.instances[2]
function segment.efield_adjust()
  local painst = (ion_px_mm < 100) and pa1 or pa2
  local ex,ey,ez = painst:field(ion_px_mm, ion_py_mm, ion_pz_mm, ion_time_of_flight)
  ion_dvoltsx_gu, ion_dvoltsy_gu, ion_dvoltsz_gu
       = painst:wb_to_pa_orient(-ex,-ey,-ez)
end

8.1 versions do not support time-dependent fields in painst:field (using ion_time_of_flight) but otherwise can do this if fast_adjust/efield_adjust segments are not used:

-- Apply E field from different PA instance (requires SIMION 8.1.2.30 or above)
simion.workbench_program()
local pa1 = simion.wb.instances[1]
local pa2 = simion.wb.instances[2]
function segment.efield_adjust()
  local painst = (ion_px_mm < 100) and pa1 or pa2
  local ex,ey,ez = painst:field(ion_px_mm, ion_py_mm, ion_pz_mm)
  ion_dvoltsx_gu, ion_dvoltsy_gu, ion_dvoltsz_gu
       = painst:wb_to_pa_orient(-ex,-ey,-ez)
end

The following code superimposes fields from two magnetic fields:

-- Superimposing two B fields (requires SIMION 8.2EA-20170127)
simion.workbench_program()
local pa1 = simion.wb.instances[1]
local pa2 = simion.wb.instances[2]
function segment.mfield_adjust()
  local bx1,by1,bz1 = pa1:bfield(ion_px_mm, ion_py_mm, ion_pz_mm, ion_time_of_flight)
  local bx2,by2,bz2 = pa2:bfield(ion_px_mm, ion_py_mm, ion_pz_mm, ion_time_of_flight)
  local bx = (bx1 or 0) + (bx2 or 0)
  local by = (by1 or 0) + (by2 or 0)
  local bz = (bz1 or 0) + (bz2 or 0)
  ion_bfieldx_gu, ion_bfieldy_gu, ion_bfieldz_gu = wb_to_pa_orient(bx,by,bz)
end
-- B field in gauss at current particle position in workbench
-- coordinates (mm) from both PA instances.  Sums them.
-- Rotates vector back to current PA instance coordinates and sets it.

8.1 versions do not support time-dependent fields in painst:field (using ion_time_of_flight), but the field can accept electrode voltages as parameters:

simion.workbench_program()
local pa1 = simion.wb.instances[1]
local pa2 = simion.wb.instances[2]
function segment.efield_adjust()
  local ex,ey,ez
  if ion_px_mm < 100 then
    ex,ey,ez = pa1:field(ion_px_mm, ion_py_mm, ion_pz_mm, {
      [1] = 200, [3] = 400
     })
     ex,ey,ez = pa1:wb_to_pa_orient(-ex,-ey,-ez)
  else
    ex,ey,ez = pa2:field(ion_px_mm, ion_py_mm, ion_pz_mm, {
      [1] = -100, [2] = 300 * sin(ion_time_of_flight)
    })
    ex,ey,ez = pa2:wb_to_pa_orient(-ex,-ey,-ez)
  end
  ion_dvoltsx_gu, ion_dvoltsy_gu, ion_dvoltsz_gu = ex,ey,ez
end

Prior to 8.1.2.30, the field function was named field_wc.

Output logging

Currently, everything written to the Log window (either via Data Recording GUI or printing in Lua) also gets written to the data recording output file (specified on the Data Recording screen).

Now, you could open your own file handle in Lua and write to it, and that would only write to that file (not the Log window):

simion.workbench_program()
local fh = assert(io.open("output.txt", "w")
function segment.other_actions()
  fh:write(ion_number, " ", ion_px_mm, "\n")
end
function segment.terminate_run()
  fh:close()
end

This useful function ensures all the arguments are strings. It concatenates them with a tab and appends a new line. So, it works quite like print but accepts as a first argument a file handle:

local function fprint(fh, ...)
  local t = {n=select('#', ...), ...}
  for i=1,t.n do t[i] = tostring(t[i]) end
  local s = table.concat(t, '\t')
  fh:write(s, "\n")
end

Example:

fprint(fh, 123.456, "test", 789, math.pi, true)

Be sure to close the file (fh:close()) or at least flush (fh:flush()) when done writing.

Something that might be useful is the simion.status() function, which prints a message to the status bar at the bottom of the SIMION window:

simion.workbench_program()
function segment.other_actions()
  simion.status("The time is " .. ion_time_of_flight)
end

Another example displaying only periodically:

simion.workbench_program()
local nions = 0
function segment.initialize() nions = ion_number end
local last_record_time = 0
local record_period = 0.5  -- microseconds
function segment.other_actions()
  if math.abs(ion_time_of_flight - last_record_time) > record_period then
    last_record_time = ion_time_of_flight
    simion.status(("%0.2g usec, ion=%d/%d"):format(ion_time_of_flight, ion_number, nions))
  end
end

You can also use the GUI Library (simion.experimental.dialog) library to display data in another window. See SIMION Example: gui (e.g. counter.lua in 2016-12-09 or above).

You can also write to Excel or gnuplot. See SIMION Example: excel, SIMION Example: gnuplot and SIMION Example: plot.

See Also