Link Search Menu Expand Document

There is an automatic object relational mapping that exists for Pic Simulations here: https://pcrumley.github.io/tristanUtils/tristanSim.html

Here is the documentation repeated here:

TristanSim is an object relational mapping (ORM) that exposes an API to access tristan-MP simulations. All of the data attributes are lazily evaluated. The first time you access it the attributes are loaded from the disk, but the next time it is cached.

Code structure & implementation

Basic Tristan Outputs

The TristanSim object takes a path to a directory upon initialization. Then it looks in that directory for any file matching names flds.tot.*, prtl.tot.*, param.* and spect.*. If all 4 files are present it creates and output point for this timestep. Here’s a quick example:

from tristanSim import TristanSim
myRun = TristanSim('/path/to/tristan/output/')

The simulation object has been built. All output points are in the simulation object. For accessing it, you can think of it as a simple list

Since you can treat the simulation object as if it just a list, we can access any of the output points using simple list operators. e.g., len(myRun) gives the number of output files, access the first one by myRun[0] or you can iterate over them.

for out in myRun:
    # Do something for each output point here

To access the data on the disk you have to look at one of the objects in an output point. For instance

myRun[4]
myRun[4]._flds

Fields is a pointer to the 5th flds.tot file in your output dir. You can see the attributes saved to the disk here by

print(myRun[4]._flds.keys())
# you can get any of those attributes by typing e.g.,
myRun[4]._flds.ex
# or more simply
myRun[4].ex

myRun[4].ex is lazily evaluated. The first time it is read from the hdf5 file, but afterwards it is in memory. If you want to delete it and reload myRun[4].reload() Then you can access it from the disk by typing myRun[4].ex. If you want to force all the field output to load and cache type myRun.loadAllFields() or myRun.loadAllPrtls() to load all the prtl outputs. In practice you shouldn’t have to worry too much about the memory access if you are ok with it being lazily evaluated.

myRun[4]._prtl, myRun[4]._spect and myRun[4]._param are similarly defined, but if there are no collisions, any attribute should be able to be accessed as myRun[4].c_omp or whatever. There are a few collisions, and then you have to which one you want to default to. For instance dens is in spect.* and flds.tot. we default to the flds value. This can be set in the __init__() funcion of TristanSim, in the self._colisionFixer dictionary. If you want to add additional hdf5 output files to look for, you can do so in the __init__() function of the TristanSim class.

Fancy Examples

If you have a suite of runs you had run with automater.py, this class comes in handy.

Building an iterable list of all your runs in a directory:

import matplotlib.pyplot as plt
import numpy as np
from tristanSim import TristanSim

# Point to the directory where the suite of runs
# In this example I have 9 different simulations
# saved in ../batchTristan
outdir = '../batchTristan'

# Let's create a list that will hold all our simulation
# instances.
runs = []

# We'll also name each run based on the directory
# it resides.
runNames = []

for elm in os.listdir(outdir):
    elm = os.path.join(outdir, elm)
    if os.path.isdir(elm):
        elm = os.path.join(elm,'output')
        if os.path.exists(elm):
            # get the directory before 'output/'
            dirName = os.path.split(elm)[0]
            dirName = os.path.split(dirName)[-1]
            runs.append(TristanSim(elm))
            runNames.append(dirName)

# Now, any run can be accessed
# e.g. runs[4], with a name runNames[4]
# to access the 3rd output time of run4: e.g. runs[4][3]

# NOTE: because runs is a 1D list of objects, treating
# it as a 2D list, e.g. runs[4,3] WILL NOT WORK.

Plotting ex of the 5th output timestep of each simulation in a 3 x 3 grid. Because of how I set up the automater, all output[i] are at the same physical time.

