The aim of a neutrino factory is to produce intense beams of neutrinos by trapping muons in a storage ring until they decay. Muons, having a halflife τ_{½}~1.56μs, are not readily available in nature and must be produced earlier on in the accelerator system. The suggested way of doing this starts with a proton beam of several GeV and fires it at a target material, which emits pions (π^{+}) that decay into the required muons (μ^{+}) soon afterwards (τ_{½}~26.0ns).
The main difficulty with this scheme is that the distribution of pions emitted from the source has an extremely large emittance: in the dataset used, 22% of the pions are actually going backwards in the lab frame (v_{z}<0), so the emittance can hardly be calculated. This means that any simulation based on conventional assumptions about the beam's general direction and size are likely to be inappropriate. There are inevitably large losses associated with capturing this beam and focussing it into a form that could e.g. be accepted into a linac.
The decay of the pions into muons adds another unusual element to the system, as the other decay products are emitted at high speed and make significant changes to the original particle's mass and velocity. The pion halflife corresponds to 7.8 lightmetres, so the decays happen on the same timescale as it takes the particles to travel through the first few magnets. This effect also makes simulation with a traditional tracking code potentially inaccurate.
This report focuses on the proposed design in which the pion source is a tantalum rod (20cm long, 1cm in radius) being hit by a 2.2GeV proton pulse, assumed to have a Gaussian longitudinal bunch shape of RMS one nanosecond. The first component in the muon capture apparatus is a channel of transversefocussing solenoids, the first of which encloses the rod (although the axis of the rod need not be parallel to the axis of the solenoid). The lack of longitudinal focussing means that the particles with larger v_{z} will work their way to the front of the bunch, which will have spread out somewhat by the end of the channel. With this in mind, the next component is a 'chicane' of speciallydesigned bending magnets, in which the faster particles will go on paths further to the outside of two opposing bends in the accelerator.
Figure 1. Schematic of the RAL muon front end design. 
The design is such that muons emitted directly down the solenoid channel with energies in range 120270 MeV are meant to take an equal amount of time to get to the end of the bends. This means that they are rebunched longitudinally at that point, ready to be accepted into the third component of the design, an 88MHz linear accelerator (not yet simulated).
These components constitute a muon front end for the Neutrino Factory, with bunched muons leaving the system with a mean kinetic energy of 400MeV. Other front end designs use different combinations of techniques to achieve the same output, for instance by using muon 'cooling' or several stages of linac acceleration.
Spacecharge effects are small in this system, due to the size of the beam (~10cm across) and low number of particles: five times fewer μ^{+} than there were p^{+} in the proton driver. In the decay channel μ^{+} and μ^{−} are present in comparable numbers, further cancelling the effect, so the simulations did not include a spacecharge solver.
This is done using a fourthorder RungeKutta (RK4) solver for each particle in the (x,y,z,v_{x},v_{y},v_{z}) phasespace, all coordinates being relative to the lab frame. A timestep was chosen by running batches of test particles through the solenoid channel and comparing the final positions with those obtained by halving the timestep. With the 10picosecond steps used in the real simulations, it is estimated that after 120 nanoseconds (approximately the total time the beam is simulated for), no errors in position of a particle will be more than 2.5mm, with RMS and mean errors 0.5mm and 0.3mm respectively. This compares favourably with the size of the beam, which is typically at least 100mm across.
The way that forces are applied to the particles is specialrelativistic and also fully integrated with the RK4 solver. RK methods in general solve firstorder differential equations, whereas equations of motion are secondorder in normal space. However, we find that the first timederivative of the position in phasespace (d/dt)(x,v) = (v,dv/dt) is a function only of the original phasespace vector (x,v)  and possibly t, in the case of varying fields such as RF cavities. The calculation of dv/dt is complicated as usual by the relativistic massincrease, but eventually one can derive the formula dv/dt = (F−c^{2}(F·v)v)/m_{0}γ, where F is the force as measured in the lab frame, i.e. q(E(x)+v×B(x)) with q=+e for our μ^{+} and π^{+}. These formulae are used in the program whenever the RK4 solver requests an evaluation of the derivative at some place in phasespace.
Currently the Bfield contribution from solenoids is calculated using an analytic solution to Maxwell's equations with cylindrical symmetry, which yields the field as the sum of an infinite powerseries. An error analysis similar to that used for determining the RK4 timestep was done, which revealed that truncating this series at the cubic term gives maximum and RMS positional errors of 0.25mm and 0.01mm respectively, so including any higher terms would slow down the simulation with no significant gain in accuracy.
The decay of pions is a random process in time, and the muon is emitted in an isotropically random direction in the rest frame of the original pion. Both of these factors would be better suited to a regime in which probabilitydensity functions f_{π}(x,v), f_{μ}(x,v) of phasespace are simulated, rather than being sampled discretely by macroparticles. While the development of such Vlasov solvers has been investigated, meaningful simulation of a 6D phasespace with reasonably high grid resolution has not been achieved.
To mitigate the noise introduced by the random decay of discrete particles, each macropion is decayed 10 times, producing 10 macromuons of onetenth the weight. Simulations without this feature are moderately accurate but the optimiser often compares the merit of quite similar results, making random variations more of a problem. When the bending chicane was added, the final muon transfer fell to the point where the optimisation was hindered by noise unless multiple decays were used.
When calculating the longitudinal beam profile, there needs to be some mechanism to simulate that the original proton pulse is not instantaneous. This is done by simulating several (usually 5) copies of the initial pions, but with the release of any one pion delayed by a random amount of time taken from the normal distribution t~N(μ=5ns,σ=1ns). In cases where only the number of particles retained is important, all the particles may leave the source simultaneously, as without spacecharge they do not interact. (They may also be simulated individually).
The simulation must take account of the various components situated in and around the beamline through which particles cannot pass. One of these is the rod source itself; others are the regions outside solenoids and the rectangular apertures of the bending magnets. Care must be taken that while stopping particles that pass outside a solenoid (for example) by a few inches, particles in a different part of the beamline intersecting the same plane (and thus technically 'outside' the solenoid) will not be removed. If this is implemented wrongly, particles can make false collisions with components.
After each run of the simulation, the result can be appended to a text file in the twoline format shown below.
zsrc=584;rotsrc=080;d1l=322;s2f=997;s3f=883;s4f=577;snf=490;alt=415;narrow=658;na1=422;na2=962;na3=474; 2.060580 (20936 particles) [v4.11] {E416DD6E}
The first line contains a list of named parameters that determine the accelerator design within the range being considered in the current optimisation. These have normalised values from 000 to 999, which can be rescaled onto the feasible range for that design parameter to find the actual value used. They may also be used to switch between different schemes of designs  for instance, the 'narrow' parameter is used to switch between four different algorithms for varying the solenoid radii along the decay channel.
The second line starts with the 'score' for the accelerator. From the optimiser's point of view, this just has to be a valid floatingpoint value with larger values corresponding to better designs (e.g. percentage efficiency). The rest of the line shows the number of particles used, the version number of the simulation software and a checksum of the beamline design file (a comment may also be added). This information is useful for distinguishing relevant results from those generated under different conditions or earlier (perhaps incorrect) versions of the software.
Optimisers are often implemented with steepest descent or hillclimbing algorithms that use the derivative of the score with respect to the parameters to determine which direction the optimiser will search along next. Unfortunately effects from pion decays add a small random variation to the score at any fixed design, making the calculation of derivatives impractical.
Instead, the optimiser uses an idea based on genetic algorithms (GAs). A conventional GA starts by running a number of random designs (the first generation) and produces a new generation using combinations of features from the best in the previous one, plus a small amount of random 'mutation'. This process is repeated until the score no longer increases. That would be far too wasteful for the purposes of this simulation, however, as it throws away the results from each generation after it generates the next (and a single simulation takes several minutes). Instead the optimiser logs all of the designs simulated (good or bad) to a file, and chooses a new design based on that data.
When the number of records in the file is small (e.g. less than 10), the program simulates a random design. Otherwise, it chooses one of the following four options with equal probability.
Since the results from independent optimisations are all still valid, this optimiser can be deployed as a distributed computing project, with nearindependent calculations being made on many different machines. This resulting project status page is currently at http://stephenbrooks.org/muon1. At times a "best results" or "results sample" file for the whole project is made available so that those who download it can base their optimisation on results found by other computers.
Before work on these simulations started, the muon front end had been developed (*) using standard modelling codes such as TRACE3D and Parmila. The basic design was as follows.
Element S = Solenoid, D = Drift space  Length (m)  Field (T)  Aperture Radius (m) 

