Using a direct simulation Monte Carlo approach to model collisions in a buffer gas cell

A direct simulation Monte Carlo (DSMC) method is applied to model collisions between He buffer gas atoms and ammonia molecules within a buffer gas cell. State-to-state cross sections, calculated as a function of the collision energy, enable the inelastic collisions between He and NH 3 to be considered explicitly. The inclusion of rotational-state-changing collisions affects the translational temperature of the beam, indicating that elastic and inelastic processes should not be considered in isolation. The properties of the cold molecular beam exiting the cell are examined as a function of the cell parameters and operating conditions; the rotational and translational energy distributions are in accord with experimental measurements. The DSMC calculations show that thermalisation occurs well within the typical 10-20 mm length of many buffer gas cells, suggesting that shorter cells could be employed in many instances—yielding a higher ﬂux of cold molecules. ' 2017 Author(s). All article content, except where otherwise noted, is licensed under a Creative Commons Attribution (CC BY) license


I. INTRODUCTION
In contrast to techniques such as laser cooling, cold buffer gas atoms can be employed to translationally cool atomic and molecular species almost independently of the internal energy structure of the species. As such, buffer gas cooling methods are very versatile, requiring only that the species of interest be generated in the gas phase and that the formation of clusters is hindered. 1,2 With a range of techniques developed to produce gas-phase species in the buffer gas cell-including capillary filling, laser ablation, beam injection, and discharge etchingbuffer gas cooling has been successfully applied to numerous atomic and molecular systems. 3,4 Within a buffer gas cell, the molecules of interest undergo elastic and inelastic collisions with a cold buffer gas such as helium. After multiple collisions, the molecular species are internally and translationally cooled. The buffer gas cell is typically kept at a density which achieves sufficient collisions for thermalisation, while preventing three-body collisions and cluster formation. To establish the optimal buffer gas cell conditions-such as the ideal cell dimensions to achieve maximal flux of cold molecules out of the cell-a methodology for simulating the rotational and translational cooling processes is required.
Gas dynamics can be described using a macroscopic or microscopic approach. The Navier-Stokes equations can be employed to provide a macroscopic description of a system, assuming a continuous medium. Alternatively, the Boltzmann equation can provide a microscopic description of a system on the single particle scale. As macroscopic properties represent averages of microscopic quantities (provided there are sufficient numbers of molecules in the system), both approaches should yield the same result. However, this is often not the case, as there are limitations with both treatments. For example, the conservation equations within the macroscopic model frequently do not form a determinate set (i.e., there is often no unique solution for a given set of macroscopic quantities). If taking the microscopic approach, the calculation of pair-wise interactions scales quadratically with the number of particles and thus becomes unfeasible for all but the smallest of systems. To overcome these difficulties, we adopt a modified microscopic approach.
In this paper, we simulate the collisions between helium buffer gas atoms and ammonia molecules within a buffer gas cell, yielding a rotationally and translationally cold molecular beam. To characterise the buffer gas cooling process, a direct simulation Monte Carlo (DSMC) approach, also known as Bird's method, is adopted. 5 The DSMC algorithm is widely utilised in systems where it is desirable to track the collisions of particles, as is the case in this work. We have adapted the methodology to enable inelastic collisions between ammonia and helium to be considered explicitly, employing state-to-state cross sections calculated at a fine grid of collision energies. This explicit treatment of inelastic collisions negates the need for models such as that proposed by Larsen and Borgnakke, 6 which require approximate inelastic collision probabilities and relaxation rates and assume a statistical distribution of energy governed by a single rotational temperature. While extensions of the Larsen-Borgnakke model have been successfully implemented, 7-9 these too have limitations in gas mixtures of non-uniform density and under strongly non-equilibrium conditions. The properties of ammonia molecules that exit the cell are examined as a function of physical parameters such as buffer gas cell length and exit aperture position, enabling optimal parameters and operating conditions to be identified. The simulated rotational and translational temperatures are in good agreement with experimental measurements, verifying the accuracy of the DSMC approach.

