Many types of accelerator error studies are highly parallel due to the fact that each run of the simulation does not depend at all on any other runs. For example, studying the effect of a random misalignment of an accelerator component by running a large number of studies with different misalignments. This type of study can be easily performed with Pynac as follows.
The strategy here will be to use the multiProcessPynac function from Pynac.Core to build a queue of different simulations, and to execute them in parallel.
First we import all the functionality that we need from general libraries, and then make sure we're working in a clean environment (i.e., no files hanging around from previous runs).
In [1]:
import numpy as np
import os
import glob
import shutil
import random
from contextlib import contextmanager
def clean():
for outputdir in glob.glob('./dynacProc_*'):
if os.path.isdir(outputdir):
shutil.rmtree(outputdir)
clean()
Then import the necessary Pynac functionality.
In [2]:
from Pynac.Core import Pynac, multiProcessPynac, makePhaseSpaceList, getNumberOfParticles
from Pynac.Plotting import PynPlt
from Pynac.Elements import Quad, CavityAnalytic
Then import some functionality for plotting results.
In [3]:
from bokeh.io import output_notebook
from bokeh.plotting import figure, show
from bokeh.models.ranges import Range1d
from bokeh.models.axes import LinearAxis
output_notebook()
In [4]:
newDir = 'workingDir'
if os.path.isdir(newDir):
shutil.rmtree(newDir)
os.mkdir(newDir)
filelist = [
'ESS_with_SC_ana.in',
'Spoke_F2F_field.txt',
'MBL_F2F_field.txt',
'HBL_F2F_field.txt',
'ESS_RFQ_out_70mA.dst',
]
for f in filelist:
shutil.copyfile('../tests/' + f, newDir + '/' + f)
os.chdir(newDir)
We need a function that adds random errors to some components, and then runs Pynac. The function applyErrorsAndRunPynac applies a random error to the fields of all the quadrupoles, and then random errors to the amplitude and phase of all the accelerating cavities. It then runs Pynac on the resulting lattice.
In [5]:
random.seed(19781216)
def getRandomScaling(scale = 1e-3):
return random.gauss(0, scale)
def applyErrorsAndRunPynac():
a = Pynac('ESS_with_SC_ana.in')
for ind in a.getXinds('QUADRUPO'):
quad = Quad(a.lattice[ind])
quad.scaleField(1 + getRandomScaling(scale = 1e-2))
a.lattice[ind] = quad.dynacRepresentation()
for ind in a.getXinds('CAVMC'):
cav = CavityAnalytic(a.lattice[ind])
cav.adjustPhase(getRandomScaling(2.0))
cav.scaleField(1 + getRandomScaling(scale = 1e-3))
a.lattice[ind] = cav.dynacRepresentation()
a.run()
Now run a large number of iterations (400 in this case) of this function. I've chosen 8 workers on my dual-threaded, quad-core, Macbook, since this gives me the best speed improvement. Your mileage my vary -- experiment a little to find the optimal number of workers for your machine.
Note that the multiProcessPynac function takes care of creating and populating the necessary new directories in order to prevent the different parallel runs clobbering each other. It needs to know which files are necessary to do the Pynac run, and so these are provided in the filelist input variable.
In [6]:
%%time
filelist =['ESS_with_SC_ana.in',
'ESS_RFQ_out_70mA.dst',
'Spoke_F2F_field.txt',
'MBL_F2F_field.txt',
'HBL_F2F_field.txt']
print(multiProcessPynac(filelist = filelist, pynacFunc = applyErrorsAndRunPynac, numIters = 400, max_workers = 8))
In [7]:
@contextmanager
def workInOtherDirectory(dirname):
presentDir = os.getcwd()
os.chdir(dirname)
try:
yield
finally:
os.chdir(presentDir)
With the help of this context manager we scan through all the directories created by the parallel Pynac run (notice the 'dynacProc_%04d' % dirnum directory signature), and extract information for plotting.
In [8]:
xemit = []
emitAtEnd = []
for dirnum in range(10000):
try:
with workInOtherDirectory('dynacProc_%04d' % dirnum):
psList = makePhaseSpaceList()
xemit.append([ps.xPhaseSpace.normEmit.val for ps in psList])
emitAtEnd.append(xemit[-1][-1])
except FileNotFoundError:
break
Now plot the data that we've just extracted.
Note that the first time I wrote this part, I had the plotting inter-twined with the data extraction but that was a mistake. The problem is that any time I wanted to change the appearance of the plots, I had to run the entire data extraction code again, and on large datasets this can take a long time. It is much better to separate the work of extracting the data from the filesystem and plotting that data.
In [10]:
p = figure(plot_width=600, plot_height=300)
for emit in xemit:
p.line(range(len(emit)), emit)
p.yaxis.axis_label = 'Normalised x emittance'
p.xaxis.axis_label = 'Emittance measurement device number'
show(p)
p = figure(plot_width=600, plot_height=300)
hist, edges = np.histogram(emitAtEnd, bins='auto')
p.quad(top=hist, bottom=0, left=edges[:-1], right=edges[1:])
p.xaxis.axis_label = 'Normalised x emittance at linac end'
p.yaxis.axis_label = 'Number of simulations'
show(p)