S1  0.4066  20.0  0.1 
D1  0.5718  
S2  0.4  3.3  0.3 
D2  0.5  
S3  0.4  4.0  0.3 
D3  0.5  
S4  0.4  3.3  0.3 
D4  0.5  
S5  0.4  3.3  0.3 
D5  0.5  
S6D23  9 more repeats of the above block  
S24  0.4  3.3  0.3 
D24  0.5  
S25  0.4  3.3  0.15 
D25  0.5  
S26  0.4  3.3  0.15 
D26  0.5  
S27D34  4 more repeats of the above block 
The tantalum rod pion source is positioned concentric with S1 and rotated by 0.1 radians about the Yaxis. This rotation reduces the number of pions lost from reimpacting the rod after one turn of a helical path in the solenoid field. Some pions are moving forward rapidly enough to clear the rod before returning to the same transverse location, but for those that are not, the shifting transverse footprint of a tilted rod with distance makes it possible to avoid a collision. Having the rod at too small an angle is also undesirable when wanting to remove the remaining proton beam, which has not interacted with the target, from the system.
The results from one preliminary simulation to determine the effect of rod angle on pion capture are graphed in figure 2.
Figure 2. How rod angle affects the performance of S1. 
This simulation was made with only the solenoid S1 in place, so the 'beam kept' figure counts the particles emerging from the correct end of this solenoid and 'wrong way' counts those from the opposite end. The graph suggests that a rod angle of ~0.45 radians would produce the lowest loss and avoid high reinteraction between the emitted particles and the rod.
For reference, the original design was tracked with the code 36 times and yielded a transmission efficiency of 3.11±0.01% μ^{+} out per π^{+} emitted from the rod, where the range quoted is a 95% confidence interval over the variation caused by random simulated pion decays. Conversions to other figures of merit for all efficiencies quoted here are found in the appendix, while units of μ^{+}/π^{+} will be used here.
Two methods of optimising the solenoid channel have been investigated. The first generalises the design given in §3.1 by varying the 12 parameters shown below.
Name(s)  Description  Range of values  Notes 