A. DSMC model
To completely describe the behaviour of an ensemble of molecules and atoms in the gas phase, the trajectories of each species must be explicitly calculated using Newton's equations of motion, considering the intermolecular forces at play. As the Newtonian trajectories are dependent on the position of a given particle and that of its neighbouring species, a complete model becomes intractable for systems comprising large numbers of particles. The DSMC method offers an alternative to the (inaccessible) exact deterministic description of a system: collisions are generated stochastically, with the collision rates and post-collision velocities derived from the kinetic theory of gases.
There are several assumptions and approximations inherent in the DSMC method. For example, each simulated particle represents a large number of physical particles. While this can lead to statistical scatter, it is only problematic in situations where there are very low densities and thus the simulated particles (which each represent 10 11 physical He or NH 3 particles) cannot adequately represent the physical conditions. As detailed in a previous publication, 10 the DSMC method can successfully describe a supersonic beam expansion. A similar approach is adopted here, where the system of interest-a buffer gas cell-is divided into subcell regions of finite size. Each cell dimension is broken into 30 subcells. Discrete time steps are employed, over which molecular motion and collisions are uncoupled. The simulation propagates the positions and velocities of particles with a step size of 1 ns, chosen to be less than the mean collision time. As the time step and cell size decrease, the DSMC treatment becomes more exact. The dimensions of each subcell are chosen such that the subcell size is on the order of the mean free path of the particles considered.

B. Rotational state-changing collision cross sections
As detailed in two recent publications, 10,11 quantum closecoupling calculations have been performed to establish inelastic cross sections for the collisions of NH 3 with He at a range of collision energies, spanning 0.01-150 cm 1 . The cross section calculations are carried out with the close-coupling method using the four-dimensional potential energy surface (PES) developed by Gubbels et al. 12 The four dimensions of the PES comprise the three spherical coordinates of He relative to NH 3 and the NH 3 umbrella (or inversion) angle. Internal ammonia bond lengths are held constant at the average ground vibrational state values. The inversion motion of NH 3 is explicitly considered in order to obtain accurate cross sections at low collision energies, as explained in Ref. 11. The surface is fit to a grid of points calculated ab initio at the CCSD(T)/aug-cc-pVQZ level of theory.
The rotational states of NH 3 are defined by the labels |J, K, ± , where J is the total angular momentum and K its projection on the NH 3 symmetry axis, with the + and − signs referring to the symmetric and antisymmetric states of the inversion doublet. The inversion splitting in NH 3 is 0.79 cm 1 . 13 Ammonia has two nuclear spin configurations; ortho-NH 3 states have A 2 symmetry, while para-NH 3 states have E symmetry. There is no interconversion between states with different symmetries (A 2 / ↔ E); hence the cross sections for symmetry-changing transitions are zero.
Considering the degeneracies of each state and assuming a thermal equilibrium, the relative number of NH 3 molecules in a given rotational state is given by where g J,K is the statistical weighting, E(J, K) is the energy of the state, k B is the Boltzmann constant, and T is the temperature. Experimentally, ammonia enters the buffer gas cell at T ≈ 210 K. Under these conditions, the probability that J ≥ 9 is 3%. Accordingly, the state-to-state cross sections are obtained by solving the close-coupling second order differential equations, considering all symmetry-and energy-permitted transitions for initial rotational states up to and including J = 8. The convergence properties of the scattering calculations are described elsewhere. 10,11

