Classical threshold law for the formation of van der Waals molecules

We study the role of pairwise long-range interactions in the formation of van der Waals molecules through direct three-body recombination processes A + B + B $\rightarrow$ AB + B, based on a classical trajectory method in hyperspherical coordinates developed in our earlier works [J. Chem. Phys. 140, 044307 (2014); J. Chem. Phys. 154, 034305 (2021)]. In particular, we find the effective long-range potential in hyperspherical coordinates with an exact expression in terms of dispersion coefficients of pairwise potentials. Exploiting this relation, we derive a classical threshold law for the total cross section and the three-body recombination rate yielding an analytical expression for the three-body recombination rate as a function of the pairwise long-range coefficients of the involved partners.


I. INTRODUCTION
Van der Waals (vdW) molecules consist of two atoms held together by the long-range dispersion interaction 1 leading to (ground state) binding energies 1 meV. Thus, setting aside ultra-long-range Rydberg molecules (with binding energies ∼ 4 neV 2-6 ), vdW molecules show the weakest gas-phase molecular bond in nature. The binding mechanism in vdW molecules is the result of the compensation between the short-range repulsion, due to the overlap of closed-shell orbitals, and the attractive vdW interaction (−C 6 /r 6 ) caused by zero point fluctuations of atomic dipole moments.
The study of vdW interactions provides a deeper understanding of crucial phenomena in physics, chemistry, and biology. For instance, these interactions play a key role in the formation and stability of gases, liquids, vdW heterostructures and biopolymers; 7-9 chemical reactions; 10-13 superfluidity of 4 He nanodroplets; 14,15 and rare gas crystals and the dynamics of impurities interacting with dense rare gas vapors. 16,17 Recent developments in cooling techniques, specifically buffer gas cooling, 18 have paved the way to new possibilities for investigating the formation of vdW molecules through three-body recombination processes. 17,[19][20][21][22][23][24] Three-body recombination is a threebody collision during which two of the particles form a bound state. These reactions play an important role in a wide range of physical and chemical phenomena, ranging from H 2 formation in star-forming regions 25,26 to loss mechanisms in ultracold dilute atomic gases [27][28][29][30][31][32][33] to the formation and trapping of cold and ultracold molecules. [34][35][36][37][38] In Ref. [39], we considered the formation of atom-rare gas vdW molecules via a direct three-body recombination mechanism at temperatures relevant for buffer gas cell experiments. As a result, we found that almost any atom a) Electronic mail: mirahmadi@fhi-berlin.mpg.de b) Electronic mail: jperezri@fhi-berlin.mpg.de in a helium buffer gas will evolve into a vdW molecule. Fueled by those results, our goal in the present work is to present a comprehensive study of A + B + B reactions and derive a classical threshold law for the formation of vdW molecules in cold environments. To investigate this problem, we use a classical approach in hyperspherical coordinates, which has previously been used to consider the three-body recombination of three neutral atoms, 33,39,40 as well as ion-neutral-neutral three-body recombination processes. 36,37,41 In order to derive a threshold law, we have obtained an effective long-range potential in hyperspherical coordinates. With it, a general expression (as a function of C AB 6 and C B2 6 ) for the corresponding dispersion coefficient is given by considering several A and B atoms chosen from alkali metals, alkaline-earth metals, transition metals, pnictogens, chalcogens, halogens, and rare gases. Furthermore, we calculated the threshold values for the three-body recombination rates at 4 K. Our results confirm that any vdW molecule AB appears with almost the same probability. This paper is organized as follows: In Sec. II, the Hamiltonian governing the classical dynamics during three-body recombination in the three-dimensional space and its counterpart in the six-dimensional space are introduced. In Sec. III, the long-range potential in hyperspherical coordinates have been obtained, and the relevant, effective dispersion coefficient as a function of dispersion coefficients of the pairwise interactions is found. In Sec. IV, a classical threshold lawà la Langevin for the total cross section and the three-body recombination rate are developed. Finally, in Sec. V, we summarize our chief results and discuss their possible applications.
with p i being the momentum vector of the i-th particle. It is more convenient to treat the three-body problem in Jacobi coordinates 42,43 defined by the following relations where R CM 12 = (m 1 r 1 + m 2 r 2 )/(m 1 + m 2 ) is the centerof-mass vector of the two-body system consisting of m 1 and m 2 . M = m 1 + m 2 + m 3 and ρ CM are the total mass and the center-of-mass vectors of the three-body system, respectively. The Jacobi vectors are illustrated in FIG. 1.
Due to conservation of the total linear momentum (i.e., ρ CM is a cyclic coordinate), we omit the degrees of freedom of the center of mass and write the Hamiltonian (1) as with reduced masses µ 12 = m 1 m 2 /(m 1 +m 2 ) and µ 3,12 = m 3 (m 1 + m 2 )/M . Here, P 1 and P 2 indicate the conjugated momenta of ρ 1 and ρ 2 , respectively. Note that, since the relations given by Eq.
(2) indicate a canonical transformation, the Hamilton's equations of motion are invariant under the transformation to Jacobi coordinates.