zsrc  Zdisplacement of the rod source  0.15 to +0.15 metres  Measured as the position of the centre of the rod relative to the centre of S1. 
rotsrc  Rotation of the rod source about the Y axis  0 to 0.5 radians  
d1l  Length of drift space D1  0.5 to 0.6 metres  The drift spaces in this area of the accelerator are primarily determined by other instrumentation, but this one can be varied slightly. 
s2f s3f s4f  The strength of solenoids S2, S3 and S4  ±2.5 to ±5.0 Tesla  The sign as well as the strength can be changed by the optimisation; very low fields are not worth considering. 
snf  Magnitude of field in solenoids S5S34  ±2.5 to ±4.0 Tesla  S5's field has the same sign as this. The fields in later solenoids will have the same magnitude but sign controlled by alt. 
alt  Controls how the solenoid fields S6S34 change sign  Discrete; two values  Either the solenoids alternate in sign (+++), or they alternate in groups of five (++++++++++). 
narrow  Method of narrowing the solenoids S10S30  Discrete; four values  The solenoids narrow from 30cm to 15cm radius in one step; two steps via 22.5cm; three steps via 25cm, 20cm; or continuously in a section of linearlydecreasing radii. 
na1 na2 na3  Solenoid narrowing points  Solenoid number; integer from 10 to 30  These control the points at which the stepnarrowings take place, or where the linear narrowing begins and ends. na2 and na3 may be unused. 
When this phase proved successful, the design range was extended to ~130 parameters by varying every component length, field and radius independently. The following practical restrictions defined the design space:
The first optimisation is a subset of this range, so the second ought to produce a yield at least as good, with the amount of improvement indicating how suitable the extra design assumptions in the first phase were. The limiting case is when no improvement is made, meaning the optimal design lies within the first subset.
The highestscoring result recorded was the following.
Parameter(s)  Value(s) 

