Observed proton beam induced disruption of a tungsten powder sample at CERN

Strathprints is designed to allow users to access the research output of the University of Strathclyde. Unless otherwise explicitly stated on the manuscript, Copyright © and Moral Rights for the papers on this site are retained by the individual authors and/or other copyright owners. Please check the manuscript for details of any other licences that may have been applied. You may not engage in further distribution of the material for any profitmaking activities or any commercial gain. You may freely distribute both the url (https://strathprints.strath.ac.uk/) and the content of this paper for research or private study, educational, or not-for-profit purposes without prior permission or charge.


I. INTRODUCTION
A continuously flowing fluidized tungsten powder jet [1] has been proposed as a potential target technology for future high energy physics facilities, such as a Neutrino Factory [2]. A gravity driven granular flowing tungsten target has been proposed for the neutron source of an accelerator driven subcritical reactor [3]. The Neutrino Factory concept requires a high power proton beam to interact with a target to generate mesons (pions) and eventually a beam of neutrinos through lepton (muon) decays. During the interaction significant heat is deposited in the target material which results in high temperatures and stress. The interest in flowing tungsten powder technology arose from its potential to accommodate very high deposited power densities while maintaining a reasonable operating temperature and low stress levels. Using small tungsten particles means they are almost uniformly heated by the beam which reduces stress, while constant circulation of the particles spreads the heat load over a large number of particles which are given enough time to cool between beam pulses. The fluidization of tungsten powder has been demonstrated off-line in a bespoke test facility [4], using helium as the carrier gas.
Two experiments [5,6], designed to investigate the response of the tungsten powder target, have been carried out at the HiRadMat facility at CERN [7,8]. The first experiment was performed in 2012, where a static tungsten powder sample of sub-250 μm crystalline particles in helium was exposed to 440 GeV=c proton beam pulses (2 mm sigma) of 1 μs duration. In 2015 a follow-up experiment was performed with a range of samples in helium at pressures of 1 bar and 20 mbar. This experiment was equipped with a sweeping mechanism which allowed the sample to be reset to its original unperturbed condition before another beam pulse was taken. In addition to the sub-250 μm crystalline particles five samples of spherical tungsten particles with mean diameters of 14, 34, 77 and 116 microns and 1 mm were also tested.
The main observation from the aforementioned inbeam experiments included: (i) powder expansion/lift proportional to the beam intensity, with a maximum observed lift velocity of 5.3 m=s for 34 μm spheres at 3 × 10 12 protons-on-target in vacuum corresponding to a peak energy density of 185 J=g; (ii) smaller tungsten particles lifted more than larger particles; (iii) lifts were observed both with helium at 1 bar and at 20 mbar but appeared smaller at 1 bar.
Following these experimental observations, three lift mechanisms have been investigated: (i) aerodynamicrapid heating of the interstitial helium gas leading to an expansion of the gas which escapes from the powder exerting a lifting drag force on the tungsten particles; (ii) thermal expansion-beam-induced thermal expansion of the tungsten particles, leading to the propagation of loads via force chains between adjacent touching particles and ultimately propulsion of some of the material; (iii) chargebeam-induced residual electrical charge on the tungsten powder leading to a force field acting on the particles as like charges repel.

II. AERODYNAMIC EFFECTS
The proton beam pulse extracted from the CERN Super Proton Synchrotron (SPS) hits the tungsten powder sample and results in significant energy deposition in the tungsten as protons are stopped, with many secondary protons and electrons released from their tungsten atoms which subsequently stop. A significant amount of the energy contained within the beam pulse ends up in the tungsten sample in the form of a sudden increase in temperature. We predict the amount of energy that will be deposited by the proton beam pulse using FLUKA [9]. Figure 1(a) shows the simple geometry, representative of the sample of tungsten powder, used for the simulations. The plane AA, which is at a distance of 11 cm from the front surface of the sample, approximately represents the axial position where maximum energy deposition occurs. Figure 2 shows a contour plot of the predicted energy deposition as a result of the beam interaction. As the energy is deposited the temperature of the tungsten increases and so does its heat capacity. Heat transfer also occurs in a very short timescale by virtue of the very short length scales between particles. Therefore, in order to calculate the expected temperature jump in the tungsten with any accuracy, temperature dependent material properties and heat transfer must be considered.
With the highest intensity beam pulse taken in the second experiment (3 × 10 12 protons) the tungsten is predicted to increase in temperature by more than 1500 K within 1 μs. The helium is also heated by the beam and thermal diffusion between the submillimeter-sized tungsten particles and the helium occurs extremely rapidly such that within the order of the beam pulse the helium will reach a very similar temperature to the adjacent tungsten particles. However, it will have had very little time to move within that 1 μs and so we could assume that the helium undergoes an almost isochoric heating process. This is expected to give rise to a large rise in the pressure of the helium. In the second HiRadMat experiment, tests were carried out in helium at 20 mbar and at a little over atmospheric pressure. Isochoric heating of helium at 1 bar by 1500 K leads to a pressure pulse of 5 bar on the beam axis, but 6 mm away at the edge of the beam spot there is very little heating and pressure rise, so a significant pressure gradient of the order 0.83 bar=mm is created. In the case of helium at 20 mbar the pressure pulse on the beam axis will be 0.1 bar, leading to a gradient of 0.017 bar=mm. The weight of tungsten powder gives rise to a downward pressure gradient of 0.00097 bar=mm. This simple comparison indicates that lift by virtue of increased gas pressure is a possibility. However, in order to assess if it is really the cause, one must consider the transient nature of the gas expansion. The high pressure helium will rapidly expand through the powder into the space above the sample and so the pressure will only act on the tungsten particles for a very short time. The escaping plume of helium may be fast and violent enough to lift some tungsten particles via aerodynamic drag. Whether or not an aerodynamic lift is possible depends on the starting helium pressure, the size of the particles and the beam intensity. To get a handle on this nontrivial transient problem we have employed a two-phase transient computational fluid dynamics (CFD) simulation [10] of a 2 mm thick slice of the tungsten sample at plane AA [ Fig. 1(b)]. The helium is treated as the continuous phase (subscript α), while the tungsten particles are treated as a dispersed solid particle phase (subscript β) within the continuous one. We now present the governing equations used, starting with the continuity equation for each phase: where r represents the volume fraction, ρ the density and U is the velocity vector of each respective phase. The momentum equations for the gas and solid phases respectively are as follows: They include the transient and advective terms on the left-hand side and the pressure gradient, dissipative viscous and gravity terms on the right-hand side, where p denotes pressure and μ denotes viscosity. The solid particle phase momentum equation includes a solid pressure term denoted by τ ij . This term depends on the solids pressure gradient as follows: where P s is determined from an empirical function of the solid phase packing fraction. The term deals with contact forces in the solid phase and ensures that the solid phase packing fraction cannot exceed unrealistic values. The terms M α and M β represent the momentum transfer between phases due to interphase drag and so M α ¼ M β . The interphase drag force is generally determined from the following equation: where c αβ depends on the flow regime, particle size and packing fraction of the solid phase. We employ the Gidaspow formulation based on the Wen Yu correlation at high void fraction and the Ergun equation at low void fraction [11,12]: The energy equation for each phase is where h tot is the total specific enthalpy and γ is the thermal conductivity.
The interphase heat transfer term is represented by Q and for the particle sizes we are considering is dominated by thermal conduction due to the short length scales of the gaps between particles. Clearly heat flow between phases has to be conserved, so The term S Eβ is the source of energy applied to the solid phase based on the FLUKA energy deposition data. Energy can be supplied to the continuum helium phase as well, but this is insignificant compared to the energy applied to the relatively dense tungsten particles. To achieve closure for the above CFD formulas, equations of state of the two phases are required. The continuous helium phase is treated as a perfect gas so that its pressure, temperature and density change physically as energy is transferred from the tungsten to the gas. The solid tungsten particles are assumed to have constant density and a temperature dependent heat capacity. The computational domain is a pseudo-2D domain at plane AA in the sample, see Fig. 1(b). The domain has an open top so that gas and powder can escape through the opening. The outside walls of the sample container are represented by no-slip walls. It is also important to note that these walls are assumed to be stationary, i.e. no movement of the sample trough. A symmetry boundary condition is used so that the computational domain only needs to be half the width of the sample container. The initial condition has a powder volume fraction of 50% from the bottom of the domain up to a fill level of 22 mm and above that the powder volume fraction is 0%. Gravity resists the powder lifting acting in the direction indicated in Fig. 1(b). The energy deposition applied to the powder corresponds to the calculated energy density at plane AA in the sample.
Within the transient simulation energy is deposited in the solid particle phase within the 1 μs beam pulse and then the resulting heat transfer to the helium and corresponding pressure spike and gas expansion can be observed. Figure 3 shows the calculated temperature of the tungsten powder and the helium gas immediately after the beam pulse of 3 × 10 12 protons, where it can be seen that the helium gas has almost reached the temperature of the tungsten by the end of the beam pulse. By 1 ms after the beam pulse the helium temperature is equal to the tungsten temperature. Figure 4 shows the simulated tungsten powder volume fraction at several times after the trigger of the beam pulse. It can be seen that immediately after the beam pulse (time ¼ 1 μs) no movement of the powder has occurred. After 1 ms the powder is starting to lift, reaching a height of 2 cm by 8 ms. The high pressure pulse results in a significant eruption, leaving a trench in the powder where the beam has passed through. Following the first HiRadMat experiment the samples were seen to have a trench in the powder along the beam axis [5]. Figure 5 shows the calculated helium velocity at various points during the lift. The highest velocity occurs soon after the beam pulse near the surface of the powder, then as the gas expands the pressure gradient and velocity decrease as expected. Figure 6 shows a summary of a series of simulation runs with varying particle diameter, proton pulse intensity and helium pressure. The main points to note are that smaller particles will lift at a lower proton pulse intensity, and that lifting will occur at a lower pulse intensity with atmospheric helium as compared to helium at 20 mbar counter to observed behavior. At this point it is useful to highlight the fact that at 20 mbar, 14 micron diameter particles were observed to lift by 30 mm within 20 ms at an intensity of 2 × 10 11 protons. The CFD simulations suggest we cannot attribute a lift of this magnitude to aerodynamic effects alone until the pulse intensity reaches at least an order of magnitude more than this. Comparison of these simulations with experiment has consistently led us to believe that while this is indeed a real physical phenomenon, powder lifting has been observed at much lower proton beam intensities than seem to be required for aerodynamic lift. This suggests that another more prominent mechanism is needed to explain the observations.

III. THERMAL EXPANSION
We again consider that the proton beam heating results in a sudden increase in temperature of the tungsten particles. The particles then expand corresponding to their coefficient of thermal expansion α and this occurs at a timescale proportional to the inverse of the speed of sound in the solid. This rapid expansion can result in the propagation of forces between particles which are often referred to as "force chains." We now consider how fast a solid cylindrical particle of height L will lift from a rigid surface if it is instantaneously heated. The timescale for expansion of the solid is simply given by where c is the speed of sound in the material. Near instantaneous uniform heating results in a uniform stress state given by the following expression: where E is Young's modulus and ΔT is the change of temperature. The net force acting on the solid (with density ρ) can be written as where A is the cross-sectional area of the cylinder in contact with the rigid surface and the second term represents the downward force on the body due to gravity (with acceleration constant g). Rearranging leads to the following expression for the vertical acceleration of a particle against gravity: FIG. 5. Helium velocities at various times after a beam pulse of 3 × 10 12 protons, a ¼ 1 μs, b ¼ 1 ms, and c ¼ 4 ms.
FIG. 6. Predicted particle lift due to pressure pulse and aerodynamic effects.
The initial upwards velocity of the particle is therefore Typically for a temperature rise of interest the second gravity term is relatively small. This leads to the interesting observation that the initial lift velocity is independent of the particle size, suggesting that any particle will jump at a velocity that depends on its temperature rise and material properties only. Figure 7 shows the initial lift velocity as a function of the instantaneous temperature jump in tungsten, where the slight departure from a straight line is due to the temperature dependent properties of tungsten. This simple expression was found to agree well with transient finite element analysis of an expanding solid on a rigid surface. We now consider a simplified model of a column of six 1 mm diameter rigid particles all resting perfectly on top of each other within the powder sample (Fig. 8). Particle 1 is on the beam axis and particle 6 is at the surface of the sample; recall that the beam axis is 6 mm below the surface of the sample.
We can calculate the initial lift velocity of each particle in the stack as a function of the Gaussian temperature rise profile. Then, assuming the column of particles is resting on a rigid surface and each particle is in direct contact with its neighbor, the velocity of the top particle (numbered 6) can be determined from a cumulative sum of the initial velocities of all the particles below it; i.e. (V 01 þ V 02 þ V 03 þ V 04 þ V 05 þ V 06 ). Note that since lift velocity is independent of particle size the lift velocity of the top particle would be significantly higher for an equivalent vertical column of much smaller particles, assuming theoretically perfect contact between neighboring particles. For the most intense beam pulse of the experiment (3 × 10 12 protons), the peak temperature rise on the beam axis is expected to be around 1500 K. Figure 9 shows the assumed temperature of each of the six particles in the column and the associated velocity, and also the cumulative velocity assuming they are all directly coupled. The vertical velocity of particle 6 at the surface of the sample is estimated to be 50 m=s. In fact, during the experiment there was no perceptible lift of the 1 mm spherical particles observed. The main reason for this is thought to be that the particles do not directly rest on each other so that there is space for a lot of the thermal expansion to occur without resulting in the propagation of vertical forces. The length scale for the attenuation of force chains in granular media has been predicted using discrete element modeling [13], where the attenuation becomes more significant as the level of disorder in a granular media increases. It was also reported that high frequency perturbations are dissipated in the order of ten to twenty particle diameters.
This suggests that force chains are less likely to propagate from the beam axis to the free surface in the case of small particles as compared to the large particles tested. In addition to this work in the literature we carried out a test using a piezocrystal and a laser Doppler velocimeter (LDV) to measure the propagation of a perturbation through tungsten powder containing particles of two different diameters of 3 mm and 150 μm. The Noliac piezocrystal expanded by 1 μm in the order of 1 μs and so could cause an initial perturbation velocity of the order 1 m=s. The piezocrystal was immersed in tungsten powder which was contained within an open topped metal container. The metal container had a hole on one side which was covered by a foil. The foil was in contact with the tungsten powder and was facing the piezocrystal surface with 8 mm of tungsten powder between the crystal and the foil. A 1.5 MHz LDV sensor was trained on the outside surface of the foil. Figure 10 shows the response measured with the LDV following a sudden expansion of the piezocrystal. It can be clearly seen that the larger 3 mm diameter particles are more efficient in propagating the perturbation to the foil window than the smaller particles. The fact that in the HiRadMat experiment the large (1 mm diameter) particles tested did not lift at all, and the smallest particles tested lifted the most, is contradictory to the results obtained with the piezocrystal. Combining this with the evidence from the literature regarding the dissipation of perturbations in granular media leads us to believe that thermal expansion and propagating force chains is not the dominant mechanism driving the observed powder eruptions.

IV. CHARGE EFFECTS
In order to consider the charge deposited in the tungsten powder sample by the incident proton beam we again use the FLUKA Monte Carlo code to simulate the collision physics in our simple model geometry shown in Fig. 1. Each beam pulse consists of a series of bunches that each travel through the tungsten powder within 0.3 ns. Many of the protons interact with the tungsten atoms resulting in many particles, including protons and electrons, being displaced from their starting positions. We use the Monte Carlo code to track the primary protons and secondary particles until their energy falls to 100 keV, where they are stopped and their kinetic energy is deposited as heat within the computational cell. Approximately 95% of the energy deposited in the tungsten by the beam results from secondary energetic electrons stopping in the tungsten.
The net effect of the deposited charged particles is a spatial charge distribution left within the tungsten powder target. The distribution shows a characteristic double layer effect (Fig. 11) with a region of positive charge in the path of the intense part of the beam, a region of negative charge in an annulus around the core of the beam where significant numbers of electrons are freed and subsequently deposited, and finally a positive region on the surface of the tungsten sample where energetic secondary electrons are able to escape the tungsten atoms leaving a positively charged layer (Fig. 12). The thickness of this positively charged layer depends on the stopping distance of the electrons produced near the surface of the sample. Figure 13 shows a higher resolution plot of the surface of the sample indicating that this layer is between 0.01 and 0.02 cm in thickness.
The Monte Carlo simulation gives us an insight into the residual charge pattern immediately after the beam interaction. However, following this period, the equalization of charge as a result of current flow in the sample is not modeled as part of the Monte Carlo simulation.
Very little charge is deposited in the atoms of the helium carrier gas. The majority of the residual charge appears in the tungsten atoms as expected. Any residual net positive or negative charge will travel to the surface of the tungsten particle in which it has been deposited in a very short time compared to the beam pulse length (in the order of 10 −19 s) [14], so by the end of the beam pulse all residual charge can be considered to be on the surface of the tungsten particles. The charge resides on the surface of the particles since Coulomb's law demands that the charges on a conductor be as far apart from each other as possible. The rate of neutralization of this charge by current flow between the metallic particles will depend on the electrical conduction between particle contacts and/or electrical breakdown across tiny gaps between particles. There is little published work on the electrical conductivity of an uncompressed metallic powder, which we suggest is due to the fact that the conductivity tends towards zero and is therefore hard to measure. Many researchers have found that the conductivity depends on compression force and some of them refer to a minimum compression pressure required for the presence of any conductivity [15,16]. Electrical breakdown between the tungsten particles depends on the potential gradient, the pressure of the helium and the size of the particles and gaps between them. We note that the maximum potential gradient occurring after the beam pulse is proportional to the proton intensity. Paschen's law states the breakdown voltage of two electrodes, as a function of the product of the gas pressure p and the gap width b between them, is given by where A ¼ σ n kT n and B ¼ AV i , η is the probability of secondary electrons being emitted, σ n is the electronneutral collision cross section, V i is the ionization potential, k is Boltzmann's constant and T n is the temperature of the neutral atoms [17]. Beam pulses were injected into the sample with helium at 1 bar and at 20 mbar. Figure 14 shows Paschen curves for helium at 1 bar and at 20 mbar. For a given gap size, if the voltage is above the curve then breakdown is expected to occur. In Fig. 15 we plot the maximum sustainable potential gradient before breakdown occurs. We also plot the peak potential gradient (electric field) due to the deposited charge determined from the solution of Poisson's equation (see Fig. 16).
It is thought that within the sample the gap size between particles will be no greater than the tungsten particle diameter, and so we expect the gap size in our samples to be always less than 1 mm which was the largest particle diameter tested. Inspection of Figs. 14 and 15 indicates that breakdown between particles with gaps smaller than 1 mm is almost impossible in 20 mbar helium. By contrast breakdown of the applied voltage gradient in 1 bar helium is expected down to a gap size of around 15 μm.
In order to simulate the mechanical effect of deposited charge on the tungsten particles we again employ the transient two-phase fluid mechanics model. By way of isolating the charge effect we now assume that there is no heating of the tungsten or helium phases and the only driver for a powder lift is the charge mechanism. We still solve the continuity and momentum equations for the two phases (gas and dispersed solid) but it is now an isothermal model so we no longer solve the energy equations. We now couple the Poisson equation with the simulation to determine the electric potential field that results from a given charge distribution; i.e.
where ε is the permittivity, φ is the potential and q is the charge density applied to the solid phase. The initial condition for the charge density is based on the FLUKA simulation results of residual net charge (Fig. 12) and thus we are assuming that there is no electrical breakdown between particles and the beam induced charge remains on the particles throughout this mechanical response simulation. The electrical field is simply determined from the gradient of the potential as follows: An additional momentum source is added to the solid phase momentum Eq. (4) and represents the Coulombic force on the tungsten particles, As the particles move the initial charge density disperses and so the electrical potential and field is recalculated as the transient simulation evolves. We choose a similar computational domain as was used for the aerodynamic lift simulation, i.e. a 2 mm thick slice of the tungsten sample normal to the beam direction at the plane AA. As before the domain has an open top and the walls are given a boundary condition of zero potential, i.e. grounded. Figure 16(a) shows the initial charge deposited on to the solid particle phase obtained from the FLUKA analysis. Figure 16(b) shows the electric field that results from the initial charge distribution. Figure 17 shows the particle volume fraction at various times following a beam pulse containing 3 × 10 12 protons incident on a sample containing 34 μm diameter tungsten particles. The predicted violent powder eruption can be clearly seen with some of the particles reaching a height of 30 mm in just 5 ms. It is interesting to note that there appears to be an eruption along the center line of the sample container but also significant amounts of powder are traveling up the walls of the container. The central eruption appears to occur in stages, with at least three distinguishable particle fronts. The first of those fronts is driven by the surface charge layer and depends on resolving the thin layer of positive charge near the surface of the powder which is a challenge from a computational perspective. To completely resolve the very steep increase in deposited charge (Fig. 13) over a narrow region of the order 0.2 mm requires a large number of small elements and prohibitively small time steps. As such there is some averaging effect of this surface charge in the simulation meaning that the peak maximum velocity of the first front could be higher than what is seen in the simulation. Having said that averaging of the surface charge layer does not affect the magnitude of the main powder eruption seen in the second and third particle fronts, which depends on the charge deposited on and around the beam axis. The second front is seen to rise by at least 20 mm in 5 ms and the particles traveling up the wall reach 25 mm in 5 ms.
If we attempt to track the leading edge of the main particle eruption in the charge induced lift simulations we get a profile similar to a true ballistic profile. Figure 18 shows the data from the simulation compared with a ballistic profile where the height depends only on the initial impulse velocity as follows: The ballistic profile requires a very large initial acceleration giving rise to the initial velocity followed by negligible acceleration. The charge induced simulations have an initial large acceleration when the residual charge has been applied but the powder remains unperturbed. Then as soon as the charged particles move apart the force and corresponding acceleration reduces rapidly. It should be noted that we would expect some error when tracking the particle front in simulations and so conclusions about deviations from a ballistic profile should be made with caution. However, it is interesting to note that after the initial acceleration some reasonably significant Coulombic forces are still present. This can be seen in Fig. 19 which shows the maximum vertical powder acceleration in the computational domain following a beam pulse. Deviations from ballistic behavior were also noted in the previously published observations [6]. Figure 20 shows the predicted lift height reached 20 ms after the beam pulse for a range of proton intensities and it can be seen that a lift is still observed with a pulse intensity of 3 × 10 10 protons. Powder lifts were observed at this lower intensity in the experiment.
If we compare the maximum (initial) lift velocity predicted by this charge model with experiment we find reasonable agreement as shown in Fig. 21. Also shown for comparison is the predicted lift due to aerodynamic forces (dotted line). We hypothesize that there is a maximum surface charge density on the tungsten particles that is possible before electrical breakdown will occur and charge will jump across to other particles with a lower potential. This limit on surface charge density means that a powder sample consisting of smaller tungsten particles can carry greater charge per unit volume and thus experience greater Coulombic eruptions than a powder consisting of larger particles. For example, if we say the maximum surface charge density possible on a spherical particle due to breakdown phenomena is 27 μC=m 2 [18] we find that the predicted charge-induced lift is now dependent on particle size and that 1 mm diameter particles, which have much less surface area and have a limited charge carrying capacity, are not predicted to lift at all compared to 100 micron diameter particles or smaller. This behavior was seen in the second experiment where a range of particle sizes was exposed to the beam pulse. Figure 22 shows a high speed video image of the multisized trough (looking in the x direction) after irradiation with 2 × 10 11 protons on target. The 14 micron average diameter particles produced the greatest response, followed by smaller responses for the larger particle sizes. The 1 mm particles did not appear to lift at all. The experiment unequivocally revealed that smaller powder spheres responded to the beam with a higher lift and velocity than larger spheres.

V. DISCUSSION
We have considered three main physical mechanisms to explain the observed proton beam induced eruptions of stationary tungsten powder samples. Thermal expansion of particles has been ruled out as it is believed that the particles are not in a direct line of contact and so much of the thermal expansion of each particle occurs without pushing on adjacent particles. It can be seen in the literature that perturbations in granular media do not propagate very far due to them relying on force chains which do not transmit force very well in an uncompacted powder sample as we have. Additionally, we have shown that perturbations can be transmitted through large particles more effectively than through small particles and yet it is the small particles that have shown the largest eruptions.
Aerodynamic effects have been considered in some detail and it is believed that with the higher intensity pulses it is possible to induce a powder lift from the resulting pressure pulse and gas expansion. However, our observations often showed powder lifts occurring at much lower pulse intensities than would be required for an aerodynamic force induced lift.
Predictions of charge-induced lift agree well with observations. The theory is that the deposited charge remains on the surface of the particles for long enough to result in significant Coulombic repulsive forces lifting the particles.
The experiments and simulations were designed to isolate and identify the relative importance of the aerodynamic and charge effects. As helium pressure is reduced the magnitude of aerodynamic effects has been shown to reduce and the chance of electrical breakdown is less likely and so charge effects are likely to become dominant. At high helium pressure the cause of powder eruptions is less clear as aerodynamic effects increase with pressure and electrical breakdown is likely so charge induced lift may be diminished. Figure 23 shows a marked difference in lift between experiments with similar proton pulse intensities but different helium pressures. The response in helium at 20 mbar was found to be significantly greater owing to less electrical breakdown occurring resulting in a significant Coulombic eruption. The response in atmospheric helium was much less, which could be attributed to some electrical breakdown occurring between the nonspherical crystalline particles within the helium environment thus reducing the magnitude of Coulombic forces.
The hypothesized maximum surface charge density on a particle before electrical breakdown provides a credible explanation as to why small particles with much higher surface area per unit volume can carry much more charge and undergo a more violent response than the larger particles.

VI. CONCLUSIONS
A stationary tungsten powder sample has been observed to erupt when exposed to a high energy proton beam pulse, with larger eruptions observed for smaller powder particles and lower helium pressure. Simulations show that this can be explained by the charge deposited on the tungsten particles by the incident beam, where the charge remains on the surface of the particles long enough to result in significant lifting Coulombic repulsive forces. Thermal expansion has little effect, owing to the gaps between neighboring powder particles limiting the propagation of force chains, while aerodynamic drag requires much higher proton beam intensities than was experimentally observed, although it may contribute to lift at higher helium gas pressures (1 bar abs) where electrical breakdown is likely.

VII. NOMENCLATURE
A area a acceleration α thermal expansion coefficient b gap width c speed of sound c αβ interphase drag parameter d p particle diameter δ ij Kronecker delta E modulus of elasticity E electric field vector ε permittivity of free space F force g gravitational acceleration γ thermal conductivity h height h enthalpy k Stefan Boltzmann constant M momentum source η probability of secondary electrons being emitted r volume fraction ρ density p pressure P s solid pressure gradient φ potential q charge density Q interphase heat transfer Re Reynolds number S heat source σ n electron neutral collision cross-section t time T temperature τ solid phase stress term U velocity vector V B breakdown voltage Vi ionization potential V 0 Initial velocity μ dynamic viscosity Subscripts α refers to continuous gas phase β refers to dispersed solid phase

ACKNOWLEDGMENTS
The authors thank the CERN HiRadMat and SPS operation teams for the technical support during the experiment. This research has received partial funding from the European Commission and Science and Technology Facilities Council (STFC) under the FP7 Capacities project EuCARD2, Grant Agreement No. 312453. We also thank Dr. Lei Yang, Director of the Laboratory for Spallation Targets at IMPCAS, and Dr. Jiangfeng Wan, Leader of The Granular Target Design group, for supplying some of the tungsten material used for the experiments.