A. Hyperspherical coordinates
In the next step, we map the independent relative coordinates of the three-body system associated with the Hamiltonian (3) in a three-dimensional (3D) space, onto the degrees of freedom of a single particle moving towards a scattering center in a six-dimensional (6D) space under the effect of the Hamiltonian H 6D . This 6D space is described by means of hyperspherical coordinates consisting of a hyperradius R, and five hyperangles α j (with j = 1, 2, 3, 4, 5), where 0 ≤ α 1 < 2π and 0 ≤ α j>1 ≤ π. 44,45 The components of a 6D vector x = (x 1 , x 2 , x 3 , x 4 , x 5 , x 6 ) in hyperspherical coordinates are given by and the volume element in this coordinate system is given by The 6D position and momentum vectors can be constructed from the Jacobi vectors and their conjugated momenta as 40,46 and respectively. Here µ = m 1 m 2 m 3 /M is the three-body reduced mass. Consequently, the Hamiltonian in the 6D space reads as

B. Total cross section and three-body recombination rate
Classically, for scattering in a 3D space, the cross section σ is defined as the area drawn in a plane perpendicular to particle's initial momentum, which the particle's trajectory should cross in order to be scattered (i.e., deviation from the uniform rectilinear motion). This concept can be extended to the 6D space by visualizing it as an area in a five-dimensional hyperplane (embedded in the 6D space) perpendicular to the initial momentum vector P 0 . Similarly, we define the impact parameter vector b as the projection of the initial position vector ρ 0 on this hyperplane, thus the necessary condition b · P 0 = 0 is satisfied. 46 Note that by treating three-body collision as a scattering problem of a single particle in a 6D space, we can uniquely define the initial conditions and the impact parameter. Hence, it is possible to obtain the probability of a three-body recombination event as a function of the impact parameter b and the initial momentum P 0 . Consequently, by averaging over different orientations of P 0 and making use of its relation with the collision energy, E c = P 2 0 /(2µ), the total cross section of the three-body recombination process will be given by with dΩ b being the differential element of the solid hyperangle associated with vector b. To obtain the second equality we made use of Ω b = 8π 2 /3. The function P in Eq. (9) is the so-called opacity function, i.e., the probability of a recombination event as a function of the impact parameter and collision energy. Note that b max represents the largest impact parameter for which three-body recombination occurs, or in other words, Finally, the energy-dependent three-body recombination rate can be achieved via the following relation

C. Potential
Throughout the present work we make use of the pairwise additive approximation which states that the total potential of a N -body system is the sum of all two-body interactions in the system. Thus, introducing the pairwise potentials U (r ij ), where r ij = | r j − r i |, we write the interaction potential V in Eq. (1) in the following form It is known that the pairwise-additive descriptions of vdW interactions provide appropriate results for the calculation of crystal binding energies 47 (showing deviations 10%), long-range coefficients of small molecules, 48 and the spectroscopy of clusters, 49 although in the latter case it is used only for its convenience. However, there are some scenarios in which a many-body interaction term is required, namely the calculation of long-range coefficients in large molecules 48 and accurate spectroscopic constants of vdW complexes. 47 On the contrary, scattering observables are accurately described at the ultracold limit without invoking many-body interaction terms in the underlying potential energy surface (see, e.g., Ref. [50]). Based on these examples and considering the nature of systems studied in this work, a pairwise approximation for the three-body potential V ( r 1 , r 2 , r 3 ) is convenient.
The relative distances r ij in Cartesian coordinate are related to the Jacobi vectors through the following equations Using these relations together with Eq. (6) we can obtain the potential V ( ρ) in Eq. (8) from Eq. (11). It is important to emphasize that, due to the relations given by Eq. (4), potential is a function of the magnitude of the 6D position vector, ρ, and the corresponding hyperangles