zsrc rotsrc  0.01006 metres, 0.03253 radians 
d1l  0.5172 metres 
s2f s3f s4f  +4.970, +4.660, +3.849 Tesla 
snf alt  +3.985 Tesla, signs alternate in groups of five solenoids (++++++++++ etc.) 
narrow na1 na2 na3  Solenoid radii decrease in three steps: to 25cm at S10, 20cm at S13 and to 15cm at S30. 
This was given a score of 7.066% μ^{+}/π^{+}, but due to statistical effects (and the fact this was the best result, so likely having an upward fluctuation), the true figure is probably closer to 6.5%. These effects were more severe than in the 'highgenerality' optimisations discussed later because an earlier version of the code, without multiple decays for each pion, was being used.
Unlike the preliminary simulations, the optimisation of the full decay channel seemed to favour small rod rotation angles, perhaps because capture of the pions by S2 and S3 is more difficult when the rotated rod produces a larger Xaxis spread of particles. The tendency towards choosing the strongest possible focussing solenoids for the majority of the beamline is not surprising given that the goal is to minimise beam loss, however the low value of s4f suggests a higher field at that point might not be beneficial.
Moving onto a larger design space, the program reached a design whose transmission was confirmed by multiple runs to be 9.65±0.01% μ^{+}/π^{+}. As before, many of the parameters were near the ends of their allowed ranges. However, not all were exactly at the extreme value, which is understandable given how difficult a yield change caused by a <0.01T field variation in one magnet is for the optimiser to detect in the presence of noise. Thus to simplify the design presented here, parameters within 10 quanta of the ends of the range were rounded to the ends. The resulting design also scored 9.65% when tracked and is listed below.
Element  Length (m)  Field (T)  Aperture Radius (m) 