fig = plt.figure()
axes = fig.subplots(3,3).flatten()
#print(axes)
j = 0
for run, name in zip(runs, runNames):
    ax = axes[j]
    istep = run[5].istep
    comp = run[5].c_omp
    ex = run[5].ex[0,:,:]
    imSize = [0, ex.shape[1], 0, ex.shape[0]]
    ax.imshow(ex,
        extent=(x*istep/comp for x in imSize),
        origin = 'lower')
    ax.set_title(name)
    #plt.colorbar()
    j += 1
    if j == len(axes):
        break
#plt.savefig('test.png')
plt.show()

Let’s plot all the total electron energy as function of time for each run, where the colors, line-styles, and marker styles depend on c_omp, ppc and ntimes.

# First get all of the unique values of c_omp, ppc and
# ntimes from our suite of runs.

c_omp_val = list(set([r[0].c_omp for r in runs]))
ppc_val = list(set([r[0].ppc0 for r in runs]))
ntimes_val = list(set([r[0].ntimes for r in runs]))

# Lists that store what the linestyles will be.
ms = ['.', 'x', '4', '8']
ls = ['-', '--', ':', '-.']
color = ['b', 'r', 'g', 'y']

fig = plt.figure()
for run in runs:
    # In this example, we have fast moving test
    # particles that have negative indices we don't
    # want to count towards this energy.
    plt.plot([o.time for o in run],
     [np.average(o.gammae[o.inde>0]-1) for o in run],
     c = color[ppc_val.index(run[0].ppc0)],
     linestyle = ls[ntimes_val.index(run[0].ntimes)],
     marker = ms[c_omp_val.index(run[0].comp)],
     markersize = 10
    )
plt.show()

Tracking Prtls.

The tristanSim class can also track particles. It requires a particular output.F90 of tristan, but if the files are saved properly, you’ll find all of the tracked particles in a trackedLecs and trackedIon object. The first call builds the database which may take awhile. It is saved afterwards.

Examples

Let’s plot a single tracked electron for these runs we didn’t track the ions.

Focus just on one run for simplicity

myRun = runs[0]

# plot t vs gamma for a random prtl

choice = np.random.randint(len(myRun.trackedLecs))
randPrtl = myRun.trackedLecs[choice]

# Each prtl has the following attributes: 'x', 'y', 'u',
# 'v', 'w', 'gamma', 'bx', 'by', 'bz', 'ex', 'ey', 'ez'

plt.plot(randPrtl.t, randPrtl.gamma)
plt.show()

This is nice, but let’s say you want to find the 10 highest energy prtls you saved. trackedLecs allows you to sort the database.

As an example, let’s sort by energy. You can pass any function here

myRun.trackedLecs.sort(lambda x: np.max(x.gamma))
# now plot the botton N
for prtl in myRun.trackedLecs[:-10]:
    plt.plot(prtl.t, prtl.gamma, 'lightgray')

for prtl in myRun.trackedLecs[-10:]:
    plt.plot(prtl.t, prtl.gamma, 'k')

plt.show()

You can also apply a mask to your particle to choose ones matching a certain criteria. You can pass any function that returns a truthy value

myRun.trackedLecs.mask(lambda x: np.max(x.gamma)>10.1)
plt.subplot(211)
for prtl in myRun.trackedLecs:
    plt.plot(prtl.t, prtl.gamma)
# Masks are applied successively. However you can unmask.
# Let's plot all the other prtls
myRun.trackedLecs.unmask()
myRun.trackedLecs.mask(lambda x: np.max(x.gamma)<10.1)
plt.subplot(212)
for prtl in myRun.trackedLecs:
    plt.plot(prtl.t, prtl.gamma)

plt.show()

All of this code can be found in examples.py in the main source directory.

Support for other pic sims.

It’s easy enough to tweak tristanSim to create an ORM for other simulations. I have done so for PICTOR, although unfortunately it uses the same keys that pictor uses. For example sim[n].ex in the TristanSim class becomes sim[n].Ex in the PictorSim class, and sim[n].xi becomes sim[n].x[sim[n].flv==1].