D. Grand angular momentum
To finalize this section, let us briefly explain the notion of grand angular momentum in hyperspherical coordinates. Further below, we will use this discussion to develop a capture model in the hyperspherical coordinate system, and with it a classical threshold law which is the main purpose of this work. In classical mechanics, angular momentum in 6D space is a bivector defined by the exterior product (also known as wedge product) of 6D position and momentum vectors as which is isomorphic to a 6 × 6 skew-symmetric matrix with elements for i, j = 1, 2, . . . , 6. This general definition applies in all higher-dimensional spaces, and for the 3D space it coincides with the familiar cross product. It is worth mentioning that, even though Λ is not equal to the ordinary total angular momentum of the three-body system in the 3D space, it contains the components of the angular momenta (associated with Jacobi vectors) among its elements. Following the original definition by Smith in Ref. [51], Λ is often referred to as the grand angular momentum.
Note that in quantum mechanics, Λ is an operator. Its square, Λ 2 , is the quadratic Casimir operator of so (6) with eigenvalues λ(λ + 4), where λ is a positive integer, and with hyperspherical harmonics as the corresponding eigenfunctions. 45,52,53 FIG. 2. A schematic illustration of the long-range vdw interaction between three particles in 3D space and its counterpart, VLR(ρ), for a single particle in the 6D space.

III. LONG-RANGE POTENTIAL IN HYPERSPHERICAL COORDINATES
Consider a three-body collision A + B + B, where A and B are neutral atoms in their ground electronic state. The interaction potential U (r ij ) consists of a short-range repulsive interaction (due to the overlap of closed-shell orbitals) and a long-range vdW tail for r ij greater than the LeRoy radius. In Ref. [39], we have shown that the formation of vdW molecules through direct three-body recombination at collision energies lower than the dissociation energy of the product molecule (E c < D AB e ) is insensitive to the short-range interaction and is dominated by the longrange tail of the potential. Note that we do not consider the contribution of higher order terms (1/r 8 , 1/r 10 , . . . ) in the long-range tail of the potential, since the effect of long-range interaction on formation of vdW complexes is mainly through the 1/r 6 term. 39 Our goal is to find a general expression for the long-range interaction potential associated with the three-body collision A + B + B, in a 6D space relevant for the classical trajectory method that we employ (see FIG. 2).
Following the discussion in Sec. II C the long-range potential in hyperspherical coordinates is obtained via the relation To fix the coefficients in this equation, we made use of different atoms chosen from alkali metals, alkaline-earth metals, transition metals, pnictogens, chalcogens, halogens, and rare gases. These atoms together with the dispersion coefficients of the pairwise interactions between atoms A and B, C AB 6 , and between two B atoms, C B2 6 , are listed in TABLE I.

