# Presentation: 556660

## Initialization of Beam Simulations using Particle Distributions Reconstructed from Experimental Data*

**Alex Friedman and David P. Grote, *****LLNL***
**Frank M. Bieniosek, Christine M. Celata and Lionel R. Prost, *****LBNL***
**Presentation BP1.110 **
**Annual Meeting of the APS-DPP****November 10-15, 2002**
**Orlando, FL**

- *	Work performed under the auspices of the U.S. Department of Energy by the University of California, Lawrence Livermore and Lawrence Berkeley National Laboratories under Contract Nos. W-7405-Eng-48 and DE-AC03-76SF00098

## Abstract

**In present-day experiments exploring the behavior of intense particle beams for Heavy Ion Fusion, a variety of diagnostics yield information on the beam particle distribution function. Good shot-to-shot reproducibility enables the acquisition of detailed 2D projections of the 4D transverse phase space using moving slits and Faraday cups**
**In trials using simulated beam data, we have observed a significant improvement in simulation fidelity when runs are initiated using “tomographically synthesized” (as opposed to simple-model) distributions based on slit scan data **
**We are beginning to employ such synthesis techniques to launch simulations with initial conditions developed from diagnostics at an “upstream” station on the High Current Experiment at LBNL**
**	****See also Thursday afternoon poster RP1.065 by C. M. Celata, *****et. al.*****, and other posters in that session**

## Introduction: problem and methods

## Simulations can be carried out at three levels

**“****ideal” beams ... straightforward**
**“****realistic” beams ... harder**
**use ensembles of integrated runs for **_**planned**_** experiments**
**“****the real machine” ... difficult!**
**need detailed machine characterization, diagnostics**
**integrated simulations tend to “drift” from **_**actual**_** experiments**
**	********  ****iterate, varying misalignments, voltage errors, etc.**
**	********  ****couple with “synthesis” approach**

## Slit-scan diagnostics yield time-dependent data

_**Crossed-Slit Intensity Map:**_
**Used to examine the distribution of beam current density**
**A scintillator may yield similar data, while a kapton foil yields a time-integrated intensity map**

_**Two-Slit Emittance Scan:**_
**Used to measure the beam emittance in direction perpendicular to the slits**

**slit 1 selects for position x**

**slit 2 selects for transverse “momentum slope” x******** ******** p****x****/p****z**

**detector**

**The second paddle typically has a shallow, reverse-biased faraday cup attached to measure the beam signal**

## A longitudinally-cold steady flow can be described by the particle distribution ƒ(x,y,x,y) in the 4D phase space, where x =  px/pz and z is the axial coordinate
	In reality, beams have a small spread in pz and transverse motion isn’t independent of pz; nonetheless, this is a useful approximation
Experimental diagnostics typically yield 2D projections of ƒ(x,y,x,y), e.g., ƒ(x,x) = ∫ ƒ(x,y,x,y) dy dy
We want to launch simulations using an initial beam state derived from experimentally measured data; a goal is to reproduce the observed beam behavior in detail 
So we are trying to solve the problem:
	Given measured ƒ(x,x) and ƒ(y,y) projections of ƒ(x,y,x,y), and perhaps ƒ(x,y) as well, synthesize a “reasonable” ƒ(x,y,x,y)
	This problem is underdetermined

**A longitudinally-cold steady flow can be described by the particle distribution ****ƒ(x,y,x********,y********)**** in the 4D phase space****, where x******** =  p****x****/p****z**** and z is the axial coordinate**
**	****In reality, beams have a small spread in p****z**** and transverse motion isn’t independent of p****z****; nonetheless, this is a useful approximation**
**Experimental diagnostics typically yield 2D projections of ****ƒ(x,y,x********,y********)****, *****e.g.*****, ****ƒ(x,x********) = ****∫**** ƒ(x,y,x********,y********) dy dy******
**We want to launch simulations using an initial beam state derived from experimentally measured data; a goal is to reproduce the observed beam behavior in detail **
**So we are trying to solve the problem:**
**	****Given measured ƒ(x,x********) and ƒ(y,y********) projections of ƒ(x,y,x********,y********), and perhaps ƒ(x,y) as well, synthesize a “reasonable” ƒ(x,y,x********,y********)**
_**	**__***This problem is underdetermined***_

**Statement of the mathematical problem**

## The problem is underdetermined as stated; to make it well posed, we add constraints