S1  0.45  20.0  [0.1] 
D1  0.5  
S2  0.6  5.0  0.3514 
D2  0.5  
S3  0.4230  5.0  0.4 
D3  0.5  
S4  0.3806  4.189  0.4 
D4  0.6612  
S5  0.5299  4.0  0.4 
D5  0.5075  
S6  0.6  3.824  0.4 
D6  0.5170  
S7  0.6  4.0  0.4 
D7  0.5  
S8D25  18 more repeats of the above block, apart from three exceptions:  
S13  0.5956  4.0  0.4 
S14  0.6  4.0  0.3862 
D20  0.5175  
S26  0.6  3.904  0.4 
D26  0.6041  
S27  0.6  4.0  0.4 
D27  0.5  
S28  0.6  4.0  0.3159 
D28  0.5  
S29  0.2  4.0  [0.15] 
Parameters not at the extremes of their range are shown in italics. Those fixed by the optimisation are shown in [brackets]. 
The tantalum rod in this design is rotated by 0.09359 radians, with centre 0.2252 metres from the beginning of S1, which makes it nearly concentric as S1 is 0.45 metres long.
Figure 3. Shape of the pion beam on entering S2. 
As before, the optimiser has chosen to have the strongest possible focussing for the majority of the beamline by using the longest solenoids with the maximum field and minimising the length of the drifts. The first few components have been optimised in other ways, perhaps to account for the unusual shape of the beam in the early stages of the channel (see figure 3).
The strength of the field in S4 is substantially lower than the maximum allowed, just as happened with the choice of s4f in the last optimisation, which indicates this is a genuine feature and not just the optimiser converging badly. As another check, the region from S5 to D28 in the above design was 'homogenised' to look entirely like the repeating block and tested with the code. This produced a yield of only 9.09%, which is significantly below the original. However, when the three 'exceptions' from the middle of the beamline were removed, a short retest produced 9.68%, meaning that they could probably be omitted.
Perhaps the most significant change of the lattice has been in the polarity of the solenoids: while the previous optimisation had the sign change every five solenoids instead of the alternation present in the original, this one has taken the process to its conclusion and does not have any changes of sign. Considering the transverse motion of particles, we see that particles travel in straight lines in drift spaces and around circular arcs (in the transverse projection) while in solenoids. For a fixed energy, there are closed transverse paths possible for both the alternating (FODO) and the samesign (FOFO) solenoidal lattices, resembling a bowtie and an oval racetrack respectively. Higher energy particles will have a larger radius of curvature in the solenoids and will turn through a smaller total angle (and lower energy particles the opposite). This means that the path of an offenergy particle will not close. In the FOFO lattice this causes the racetrack orbit to precess because of the extra rotation on each circuit. However, due to the opposing bends in the FODO bowtie, this orbit suffers a translation each time it fails to close, shifting offenergy particles ever more offaxis and dispersing the beam. The precession process, by contrast, remains bounded in transverse space since any translation that may be present in the orbit defect will eventually be cancelled when it has rotated by 180°.
The way this was found illustrates the ability of a highdimensional optimiser to generate solutions qualitatively different to the reference design, sometimes even incorporating design principles that (as in this case) were not foreseen by the authors.
The second optimisation presented here varies the solenoid channel in the same way as before but also tracks the particles through the second segment of the design. The score for each configuration is calculated from the number of particles retained at the end of this second part, meaning that although the chicane itself cannot be modified due to beam optics considerations, the decay channel can be refined to be optimally compatible with it.
The bend consists of four identical 'triplets' of magnets (B1, B2, B1) with fields such that a pencil beam of muons with varying energies emerges from a triplet as a dispersed beam and rotated by 101°. When that dispersed beam enters a second triplet at the correct angle, it is turned back into a pencil beam (up to the correctness of the field within the magnets) but by then the higherenergy muons will have travelled further, giving the slower ones time to catch up. To cancel the dispersion produced by the decay channel, 404° of bends are needed, so to avoid layout problems the second pair of triplets bends in the opposite direction.
Figure 4, produced using the same tracking code as used for the optimisation, shows how the magnets will bend the paths of muons with energies of 120 MeV (red) up to 270 MeV (violet) in steps of 10 MeV.
Figure 4. The reverse phase slip bending chicane in action. 
Not all the particles make it through the chicane: the iron poles of the magnets are at roughly ±15cm from the midplane, close enough to create losses, especially since vertical focussing is quite weak. This is modelled by placing aperture restrictions at either end of each sector magnet.
The reference design for the solenoid channel (§3.1) was tested with the chicane added and scored 1.127±0.006% μ^{+}/π^{+}, giving the chicane an apparent efficiency of 36%.
A preliminary simulation calculated the distribution of muon arrival times at the end of the bending chicane, to see if the beam will be sufficiently bunched to be accepted into the 88MHz linac.
Figure 5. Longitudinal beam profile upon rebunching by the chicane. 
As figure 5 shows, the bunch emerging from the chicane is approximately the same length as the linac bucket (an 88MHz wavelength is shown for comparison), though some loss of the distribution tails will occur.
A lot of the apparent bunch length reduction in the chicane is due to energy selection by the horizontal aperture of the magnets when the beam is dispersed, since the particles with extreme energies tend to be at the ends of the bunch. For example, in one configuration the raw output of the solenoid channel has an RMS duration of 9.17ns, which reduces to 2.88ns for the 59.9% of the beam that is within the chicane's design acceptance of 120 to 260 MeV. If only the sub120 MeV stragglers are removed, the figure is just 3.40ns for 90.7% of the beam.
On leaving the chicane, the bunch has an RMS duration of 2.11ns in which 93.1% of the particles are of the design energy. If only particles in that range are considered, the bunch extent has reduced by 30% from 2.88 to 2.01ns, which is due mainly to the dynamics of the chicane rather than losses.
The design from the 12parameter optimisation of the decay channel (§3.3) was run with the bending chicane appended, to provide a comparison between separate and joint optimisation, and produced a yield of 1.878±0.027% μ^{+}/π^{+}. Then the full optimisation, scoring by transmission through the end of the chicane instead of just the solenoid channel, generated a design with a transmission of 1.930±0.013%. This is higher by a statistically significant amount, but the code without multiple pion decays was still being used at this point, and with the lower yields coming from the chicane, noise started to impair the optimiser, so these figures may not be quite optimal.
Increasing the number of parameters improved the transmission to 2.704±0.004% μ^{+}/π^{+}. One comparable result for separate optimisation is the design of §3.4 with the chicane added, which scores only 2.004±0.005%. Closer examination reveals that the restriction of fixing the final solenoid's radius in this lattice (which was meant to be a substitute for the finite aperture of the chicane) actually produces a betatron focus there when optimised, in which the transverse velocities of the particles are still large. Removing this restriction gives a decay channel transmitting 15.80±0.01% on its own and 2.238±0.004% through the chicane. Still, jointly optimising the two sections has produced an improvement of just over 20%: not insignificant given that the full neutrino factory has many sections and a compounded 20% improvement at each interface would correspond to an overall factor of 2× or 3×. This performance would only be accessible via holistic approaches such as optimising using a (currently hypothetical) code capable of tracking through the whole machine, or tuning the magnet strengths on the real machine after it is built, the latter being able to recover only the fraction not requiring machine geometry changes.
Only a few of the components in the new design deviate from the pattern of using the highestpossible solenoid fields and the shortestpossible drifts for maximum average field onaxis. There are some weaker solenoids near both ends of the channel, as has happened before for matching with the input and output. The only nonminimal drift could well be an artefact of the optimiser having not quite converged.
Solenoid  Field (T)  Drift  Length (m) 