A. ρ-distribution
To find the radial dependence of the potential V LR (ρ, α), referred to as V LR (ρ) in FIG. 2, we set the right-hand-side of Eq. (16) equal to a constant, −Q, and solve the equation, for randomly sampled hyperangles. Different hyperangles α j are generated by means of the probability density function associated with dΩ given in Eq. (5) to generate random points uniformly distributed on the 6-sphere (in the geometrical sense). This procedure implies that the solution of Eq. (17) will be obtained as a distribution of ρ values, f (ρ), for each particular Q. Finally, we choose ρ = ρ m with the maximum likelihood in the ρdistribution, as a single value solution of Eq. (17). We select Q from interval [0.001, 1000] K and solve the equation for 10 4 randomly generated sets of hyperangles. Note that, even though for each set of α, Eq. (17) will be transformed into a sixth-degree equation of the variable ρ, we have seen that (regardless the value of Q) four of the roots are complex and from the two remaining real roots only one is positive. Interestingly enough, for all chosen values of Q and dispersion coefficients, the probability density function (PDF) of the ρ-distribution, F f (ρ), is described by the PDF of the generalized extreme value (GEV) distribution, i.e., The GEV distribution is a family of continuous probability distributions developed within the extreme value theory. 65,66 It is parametrized with a shape parameter ξ = 0, a location parameter β, and a scale parameter δ. It is worth mentioning that here ξ > 0,  [56], and for Ti-He from Ref. [39]. b Dispersion coefficient C 6 is taken from Ref. [57]. c Dispersion coefficient C 6 is taken from Refs. [58][59][60]. d Dispersion coefficient C 6 for Sr-B is taken from Ref. [61] and for pnictogens-He from Ref. [62]. e Dispersion coefficient C 6 for O-B is taken from Ref. [63] and for F-B from Ref. [64]. f Dispersion coefficient C 6 obtained through the London-Drude theory of dispersion interactions and Taken from Ref. [8].
which corresponds to the type II (also known as Fréchet distribution) of GEV distributions. 66 The parameters of GEV in Eq. (18) are fitted yielding an uncertainty below 0.05%. As an illustration, the distribution f (ρ) obtained by setting Q = 1 mK in Eq. (17) for the three-body collision As + He + He is shown in FIG. 3. As it can be seen in this figure, the probability density is positivelyskewed. Therefore, due to this skewness, the maximum ρ m is achieved by using the mode of the GEV distribution, instead of its mean value. It is worth mentioning that we have observed the same (hyper)radial distribution f (ρ) also in the case of pairwise interactions proportional to 1/r 4 . In other words, if we change an atom by an ion, the distribution of the hyperradius ρ is of the same kind. However, the results of such a study are beyond the scope of this work and will be considered elsewhere.

B. Effective dispersion coefficient
The general form of the long-range potential V LR (ρ) in the 6D space can be derived from the power-law relation between Q and ρ (= ρ m ). An example has been shown in FIG. 4. It can be seen that, as expected, Q is proportional to 1/ρ 6 and the corresponding effective dispersion coefficient, C eff , is the slope of the fitted line in the log- log scale. Therefore, C eff can be obtained as a function of C AB 6 and C B2 6 of the pairwise vdW interactions. Finally, utilizing Eqs. (16) and (17), we have Through the same procedure for all A-B-B systems mentioned in TABLE I and calculating the corresponding vdW coefficients C eff , we have found the general expression and C B2 6 .

IV. CLASSICAL THRESHOLD LAW
A threshold law for the three-body recombination cross section and rate can be established based on the pioneering capture theory of Langevin. 67 In the framework of this classical capture model, every trajectory associated with the collision energy (E c ) above the potential barrier leads, with unit probability, to a reaction event. To implement the same idea in 6D space, we first need to define the effective potential.
In a three-body collision, after including the centrifugal energy, the effective long-range potential in hyperspheri- cal coordinates reads as 51 (in atomic units) with a maximum (the so-called centrifugal barrier) at can be obtained from the components given by Eq. (14) and after applying the (algebraic) Lagrange's identity 68 one finds Utilizing the relation between the impact parameter vector b and the initial position and momentum vectors, 40 we have ρ 0 ∧ P 0 = b ∧ P 0 , where we used the fact that P 0 ∧ P 0 = 0. Finally, taking into account that b ⊥ P 0 , one finds which establishes the intimate relation between the grand angular momentum and the magnitude of the impact parameter.
Knowing that a reaction occurs only if E c ∼ V eff (ρ 0 ), we can find a threshold value for the impact parameter, b max , which is assigned to E c = V eff (ρ 0 ). Upon substituting for Λ 2 obtained from Eq. (24) into V eff (ρ 0 ), we derive the relation Inserting Eq. (25) into Eq. (9) we obtain the geometric cross section (P(E c , b) = 1 for b ≤ b max ) in the following form where in the last line we made use of Eq. (21). Employing Eq. (10), the three-body recombination rate can be calculated as a function of collision energy and dispersion coefficients of the pairwise interactions as where k B is the Boltzmann constant and Γ(x) is the gamma function of argument x. We should emphasize that these relations are best valid for collision energies smaller than the dissociation energy of the vdW molecules, which is typically below 100 K (≈ 10 meV). Moreover, one should also verify the validity of the classical approach based on the number of involving partial waves, which will be discussed in what follows.
A. Estimated number of contributing three-body partial waves From the perspective of a scattering problem of a single particle in a 6D space, each (grand) angular momentum quantum number λ is a so-called partial wave. Following this fact, we introduce λ as the partial wave associated with a three-body recombination in 3D space.
It is known that the reliability of the classical approach depends on the number of partial waves contributing to the scattering observables. 46,69 In other words, the large number of partial waves (≈ 20) contributing to the scattering washes out quantum effects such as resonances. The number of partial waves that impart a significant effect to the scattering problem can be estimated from the strength of interaction, i.e., the collision energy E c . 39,46 Setting Λ 2 equal to the eigenvalues of its quantum mechanical counterpart, λ(λ + 4), from V eff (ρ 0 ) = E c we can establish the following relation between the maximum three-body partial wave λ max and a given collision energy where we made use of the fact that for λ 4, λ(λ+4) → λ 2 . Equation (29) provides a measure to check the validity of classical calculations for different A-B-B systems based on the collision energy, reduced mass, and pairwise dispersion coefficients. For a more detailed comparison between the quantum and classical results obtained by hyperspherical classical trajectory method in three-body recombination see Ref. [40].