C. NH 3 -buffer gas collisions
The probability distribution of the collision energy between 6 K He and 210 K NH 3 is provided in Figure 1 is calculated using the Boltzmann distribution and considering the relative velocity and temperature-weighted centre-ofmass velocity of the colliding particles. While the distribution in Figure 1 is only applicable to the initial collisions, as it applies to an unphysical homogeneous mixture of molecular gas and buffer gas at two very different temperatures, it provides an indicative experimental collision energy range of 0-150 cm 1 . As such, the collisional energy transfer between all symmetry-and energy-allowed transitions is considered at collision energies from 0.1 cm 1 to 150 cm 1 , in accordance with the experimental parameters.
As each simulated NH 3 particle represents a large number of physical molecules, it is assigned a rotational state distribution, defined by a normalised vector Π. The rotational distribution is again limited to J < 9; given the initial NH 3 rotational population distribution and the low collision energy, the de-excitation of rotational population to lower energy levels is the dominant outcome of a state-changing collision. A Markov chain description is used to evolve the rotational state distribution after each collision. A matrix P establishes the transition probabilities, with the final rotational population distribution calculated from the product of the initial rotational population distribution vector and the probability matrix, Π f = Π i P. Each element in P is calculated from the energy-dependent ratio of the inelastic cross section for a given collision-induced transition and the total collision cross section, P ij (E) = σ ij (E)/σ tot (E). The principle of detailed balance is satisfied, and the total energy is conserved at each step, as described in previous work. 10 Collisions are only considered between pairs of particles within the same simulated subcell. The probability P coll of a collision between two simulated particles (which each represent F N real particles) during time interval ∆t is equal to the ratio of the volume swept out by the total collision cross section at the relative velocity of the collision, divided by the volume of the subcell (V cell ), As such, P coll is proportional to c r , the relative speed of the pair of particles. When a collision is deemed to have occurred between a NH 3 molecule and a He buffer gas atom, the collision energy is calculated and the appropriate rotationally inelastic cross section is selected. 10 Post-collision velocities are randomly oriented, with the magnitude of the velocity determined by the energy transferred in the collision.
Calculating collisions between all M pairs of simulated particles in a subcell is inefficient, given the frequently large number of simulated particles and the low collision probability for any given pair. As such, only M possible collision pairs are explicitly considered in a given time step. Bird's modified No Time Counter (NTC) method 14 was designed to aid the selection of representative collision pairs, where σ tot c r in Equation (2) is replaced by (σ tot c r ) max to yield the fraction of collision pairs explicitly considered, In Equation (3), (σ tot c r ) max is a simulation parameter, set to an initially low value and adjusted whenever a higher value of (σ tot c r ) is encountered. (Note that the value of (σ tot c r ) max is not critical; as Bird explains, the term later cancels out of the calculation. 5 ) Collisions thus occur between the selected pairs with probability P coll = (σ tot c r )/(σ tot c r ) max . In the DSMC simulation, M is calculated using Equation (3), and potential collision partners are selected using a random number generator. If P coll calculated for a given collision pair is higher than a random number between 0 and 1 (R), the collision is accepted; if P coll < R, the collision is rejected. The overall collision probability in a subcell is subsequently calculated by weighting the sampled collision probability by the total number of simulated particles in the subcell. Wall surfaces are modelled thermally: collisions between molecules and the walls are treated elastically. Post-collision velocities are determined from the magnitude of the particle's pre-collision velocity, the collision angle, and the wall temperature, with linear momentum conserved. It is assumed that collisions between pairs of NH 3 molecules are purely elastic, with the hard-sphere collision model employed for collisions between pairs of He atoms. The validity of these assumptions has been previously established. 10