S1  20  D1  0.5 
S2  5  D2  0.5 
S3  3.088  D3  0.5 
S4  4.429  D4  0.5 
S5S25  4  D5D25  0.5 
S26  3.752  D26  0.5 
S27  4  D27  0.5345 
S28  3.760  D28  0.5 
S29S31  4  D29D31  0.5 
S32  3.431  D32  0.5 
The other parameters in this design are of more interest. In the decaychannelonly simulations, the program predominantly opted for the (logical) approach of using the largest possible solenoid radii. By contrast, the new design has an evident solenoid radius structure (see figure 6).
Figure 6. Solenoid radii in the highgenerality optimisations. 
The decreasing trend is fit well by the linear formula 0.402−0.003n for the radius of Sn from S7 to S28. Note that the optimisation algorithm itself does not treat the parameters as being in any particular order, so trends such as this cannot be artificial and must come from the beam dynamics. Secondly, note that if hardended solenoids had been used in the simulation, narrowing the radii in this way would have produced a result no better than keeping them at maximum size until the very last solenoid, because particles would experience the same fields in such a linear model in either case until they are lost. This suggests that the nonlinear part of the solenoid endfields is being used to achieve some sort of additional focussing.
The fact that the solenoids have not narrowed all the way to 15cm is another indication that the use of that design constraint in §3.4 was not appropriate.
Figure 7. Solenoid lengths in the highgenerality optimisations. 
Figure 7 shows that the optimiser has settled on an intermediate value for the solenoid lengths in the main part of the channel, in contrast to how it made them all as long as possible without the chicane. The median of the lengths used is 0.463m.
Some analysis of regular solenoid focussing channels using linear optics was done in order to explain this figure. The Larmor rotation within each solenoid is given by
θ = 
 ,  q, m_{0} = Particle charge, rest mass B_{z}, L = Solenoid field, length γ, v_{z} depend on particle velocity 