B. Low-energy limit: s-wave collisions
In the final part of this section we derive a classical threshold law associated with the quantum s-wave scattering, i.e., λ = 0. In this case one may define the parameter b max as the distance at which the collision energy is comparable to the strength of the interaction potential, i.e., E c = C B2 6 r −6 12 + C AB 6 r −6 23 + C AB 6 r −6 31 in 3D space or equivalently E c = C eff ρ −6 in 6D space. Therefore, the maximum impact parameter in the hyperspherical coordinate system reads as By a similar argument as above, we obtain the cross section from Eq. (9), which yields Consequently, the three-body recombination k sw 3 (E c ) and its thermal average k sw 3 (T ) are given by and respectively.

C. Results
The results derived by performing the thermal average (28) for different A+B+B reactive collisions for T = 4 K (relevant for buffer gas cells) are shown in TA-BLE I. To calculate these values we used the mass of the most abundant isotopes of A and B atoms. Note that the calculated recombination rates account for both AB and B 2 products of the three-body process. However, based on the relative values of the dispersion coefficients, AB molecules are formed more often than B 2 ones, unless the dispersion coefficient for B 2 is larger than that of AB. It is important to notice that all calculated three-body recombination rates are of the same order of magnitude. One reason is the very close values of C eff obtained for different systems, which almost neutralizes the effect of the three-body reduced mass µ.
TABLE II shows the three-body recombination rates k 3 (T ) given by Eq. (28) for six different A+He+He collisions at T = 4 K, together with values of k num 3 (T ) taken from Ref. [39]. k num 3 (T ) are the numerical values calculated via the classical trajectory method introduced in Ref. [40]. The recombination ratesk 3 (T ) in this table are obtained from a capture model that only takes into account the pairwise interaction of the stronger long-range tail, i.e., 39k It can be seen that the trend of the rates calculated with a capture modelà la Langevin (k 3 ) is in reasonably good agreement with the trend of k num 3 except for Li and Na. This is because the dissociation energy of LiHe and NaHe are below 2 K. Hence, the collision energies considered to calculate k 3 (T = 4 K) have reached the high-energy regime (E c > D e ). Since this regime is sensitive to the short region of the potential, the capture model is not as accurate as in the other cases. Finally, we must highlight the considerable improvement (about one order of magnitude) in the threshold values k 3 obtained from Eq. (28), over those values derived from the threshold law given by Eq. (34).

V. CONCLUSIONS AND PROSPECTS
After developing a clear picture of the (hyper)radial dependence of a three-body potential in a 6D space and studying more than 40 three-body systems relevant for vdW molecule formation, we have found how the longrange interaction of three-body systems depends on pairwise interactions between the colliding partners. Then, employing a classical trajectory method in hyperspherical coordinates, we have established a classical threshold law for the formation of vdW molecules through direct three-body recombination processes relevant for buffer gas cooling experiments. In addition, we have shown that at a given temperature, the three-body recombination rate is of the same order of magnitude independently of the atomic species under consideration, which corroborates our previous studies on the matter.
The most valuable achievement of this work is to offer a simple expression which makes it possible to obtain the three-body recombination rate by only using the long-range dispersion coefficients and masses of the colliding atoms. This result helps to quickly estimate the role of three-body recombination in a given scenario, and with it, provides a new avenue for the calculation of three-body recombination rates avoiding costly computations. Finally, we hope that our findings help to make three-body collisions more approachable for the chemical physics community.