D. Modelling an experimental buffer gas cell
The simulations model an experimental buffer gas cell, illustrated in Figure 2 and based on the design of Sommer et al. 15 As described in Ref. 16, the 20 × 40 × 40 mm (length × height × width) buffer gas cell in the experimental setup is attached to the second stage of a two stage pulse-tube cryocooler.The buffer gas inlet pressure is typically maintained at 0.6 mbar, with an ammonia inlet flow rate of 1.0 SCCM. Helium enters the buffer gas cell at 6 K, as the buffer gas line thermalises with each cryocooler stage prior to entering the buffer gas cell. As the molecular gas line is (largely) thermally disconnected from the cryogenic environment, the temperature at the inlet can be controlled by a small heater block. A 3 mm section of the molecular gas line extends into the buffer gas cell (inner diameter 1 mm, outer diameter 2 mm) FIG. 2. Schematic illustration of the experimental buffer gas cell. Ammonia molecules (indicated as red spheres) enter the front face of the cell through the molecule inlet; helium buffer gas atoms (depicted as blue spheres) can be seen entering the cell from below. NH 3 molecules that exit the cell with sufficiently low velocities are electrostatically guided into a differentially pumped reaction chamber. and is simulated as a cylinder. Ammonia molecules are simulated with a Boltzmann distribution of velocities centred on the mean velocity (as is typical for the treatment of a moving gas at equilibrium). 5 The inlet buffer gas line and the exit aperture have an inner diameter of 1 mm. While the exit aperture is modelled in 2D, the simulated particles complete their trajectories during each discrete time step. As such, there is some propagation of trajectories through the exit aperture and thus a spatial distribution in three dimensions at the outlet.
Experimentally, the buffer-gas-cooled ammonia molecules effuse out of the exit aperture into a region under high vacuum, where they enter an electrostatic quadrupole guide positioned 1 mm from the cell exit. 16 The quadrupole is assembled from hand-polished stainless steel rods with a 2 mm diameter circular cross section. Voltages of up to ±5 kV are applied to the guide electrodes, achieving maximal field strengths of up to 90 kV/cm at the electrode surfaces. Translationally cold molecules in low-field seeking states are guided around three 90 • bends of 100 mm curvature into a differentially pumped reaction chamber, where the rotational state population is probed by resonance-enhanced multiphoton ionisation (REMPI). The guide is not included in the simulations but is sufficiently far from the exit aperture not to affect the gas flow out of the cell.