**The available data is insufficient to fully specify the 4D distribution**
**For ƒ(x,y,x********,y********) on a 4D grid with, e.g., 20 points along each Cartesian axis, need 20****4**** = 1.6********10****5**** points**
**The measured (x,x********) and (y,y********) projections each consist of 400 data points, so we have a total of 800 data points, far fewer than is demanded by a complete description**
**We may also know the intensity ƒ(x,y); then we have 1200 points**
**We usually don’t know the other “Cartesian” projections (x,y********),     (y,x********), or (x********,y********); but even if we knew them, plus a few other projections at “angles” to the axes, it would not suffice**
**We impose additional constraints to complete the specification**
**We use a Monte-Carlo approach related to “maximum entropy” **
**This is augmented (optionally) by imposing a sampling region with “rounded corners”**

## This “2-plane” method reproduces, in the limit of many particles, the measured (x,x) and (y,y) data

**This “2-plane” method reproduces, in the limit of many particles, the measured (x,x********) and (y,y********) data**

**Assign target “counts” N(x,****x********) and N(y,y********)**** — the numbers of particles to be loaded into each “bin” — proportional to the measured ƒ in each bin**** **
**	****These bin counts should sum to the desired total number of particles**
**Repeat the following until all particles have been set: **
**Generate a random point in a 4-box, and accept it only if it ...**
**(a) falls in a pair of bins whose counts are nonzero, and **
**(b) falls within a user-specified “sampling” region**
**If it is accepted, decrement the counts in the corresponding bins**
**The sampling region is the intersection of:**
**four 4-cylinders, the extrusions of circles in (x,y), (x********,y********), (x,y********), & (y,x********) (nominally of radius 1 relative to the data extent, but adjustable)**
**a 4D hyperellipsoid (nominally of radius √2, but adjustable)**
**Method can “paint itself into a corner” (fail to load all particles desired)**
**One 3-plane method proceeds similarly, with addition of counts N(x,****y****) **

## An attractive 3-plane version uses bin counts N(x,y) only in the spatial plane

**An attractive 3-plane version uses bin counts N(x,y) only in the spatial plane**

**Compute bin counts N(x,****y****) ******** ƒ(x,y) mapped onto working grid**
**Correct probabilities to account for “rounded” 4D sampling region:**
**Compute V(x,x********) = number of (y,y********) grid points in sampling region**
**Compute ƒ*(x,x********) = ƒ(x,x********) N****y ****N****y******** H[V(x,x********) ] / [V(x,x********)+10****-99****] ;**
**	*****i.e.*****, divide the acceptance probabilities in ****(x,x********)**** by the area of the “slice” of the sampling region in the (y,y********) plane at that (x,x********)**
**Similarly for ƒ*(y,y********)**
**Until done, generate random points in the four-space; accept each if:**
**It falls inside the sampling region**
**ƒ*(x,x********)/Max(ƒ*(x,x********)) > Random(0,1) and similarly for (y,y********)**
**The “remaining” bin count N(x,y) > 0**
**	****If a point is accepted, decrement the corresponding N(x,y)**
**This algorithm avoids the “painting into a corner” problem**

**These algorithms, and all other processing, were implemented in the “very high level” scripting (interpreted) language Yorick, using vector constructs for efficiency (very few element-at-a-time loops are used) **

## Tests using simulated data

## Tests of synthesis at the exit of the HCX ESQ injector are based on a WARP3d run (so we know the 4D distribution)

**Tests of synthesis at the exit of the HCX ESQ injector are based on a WARP3d run (so we know the 4D distribution)**

- (frame from a WARP movie by J-L. Vay)

## Matching section compresses beam significantly    before it enters the HCX transport line

**Matching section compresses beam significantly    before it enters the HCX transport line**

- (frame from a WARP movie by J-L. Vay)

## Self-consistent simulations of HCX use WARP3d in injector, and WARPxy in matching section and transport line

**Self-consistent simulations of HCX use WARP3d in injector, and WARPxy in matching section and transport line**

**End of injector**

**End of matching section**

**50 lattice periods (100 quads)**
**past matching section**

- x  px/pz
- -0.002  0   .002

- -0.01        0          .01

- -0.01        0          .01

## Reconstructions of self-consistent simulated ƒ at injector exit

**self-consistent simulation (E. Henestroza)**

**2-plane synthesis using (x,x********) & (y,y********) data **

**3-plane synthesis using (x,x********), (y,y********), & (x,y) data **

## Comparison of HCX simulations shows benefit of initiating the run with a good approximation to ƒ(x,y,px,py)

**Comparison of HCX simulations shows benefit of initiating the run with a good approximation to ƒ(x,y,p****x****,p****y****)**

**HCX simulations initialized using synthesized 4D particle distributions**** based on ESQ-exit slit-scan data (simulated, here) are far better than those using an ****ideal distribution**
**An (x,y) scan****, in addition to (x,x********) and (y,y********), is important**

**z (m)**

**Normalized x emittance (********-mm-mr)**

**semi-Gaussian**

**Integrated simulation, starting at source**

**3-plane synthesis**
**adding (x,y) data**