which equals 160° in the case of 120MeV muons down to 90° in the case of 260MeV muons. Angles that are too small (or more generally, too close to a multiple of 360°) will mean particles will only have rotated slightly each time the come to a drift space, leading to a large transverse orbit. The FOFO lattice is also not unconditionally stable, as if the drift length ratio D:= L_{drift} / L_{solenoid} is greater than a certain value
D_{critical} = 

 , 
the transfer matrix becomes unbounded when iterated. For this lattice, D=1.08 so we must have θ<192° for stability (discounting solutions θ>360°). The range 90°160° can be seen as avoiding both of these problems by a suitable margin.
Finally, to check how much of the optimiser output is "noise" surrounding the trends already identified, rounded versions of the lattice were checked, described by the table below.
Solenoid  Length (m)  Field (T)  Aperture Radius (m) 

S1  0.3949  20  [0.1] 
S2  0.6  5  0.3093 
S3  0.5235  3.088  0.4020.003n 
S4  0.463  4.429  
S5  4  
S6  ^{1}0.463, ^{2}0.5584  
S7S30  0.463  
S31  0.3461  
S32  0.4651  3.431 
All drifts were set to 0.5m long. There was some doubt as to whether the particularly long length of S6 (see figure 6) was genuine or an optimiser error. In the first simulation, the solenoid lengths were rounded as far back as S4 but this decreased the transmission to 2.678±0.008%, which is lower than the original with a statistical significance of 5.7σ. In the second rounded design the value of S6 (but not S5 or S4) was restored, which gave an efficiency of 2.704±0.006%, just as good as the nonrounded design.
Recently interest has moved onto a new sequence of components after the decay channel (*), consisting of a muon cooling ring with absorbers to reduce the emittance of the beam over several turns, preceded by an RF system to phaserotate the decay channel output into a smaller energy spread. At the other end of the beamline, it has become clear that various practicalities demand the target area be changed from the model used here. These and other work to be completed with this code are listed in detail below.
The simulations described here used a dataset generated for 10^{5} 2.2GeV protons hitting the rod. The total π^{+} weight produced was 21285.3, so conversion is done at 0.212853 π^{+}/p^{+}.
Beamline  Efficiency in μ^{+}/π^{+} 95% C.I. uncertainties for last digits in brackets  μ^{+} per proton on rod  μ^{+} per proton.GeV 
Original decay channel  3.114(10)%  0.6627(22)%  3.012(10)×10^{3} 
12parameter optimised  ~6.5% Simulated with earlier, lessaccurate version of the code  ~1.4%  ~6.3×10^{3} 
Highgenerality optimised  9.645(13)%  2.053(3)%  9.332(12)×10^{3} 
Original to end of chicane  1.127(6)%  0.2399(13)%  1.090(6)×10^{3} 
12parameter optimised  1.930(13)%  0.4107(28)%  1.867(13)×10^{3} 
12parameter for decay channel through chicane  1.878(27)%  0.3998(57)%  1.817(26)×10^{3} 
Highgenerality optimised  2.704(4)%  0.5756(8)%  2.616(4)×10^{3} 
Highgenerality for decay channel through chicane  2.004(5)%  0.4267(11)%  1.939(5)×10^{3} 