A. Rotational and translational energy distributions
For each simulated particle leaving the cell through the exit aperture, the molecular properties-such as the rotational state and velocity-are recorded. A range of cell lengths is considered, with the goal of establishing the minimum length required for full equilibration of the internal and external degrees of freedom-yielding the maximal flux of cold molecules out of the cell. As anticipated, the rotational and translational temperatures of the NH 3 molecules initially decrease with increasing cell length. Convergence of the molecular parameters is observed with cells of length approximately 3 mm beyond the molecular inlet, well before the position of the experimental exit aperture (located 17 mm beyond the molecular inlet tube). Given this effective equilibration of molecular properties, and the computational effort required, longer cells have not been simulated. It is worth noting that the equilibration of rotational and translational temperatures occurs at roughly similar rates, implying that the energy exchange per collision is approximately equal for translation and for rotation (bearing in mind that these two are not wholly separable in our treatment).
The translational temperature obtained from fitting the simulated velocity distribution with an effusive beam velocity distribution function is plotted against the cell length in Figure 3. (For the purpose of this section and the figures appearing therein, "cell length" refers to the distance between the end of the molecular inlet tube and the exit aperture.) For an effusive gas flow, the exit aperture thickness and diameter are both much smaller than the mean free path of the gas (under steady-state conditions). Thus the velocity distribution of the beam exiting the buffer gas cell is given by wherev is the mean thermal velocity. In this way, the mean forward velocity and full width at half maximum of the characteristic angular spread can be established. 4 The simulated translational temperature of NH 3 equilibrates to 5 ± 1 K at cell lengths of 2 mm and longer. This is in agreement with the experimental buffer-gas-cooled ND 3 translational temperature of 5.96 ± 0.05 K, measured at the end of the quadrupole guide, 16 and with the 6.0 K He buffer gas. The uncertainty in the NH 3 translational temperature primarily arises from the uncertainty in the fit of the velocity distribution to the simulated particles. Statistical scatter may also influence the simulated translational temperature, as only a fraction of simulated molecules reach the exit aperture. At cell lengths ≤1 mm, a bimodal rotational temperature distribution is obtained. The presence of two distinct temperature regimes is attributed to the rotational state dependence of the cross sections: high-J states typically require more collisions for de-excitation than low-J states, arising from the reduced inelastic cross sections. This is consistent with the formula proposed by Polanyi and Woodall some 44 years ago, where they describe the probability of de-excitation based on the energy difference between the two states, ∆E, where C is a constant. 17 As the energy difference between consecutive rotational states increases with increasing J, the transition probability decreases. It is not possible to assign a single rotational temperature to the bimodal distribution observed at short cell lengths. Accordingly, only the "cold" part of the distribution is considered when fitting a rotational temperature for cell lengths ≤1 mm. Thus the temperatures plotted in Figure 4 for short cell lengths represent the lower limit of the NH 3 rotational temperature. For cells longer than 1 mm, there is only one component to the distribution and thus the entire population is FIG. 4. Rotational temperature of the NH 3 beam, at a range of simulated buffer gas cell lengths, ascertained from a Boltzmann fit to the simulated rotational state distributions. As states with A and E parities do not interconvert, these are considered (and plotted) separately. The shaded region indicates the uncertainty in the temperatures derived from the Boltzmann distribution fit, with the solid lines illustrating the trend in rotational temperature. considered when assigning a rotational temperature. The NH 3 rotational temperature equilibrates to a constant value (within the reported error) beyond cell lengths of ∼2.0 mm. Plots of the simulated population in each J state, provided in Figure 5 and Table I, confirm that the rotational population distribution does not change (within the uncertainty of the simulations) for cell lengths beyond 2.0 mm.
The equilibrated rotational temperature of buffer-gascooled NH 3 molecules established from these simulations, as presented in Figure 4, is 6.2 K. As states with A and E symmetry do not interconvert, the rate of rotational cooling depends on the relative collision cross sections of the states within each symmetry group. The experimentally measured rotational temperature of buffer-gas-cooled and electrostatically-guided ND 3 under comparable conditions is 10 K, as ascertained from PGOPHER 18 fits to the experimental REMPI spectra. 16 The influence of alterations in the buffer gas and ammonia densities has been experimentally examined 16 and is also explored with the simulations in this work. When the experimental flow rate of ammonia is increased from 1.0 SCCM to 3.0 SCCM, with the inlet pressure of helium held constant at 0.6 mbar, the rotational temperature of the buffer-gas-cooled and electrostatically guided beam increases by (35 ± 5)%, from 10 K to 13.5 K. A similar trend occurs in the simulations: as the density of NH 3 is increased in the same manner (with He held constant), the rotational temperature of the beam at the exit aperture increases by (21 ± 12)% or 1.2 ± 0.7 K in the 3 mm cell.
The rotational (and translational) temperatures obtained from the DSMC simulations are consistently lower than the experimental measurements, although the experimentally observed temperatures do fall within the uncertainty of the DSMC values. This could be due to a number of effects. First, it is not possible to experimentally probe the buffer-gas-cooled beam as it exits the cell (in our set-up). Only molecules in lowfield seeking states are transmitted through the guide, and thus no molecules in J = 0 can be directly detected. (It should be noted that in assigning a rotational temperature to the experimental REMPI spectra, the probability that molecules in a given state are transmitted through the guide is accounted for.) Furthermore, any imperfections in the 2 m three-bend electrostatic guide (such as misalignment of adjacent sections) could hinder the transmittance of, for example, the slowest-moving molecules. Black body radiation (BBR)-induced transitions could also increase the rotational temperature of molecules as they travel through the guide. The electrostatic guide is disconnected from the cryogenic environment of the buffer gas cell, and thus ammonia molecules will be exposed to ambient (298 K) BBR while they traverse the 2 m length of the guide, which takes approximately 2.6 ms.
There will also likely be some subtle differences between the behaviour of NH 3 and ND 3 , owing to the different rotational constants and inelastic cross sections. For example, more rotational levels are populated in ND 3 at 210 K, meaning that rotational states up to J = 12 or 13 must be included in the calculations. Collisions between ND 3 and He need to be explicitly modelled (or, alternatively, further experimental measurements undertaken with NH 3 and He) in future work to establish the differences in the behaviour of the ammonia isotopologues.
The DSMC simulations also make several assumptions about the buffer gas cell conditions. For example, molecules that diffuse to the cell walls before they are fully thermalised will freeze in the experiments; collisions between simulated particles and the wall are modelled thermally in the simulations. As such, there will likely be some experimental deviations from the assumed buffer gas cell conditions. Agreement within the uncertainty of the calculations, and reproduction of the experimental trends observed with changes in the buffer gas and molecular gas densities, serves to validate the approach taken in this work.
An important observation is the inseparability of rotational and translational energy. DSMC simulations performed both including and excluding inelastic collisions confirm that the final equilibrated translational temperature of the molecules is influenced by inelastic collisions-with the translational temperature thermalising to a value some 10% higher when rotational-state-changing collisions are included in the simulations. This result indicates that elastic and inelastic processes cannot be considered in isolation, as occurs in some approximate treatments of buffer gas cell behaviour.