**2-plane synthesis**
**using (x,x********) and (y,y********)**

**Notes:**
(This slide is a slightly updated version of the one shown at HIIF2002; the committee saw a version of this at the February meeting.) 

A simulation of HCX is used as the reference case since we know the "right answer" (black curve). Launching a simulation with semi-gaussian initial conditions (with a,a',b,b',emittances, etc. matched to the integrated simulation's parameters at injector exit) fails utterly (magenta curve). Use of a synthesized initial distribution based on x,x' and y,y' slit-scan data does better, but produces a spurious "hump" in the matching section and too slow a rise in the 40-quad HCX line (blue curve).  When a synthesized disttribution based on both slit-scan and spatial image x,y data is used to start the run, the behavior is quite close to correct, especially in the shape of the "wiggles" in the HCX line and in the asymptotic flatness; it fails primarily in that the jump at the exit of the matching section is too large (red curve).

Simulations using initial conditions synthesized from the actual HCX data are in progress.

## Synthesized ƒ leads to WARPxy simulations with higher fidelity than RMS-equivalent semi-Gaussian (y,y shown)

**Synthesized ƒ leads to WARPxy simulations with higher fidelity than RMS-equivalent semi-Gaussian (y,y******** shown)**

**semi-Gaussian**

**Integrated simulation, starting at source**

**3-plane synthesis**
**adding (x,y) data**

**2-plane synthesis**
**using (x,x********) and (y,y********)**

**Normalized y emittance (********-mm-mr)**

## Initial synthesis and simulation using experimental data

## High Current Experiment (HCX) layout

**Detail of diagnostics at transport line entrance (QD1)**

**beam**

**kapton**

**x  y**

**upstream**
**slits**

**x**

**downstream slits**** would collide, can’t be used simultaneously**
**crossed-slit map uses upstream       x slit and downstream y slit**

**y**

**(x slits are vertical, move horizontally from shot to shot; y slits are horizontal )**

**9.1 cm**

**start - simulation - end **

## Quite a few manipulations are needed to put the experimental data in a usable form

**Quite a few manipulations are needed to put the experimental data in a usable form**

**Preliminaries: **
**Read in raw data and transform to “HCX Global Coordinates”**
**Plot beam moments vs. time, and select a beam-frame slice for analysis / synthesis (WARPxy follows a single slice per run) **
**Apply “threshold” to experimental ****ƒ’s**** (to remove noise & bias)**
**Map emittance scan data to “working grid” (finer than the raw data grid)**
**Compute and then remove ********x********, ********x********, ********y********, ********y******
**Compute and then remove net beam envelope “convergence” in each plane, proportional to ********xx******** and ********yy******
**For each data column at constant x (which may start at its own unique x********),   fit a cubic spline and compute ƒ on working-grid x******** nodes**
**Then use a spline fit along x, for each row in x********, to get a 2D ƒ on working grid**
**Zero-out nodes outside the raw-data domain (spline fits can do weird things)**
**Normalize the data to range from 0 to 1, yielding ****ƒ****(x,x********) and ****ƒ****(y,y********)**
**Map crossed-slit data to working grid**
**Use the net convergence from the emittance scan (for each time slice) to rescale the y values, since the synthesis is done at a plane well upstream of the y slit**
**Threshold and normalize**
**Map kapton data from its native very-fine grid to working grid, for comparison**

## HCX data at exit of matching section (station “QD1”) at one time-slice is the preferred input to WARPxy simulations

**HCX data at exit of matching section (station “QD1”) at one time-slice is the preferred input to WARPxy simulations**

**emittance scan data**
**(after removal of net envelope convergence / divergence)**

**raw crossed-slit        ****(distorted because y slit is far downstream)**

**Time-integrated spatial intensity is “smeared”**

- Crossed-slit plots at time slice clipped to 85% of peak

- Dots denote raw data locations

**kapton**

**remapped crossed-slit**

**remapped crossed-slit**

**Data at a single time slice**

## First 3-plane synthesis using HCX data at entrance to    transport line was based on emittance scan and kapton data

**x,x ́ slit data**

**y,y ́ slit data**

**n(x)**

**x,y kapton data**

**x**

**x,y as synthesized**

**The kapton data is time integrated; since I****beam**** varied with time, that image is smeared and results were marginal**

**Notes:**
Measured f(x,x'), f(y,y'), and f(x,y) at the upstream end of the HCX transport line are being used to reconstruct an approximation to the 4D f(x,y,x',y').  Simulations using initial conditions so derived will help interpret the observed distribution downstream.  This reconstruction/simulation effort started immediately upon data acquisition, on the basis of much advance preparation.
The leftmost column shows experimental data from the parallel-slit "emittance" scanners (upper 2 images) and from the kapton film (lower image), all at the exit of the matching section (= entrance to the 10-quad HCX transport line).  The slit-scan data is shown at a single time slice near mid-pulse, as used in the synthesis; the mean convergence from the AG focusing has been subtracted.  The kapton data is time-integrated, with the low-density fringe outside the beam (especially at top and bottom) due to the head and tail which are mismatched to the channel. The spurious image at the left results from the kapton sharing a diagnostic "paddle' with a slit and so suffering damage when the slit is used.
The middle column, upper frame, shows the integrated quantity n(x) obtained by: integrating the kapton x,y data over y (black curve); integrating the slit-scan x,x' data over x' (blue curve); and computing the spatial density of the synthesized set of particles (magenta curve).  The lower frame shows the x,y density as reconstructed from the experimental data.  The top, bottom, and sides of the image are cut off because the slit-scan data indicates a smaller beam than that in the kapton image, because the slit is located 1 cm upstream of the kapton (see plots on upper left of page).
The rightmost two columns show a subset of the particles that are generated by the synthesis procedure, projected onto the six principal Cartesian axes). The scales differ from those used in the x,y plots of the other columns. Note that the synthesis uses a Monte-Carlo procedure that actually loads particles; densities such as those in the middle column are derived.