B. Collisions in DSMC simulations
The number of collisions in the DSMC simulations before a given molecule reaches the exit aperture is orders of magnitude higher than that predicted by considering the mean free path of a molecule in the buffer gas cell and the shortest route to the exit aperture. Invoking the hard sphere assumption and taking a typical thermally averaged collision cross section of 10 14 cm 2 at an approximate buffer gas density (from our experimental cell design) of 5×10 15 atoms cm 3 yields a mean free path of ∼ 8 × 10 −5 m. 4 Hutzler et al. use this mean free path value to suggest that, as up to ∼100 collisions are required for translational and rotational "quenching", the thermalisation length for the molecular species in the buffer gas cell is not more than 1 cm. 4 The DSMC simulations indicate that, under our conditions, no molecules take a direct path to the exit aperture. As such, molecules undergo far more collisions than required for the thermalisation of the internal and external degrees of freedom by the time they reach the exit aperture. By the time ammonia molecules travel a net distance of 3 mm from the molecular inlet, they undergo an average of 10 6 collisions. While this is a very high number of collisions, it is consistent with the diffusion properties of gases. Gaseous molecules in atmospheric conditions, with velocities on the order of 10 2 m s 1 , exhibit diffusion velocities of around 10 2 m s 1 . This is because gaseous species take a circuitous path due to collisions, typically travelling 10 4 times the net distance between their starting and ending points. 19 Combining the approximate mean free path with the typical random walk properties suggests that molecules undergo roughly 3 × 10 5 collisions when travelling a net distance of 3 mm in the buffer gas cell-in reasonable agreement with the 10 6 collisions calculated by the DSMC simulations.
As set out in Equation (2), the probability of a collision between two simulated particles is proportional to ∆t: given a longer time, a larger volume will be swept out by the cross section and thus the collision probability will increase. At some point, this direct relationship breaks down; P coll (and P coll ) must be ≤1, whereas there are no restrictions on ∆t, although the time step should be set to a fraction of the local mean collision time. Figure 6 illustrates this breakdown, with the DSMC collision probability exceeding the collision probability calculated analytically for P coll ≥ 0.1.
Previous applications have found a dependence of certain physical properties on the size of the time step, for example, in DSMC treatments of ultra-thin gas film lubrication systems. 20 Ng et al. state that limiting the time interval to be a fraction of the mean collision time is insufficient, proposing a further requirement that ∆t < cell size mean velocity for DSMC simulations under such conditions. To ensure that we are operating within the region where the DSMC collision probability and that calculated analytically are in agreement, we select ∆t to be 1 ns. Simulations performed with ∆t = 0.5 and 10 ns yield the same results as achieved with ∆t = 1 ns. Thus, in this work, we find no dependence of physical properties on the size of the time step.
In sophisticated DSMC (DSMC07), the time step and cell size vary dynamically, with preferential selection of the nearest particle in the cell as the collision partner (unless those particles recently collided). 21 This is in contrast with the traditional DSMC (DSMC94) treatment-adopted in this work-where collision partners are selected at random within a cell, with ∆t and cell size kept constant throughout the simulation. Employment of the traditional DSMC approach avoids the underestimation of average collision numbers that can occur when there is correlation between ∆t, the number of particles, and the cell size 22 and reduces the sensitivity to ∆t and the number of particles within a cell. 21

C. Influence of buffer gas cell design
Experimentally, the exit aperture is located in the centre of the far wall of the buffer gas cell, opposite the molecular inlet line (see Figure 2). The influence of the position of the exit aperture on the properties of the molecules that exit the FIG. 6. The ratio between the DSMC collision probability and that calculated analytically is plotted as a function of P coll . It can be seen that the collision probability calculated using the DSMC approach deviates from the analytical result by more than 1% for values of P coll ≥ 0.01. cell is examined by comparing simulations featuring different exit aperture locations on the far wall. The position of the exit aperture has a negligible effect on the properties of ammonia molecules that exit the cell. (Similarly, there is no dependence of the molecular properties on the position of the buffer gas inlet line on the bottom wall.) These findings are attributed to the high number of collisions that molecules undergo within the cell, resulting in complete thermalisation well before any position on the far wall is reached.
The high number of collisions observed in the DSMC simulations suggests that the optimal buffer gas cell length might be shorter than the 1-2 cm selected in many current designs. As Figure 7 illustrates, the flux of ammonia molecules exiting the cell decreases exponentially as the cell length increases. From a 2 mm cell to a 3 mm cell, there is a 66% decrease in the molecular flux out of the cell. This behaviour is in accordance with the extraction model proposed by Hutzler et al., where the rate at which the buffer gas exits the cell,Ṅ b , is given byṄ where N b is the total number of buffer gas atoms in the cell, V cell is the cell volume, and A aperture is the area of the exit aperture. 4 Given the computational challenges associated with simulating longer cells, we are not in a position to quantitatively compare the relative molecular flux from a 3 mm cell with a 17 mm cell. However, the trend observed in Figure 7 and the model proposed by Hutzler et al. indicate that the flux will continue to decrease as the cell dimensions are increased. The buffer gas cell adopted in this work features a molecular inlet line that extends 3 mm into the cell. A cell length extending 3 or 4 mm beyond the molecular line inlet-for a total cell length of up to 7 mm-would serve to achieve higher molecular density at the exit aperture whilst still achieving the full thermalisation of NH 3 molecules. Furthermore, a shorter cell will limit diffusion to the cell walls. It should be noted that longer cells will still be required FIG. 7. The relative number of ammonia molecules that exit the buffer gas cell is shown as a function of the cell length. in circumstances where vibrational relaxation is required, as vibrational quenching cross sections are often several orders of magnitude lower than rotational state-changing cross sections. 23

IV. CONCLUSIONS
Collisions between He atoms and NH 3 molecules in a buffer gas cell are modelled using a DSMC method. The inclusion of state-to-state cross sections enable inelastic collisions to be explicitly considered. This allows the rotational state population distribution to be tracked and also facilitates an accurate calculation of the translational temperature of the beam. The final translational temperature of the beam is affected by inelastic collisions-demonstrating the inseparability of elastic and inelastic processes in models that describe buffer gas cell beams. The rotational and translational energy distributions of the simulated buffer-gas-cooled NH 3 beams are in accord with experimental measurements on a similar system. 16 Interestingly, the DSMC simulations indicate that thermalisation occurs well within the typical length of experimental buffer gas cells, which are frequently 10-20 mm long. While the ideal length of the buffer gas cell will depend on the properties of a given system, our results indicate that shorter buffer gas cell designs could be employed in many applications. This would serve to yield a higher flux of cold molecules, without compromising the rotational and translational temperature of the resulting beam.