## 3-plane syntheses are now based on time-resolved emittance scanner and crossed-slit data, mapped to a common station

**3-plane syntheses are now based on time-resolved emittance scanner and crossed-slit data, mapped to a common station**

**effects of round-edge constraint**

**Projections of synthesized distribution agree well with input data; here, ƒ(x,x********)**

## This early “D-end” data shows the corresponding time-slice of the beam, at the exit of the 10-quad HCX transport line

- (data has been recentered; beam was sent off-axis by misaligned quadrupole lenses)

## Initial simulations are only in coarse agreement with early data at “D-end” of HCX 

**Initial simulations are only in coarse agreement with early data at “D-end” of HCX **

**remapped crossed-slit**

**emittance scan y,y******

**Spatial hollowing is a common feature**

**Notes:**
                              RMS Epsilon x     RMS Epsilon y 
Data (thresh 10%) 7.66639e-06      1.03185e-05 
Simulation             1.40659e-05      1.24266e-05

## Discussion and plans

## Future Plans

**We plan to track down the sources of discrepancy between HCX data and simulation output, in search of quantitative agreement**
**In addition to this “inverse” approach, we are also working the “forward” problem via integrated simulations starting from the source **
**To achieve agreement with experimental data at downstream stations, it may be necessary to vary (within the error bars) voltages, alignment, and/or source inhomogeneities.  **
**While the “forward” approach is more consistent, it is laborious, and we anticipate that much of our day-to-day work will be done using the mid-stream initialization procedures described here**
**Recently, slits have been combined with scintillator-based imaging methods, and the use of multi-aperture hole plates is under study **
**These methods yield 3D (projectional) and 4D (full) phase space information, but with relatively coarse sampling. **
**We plan to develop the capability to make use of such data to initialize simulations, and in general, to fold any extra information, *****e.g.*****, localized data at high resolution, into the synthesis**

## Extras

## Integrated WARP runs, compared w/ data

|  | Simulation | Pg. 108 data |
| --- | --- | --- |
| current | 555 mA | 175 mA |
| energy | 1.79 MeV | 1.045 MeV |
| beta | 9.9  10-3 | 7.6  10-3 |
| emittance (RMS) |  | 1.2  10-5 |
| emittance (edge) |  | 4.8  10-5 |
| norm emit (edge) | 5 - 12  10-7 | 3.6  10-7 |
| norm emit (edge) scaled by square root of current ratio |  | 6.5  10-7 |

## Spatial hollowing is a common feature

**spatial intensity x,y**

**phase space x,x******** view**

**phase space y,y******** view**

**WARPxy simulations initiated using a synthesis based on data at QD1: views of simulation phase space at D-end**

**In this run, the simulated distribution is only in coarse agreement with early HCX data**
**This work is quite new (Nov. 2002); we plan to track down the source of the discrepancy**
**Ongoing experiments are yielding well-centered beams with improved temporal “flat-tops” and better diagnostic resolution**

## HCX transverse phase space data at exit of matching section (station “QD1”) is input to WARPxy PIC simulations

**HCX transverse phase space data at exit of matching section (station “QD1”) is input to WARPxy PIC simulations**

**Emittance scan data**
**(after removal of net envelope convergence / divergence)**

**Raw crossed-slit map**

**Kapton foil damage pattern**

- Crossed-slit plots at time slice clipped to 85% of peak

- At time slice

- At time slice

- Time integrated

- Time integrated

**Crossed-slit data remapped to z(kapton)**

- At time slice

- At time slice

- Dots denote raw data locations