A rigorous adiabatic approach to ultracold atom–molecule collisions in a magnetic field

Abstract

We extend the rigorous adiabatic coupled-channel formalism to ultracold nonreactive atom–molecule collisions in the presence of an external magnetic field. The wavefunction of the collision complex is expanded in adiabatic basis states obtained by solving the eigenvalue problem for the adiabatic Hamiltonian (the total Hamiltonian of the collision complex minus the radial kinetic energy) on a grid of atom–molecule distances R. The resulting coupled-channel equations are solved using the diabatic-by-sector method. We show that the adiabatic approach provides accurate cross sections for cold and ultracold Mg (¹S) + NH (³Σ⁻) collisions in a magnetic field with ≃2 times fewer channels than the standard diabatic basis. We further develop an efficient R-dependent basis truncation protocol (RBT), in which the elements of the log-derivative matrix are sampled and discarded as it is propagated from small to large R. When implemented in the standard diabatic basis, RBT affords a computational gain of more than one order of magnitude. The adiabatic basis can be truncated even more aggressively as a function of R to include just the open channels at long range. This leads to an overall computational gain of ≃15–30 for the propagation part of the calculation, with direct CPU–time measurements indicating reductions exceeding a factor of 50. The gain is particularly significant in situations where substantial errors in the calculated cross sections (<50%) can be tolerated, making the adiabatic basis formulation a promising approach to strongly anisotropic collisions and chemical reactions in the presence of an external magnetic field.

---

I. Introduction

Ultracold molecular gases are a promising platform for quantum science,^1–4 quantum control of chemical reactions,^5–9 and precision tests of fundamental physics.^10–12 These applications demand precise quantum control over molecular degrees of freedom, including their translational motion and internal electronic, rovibrational and hyperfine states. Major experimental advances have been made over the last decade in developing the tools for such control.^13,14 In particular, cooling, trapping, and state-selection of molecules^13,14 have enabled state-to-state tracking of reactants and products in the ultracold KRb + KRb → K₂ + Rb₂ chemical reaction,¹⁵ observation of a Feshbach resonance in ultracold molecule–molecule collisions,¹⁶ association of ultracold triatomic molecules via atom–molecule Feshbach resonances,¹⁷ magnetic control of three-body recombination products,¹⁸ energy transfer from hyperfine to rotational degrees of freedom,¹⁹ and the realization of a Bose–Einstein condensate of ground-state molecules.²⁰

Molecular collisions are central to the ability to control intermolecular interactions in the ultracold regime.^5,6,10,11 Ultracold molecular gases typically exhibit universal loss, in which detrimental inelastic collisions occur with near-unit probability once the molecules approach each other at close range.²¹ However, elastic collisions are necessary for collisional cooling techniques such as sympathetic and evaporative cooling.^{10,11,22–24} Ideally, collisions are expected to be of use as a mechanism to control ultracold molecular gases through Feshbach resonances that are tunable with external fields,¹⁶ as they have been for ultracold atomic gases.²⁵ Although magnetic Feshbach resonances have been observed in ultracold atom–molecule and molecule–molecule collisions of K + NaK,^17,26,27 Na + NaLi,^28,29 and NaLi + NaLi,¹⁶ the mechanisms behind these resonances have yet to be fully elucidated because of enormous computational challenges that generally prohibit rigorous quantum dynamics calculations.^21,30,31 These challenges stem from the large number of molecular states (rotational, vibrational, fine, and hyperfine) coupled by strongly anisotropic intermolecular interactions,^30,32 and the steep scaling of numerically exact coupled-channel calculations (see below), which currently limit the maximum number of scattering channels to about 18 500,³³ with most practical calculations becoming computationally unfeasible already for >10 000 channels. Advancements in molecular quantum scattering theory and computational methods are therefore needed to analyze current and past experiments on ultracold molecular collisions and to provide guidance for future experiments.

The most complete and rigorous theoretical description of molecular collisions is based on the quantum coupled-channel (CC) approach that solves the Schrödinger equation numerically exactly for a given Hamiltonian.^34,35 These calculations rigorously account for the coupling of molecular rotational and orbital angular momenta, as well as electronic and nuclear spins. In addition, the inclusion of excited vibrational and electronic states is sometimes necessary.¹⁹ Highly anisotropic interactions between the collision partners can couple hundreds of rotational states,^19,31,36 necessitating enormous basis sets to achieve convergence. Because CC calculations scale cubically with the number of basis states (or scattering channels), they tend to be computationally intensive and quickly become intractable as more basis states and/or degrees of freedom are added.

Fortunately, the efficiency of CC calculations on molecular collisions in external fields can be drastically improved by using optimized channel basis sets, such as those based on the total angular momentum (TAM)^37,38 or the total rotational angular momentum (TRAM) representations.^{19,32,36,39–41} The TRAM is the vector sum of all the angular momenta for mechanical rotation (i.e., the rotation of the diatomic fragment and the orbital motion of the atom in the atom–molecule collision complex). It is rigorously conserved in the absence of anisotropic spin-dependent interactions, which are often weak compared to the rotational energy splitting and short-range forces.^19,32,36 By leveraging this property, Morita et al.³⁶ have recently succeeded in obtaining converged cross sections for ultracold Rb + SrF collisions in a magnetic field based on ab initio potential energy surfaces (PESs) in the rigid-rotor approximation. Their CC calculations fully included the hyperfine structure of both Rb and SrF in an external magnetic field, providing the first rigorous insights into magnetic Feshbach resonance spectra in strongly anisotropic atom–molecule collisions.³⁶ By contrast, evidence of TRAM non-conservation has been reported in a recent experimental and theoretical study of hyperfine-to-rotational energy transfer in ultracold Rb + KRb collisions.¹⁹ The breakdown of TRAM conservation is likely induced by short-range spin-dependent interactions and conical intersections between the ground and first excited PESs of the Rb–KRb collision complex.¹⁹

The CC calculations can also be made more efficient by neglecting certain weak interactions in the short-range region where the anisotropic atom–molecule interactions dominate. Recent studies of ultracold Mg + NH collisions have found that intramolecular hyperfine⁴² and Zeeman^42,43 interactions can be neglected at short-range, which can be understood using multichannel quantum defect theory with a frame transformation (MQDT-FT).⁴² This allows for a several orders-of-magnitude reduction in computational effort compared to full CC calculations including all degrees of freedom.⁴² However, even though MQDT-FT can give quantitatively accurate results for ultracold atom–molecule collisions in a magnetic field without explicitly including the hyperfine structure, its limits remain to be explored. For example, the strong short-range coupling between the spin and rotational degrees of freedom in ultracold Rb + KRb collisions¹⁹ can pose challenges for the MQDT-FT approach, which benefits from neglecting this coupling.⁴²

All previous CC calculations on ultracold atom–molecule and molecule–molecule collisions in magnetic fields have used diabatic basis functions, which are independent of the atom–molecule distance R. Examples include the uncoupled space-fixed basis^44–51 and the computationally efficient TAM^{22,23,38,52,53} and TRAM^19,32 representations. A disadvantage of the diabatic basis is that it does not account for the peculiar properties of atom–molecule interactions, which typically vary dramatically with R. Specifically, at small R these interactions are often strongly anisotropic, coupling a large number of diabatic rotational states.^{22,23,38,42,52} At intermediate-to-long range, the anisotropy decreases and the couplings between the basis states decline sharply. A single R-independent diabatic basis lacks the flexibility to capture these variations, necessitating the use of extensive diabatic basis sets to achieve convergence of scattering observables, even when the computationally efficient TAM or TRAM representations are used.^22,23,42

An appealing alternative is offered by the adiabatic basis composed of the eigenstates of the full atom–molecule Hamiltonian evaluated at a fixed value of R. The corresponding eigenvalues, known as the adiabatic potentials, describe how the atom–molecule interaction energy varies as a function of R. Because the adiabatic basis states explicitly account for the atom–molecule interaction, they are expected to provide a better, more compact description of the quantum states of the collision complex than their diabatic counterparts. Accordingly, the adiabatic approach is widely used in the theoretical studies of ultracold few-body recombination⁵⁴ and reactive scattering.^35,55,56 The adiabatic hyperspherical method has enabled rigorous insights into the quantum dynamics of four-body and five-body recombination^57,58 and ultracold chemical reactions Rb + Rb + Rb → Rb₂ + Rb¹⁸ and K + KRb → K₂ + Rb.⁵⁹ However, the adiabatic basis has not yet been used in rigorous calculations of ultracold nonreactive atom–molecule collisions in the presence of an external magnetic field.

Early implementations of the adiabatic formulation for model atom–molecule collisions in the absence of external fields were based on the R-matrix propagation method.^60,61 More recently, log-derivative propagation in the (quasi)-adiabatic basis was implemented in the quantum scattering package MOLSCAT⁶² even though, to our knowledge, it has not been used in practical calculations. Instead, adiabatic potentials are often used as a valuable tool to analyze the results of CC calculations, as in exploring the density of magnetic Feshbach resonances in ultracold Rb + SrF collisions,³⁶ explaining the strong dependence of ultracold reaction rates on the initial states of the reactants,⁶³ and elucidating the mechanisms of rovibrational energy transfer in H₂O + H₂ collisions.⁶⁴

A variety of classical and quantum adiabatic capture theories rely on adiabatic potentials to determine the flux from the initial channel that is transmitted from long-range to short-range.^65–72 One can then approximate reaction rates by assuming that the flux is absorbed with unit probability due to the formation of a long-lived complex at short range. Statistical adiabatic models are similar to adiabatic capture theories, but assign channel-dependent probabilities for reactive or non-reactive transitions that can take a range of values instead of assuming unit probabilities for all channels.^73–75 They are also instrumental in analyzing the mechanisms of microwave and electric field shielding.^20,76–80 Notwithstanding their excellent performance for ultracold collisions of alkali-dimer molecules,^20,76,77 these methods rely on the fundamentally unphysical absorbing boundary condition at short range,⁸¹ making them unsuitable for rigorous quantum dynamics calculations.

In this paper, we develop a rigorous adiabatic methodology for CC calculations on ultracold atom–molecule collisions in a magnetic field. We apply the methodology to calculate the cross sections for ultracold Mg + NH collisions in a magnetic field, using a realistic ab initio PES. Significantly, we find that converged scattering cross sections can be obtained with far fewer adiabatic channels. However, this advantage is somewhat offset by the additional computational effort needed to generate the adiabatic basis states, and to calculate the overlap matrices between the adjacent sectors. We also quantify the numerical performance of the adiabatic vs. diabatic basis sets.

In addition, we explore the possibility of truncating the CC basis set “dynamically”, as the log-derivative matrix is propagated from small to large R. We found that doing so in the adiabatic basis leads to a 15–30 fold reduction in computational cost compared to the diabatic basis for the propagation part of the calculations. Our R-dependent basis set truncation (RBT) procedure is an extension of the approach proposed by Light and co-workers.⁶¹ The key element of our RBT protocol is that, instead of truncating the intersector overlap matrix, as done in ref. 61, we truncate the log-derivative matrix. As a result, our method is more general and can be applied in the adiabatic as well as diabatic formulations.

The structure of this paper is as follows. In Section II, we formulate the adiabatic CC theory for ultracold atom–molecule collisions in a magnetic field and describe the RBT approach. In Section III, we apply RBT to ultracold Mg + NH collisions in a magnetic field and discuss the computational gains achieved over the standard diabatic CC approach. Section IV concludes by summarizing the advantages the adiabatic approach has over the diabatic method, and the regimes in which the adiabatic basis is expected to provide the maximal advantage.

II. Theory

We begin by presenting the adiabatic coupled-channel (ACC) approach to ultracold atom–molecule collisions in a magnetic field before describing our R-dependent basis truncation (RBT) algorithm. In the next section, the ACC/RBT approach will be applied to ultracold Mg–NH collisions in a magnetic field. Throughout this paper, we restrict attention to nonreactive atom–molecule collisions in a magnetic field. The reader is referred to ref. 82 for an extension of the adiabatic approach to chemically reactive scattering in an electric field using hyperspherical coordinates.⁵⁶

A. Hamiltonian and adiabatic coupled-channel (ACC) equations

The Hamiltonian for an atom–molecule collision complex A + BC can be written as (in atomic units)^34,44,45where μ is the reduced mass of the complex, and R and r are the standard Jacobi vectors: R extends from the center of mass of BC to A and r joins the nuclei of BC. The adiabatic Hamiltonianincludes, in order, the centrifugal kinetic energy due to the orbital motion of A around BC, the interaction potential describing the interaction between the atom and molecule, and the asymptotic Hamiltonian that describes the isolated (i.e., non-interacting) atom and molecule, defined asIn eqn (2), L̂ is the orbital angular momentum that describes rotations of R, R = |R|, r = |r|, and θ is the angle between the vectors R and r. We assume that the A–BC interaction is adequately described by a single adiabatic potential energy surface (PES), and freeze the internuclear distance of the diatomic molecule at its equilibrium value, r = r_e. The rigid-rotor approximation is known to be reliable for rigid diatomic molecules such as NH colliding in their ground vibrational states.^44,45

Here, we consider the prototypical example of a binary collision between a ³Σ diatomic molecules and a structureless ¹S-state atom in an external magnetic field B. The asymptotic Hamiltonian then reduces to that of an isolated ³Σ molecule such as the NH radical⁴⁵

Ĥ_as(r) = B_eN̂² + γ_srN̂·Ŝ + 2μ₀B·Ŝ + V̂_SS (4)

where N̂ and Ŝ are the angular momentum operators that describe molecular rotation and electronic spin, respectively, B_e is the rotational constant, γ_sr the spin–rotation constant, μ₀ is the Bohr magneton, and V̂_SS is the spin–spin interaction.⁴⁵

Following the general strategy of the adiabatic approach,^34,54,57,83 we expand the wavefunction of the triatomic atom–molecule collision complex at a total energy E in the adiabatic basis as

Here, F_ii0(R) are the radial solutions of the time-independent Schrödinger equation Ĥψ_i0 = Eψ_i0 which form a square M × M matrix F(R), whose columns represents linearly independent radial solutions for a given initial state i₀.³⁴ The fact that M can be R-dependent will be leveraged in the development of the RBT protocol in Section IIC.

The adiabatic basis functions Φ_i are eigenstates of the adiabatic Hamiltonian

Ĥ_adΦ_i(r,R̂;R) = ε_i(R)Φ_i(r,R̂;R), (6)

where ε_i(R) are the eigenvalues of the adiabatic Hamiltonian, also known as adiabatic potentials. The eigenfunctions Φ_i(r,R̂;R) are also known as adiabatic surface functions.^55,56 The adiabatic eigenvalues and eigenfunctions depend on R only parametrically, and approach those of the asymptotic Hamiltonian as R → ∞, which is a direct consequence of eqn (3). The adiabatic eigenvalue problem (6) is solved at a fixed value of the atom–molecule distance R by expanding the adiabatic eigenfunctions in a primitive R-independent basis setwhere is the number of primitive basis functions used in the expansion. The expansion coefficients T_ji(R) are found by solving the matrix eigenvalue problemobtained by substituting eqn (7) into eqn (6). The expansion coefficients T_ji(R) form columns of the eigenvector matrix T(R), which diagonalizes the matrix of the adiabatic Hamiltonian in the primitive basis with elements 〈χ_i|Ĥ_ad|χ_j〉.

To solve the adiabatic eigenvalue problem for our triatomic complex composed of a rigid-rotor ³Σ molecule colliding with a ¹S atom, we use the uncoupled space-fixed (SF) representation^44,45 as our primitive basis

|χ_j〉 = |NM_N〉 |SM_S〉 |LM_L〉. (9)

These basis states are simultaneous eigenstates of N̂², N̂_z, Ŝ², Ŝ_z, L̂², and L̂_z. The SF quantization axis points along the magnetic field vector such that B = Bẑ. The uncoupled SF basis is well suited for weakly and moderately anisotropic interaction potentials. For strongly anisotropic interactions, the total angular momentum³⁷ or total rotational angular momentum³² bases are much more computationally efficient.

The matrix elements of the adiabatic Hamiltonian in the primitive basis 〈χ_i|Ĥ_ad|χ_j〉 required to set up the matrix eigenvalue problem (8) are obtained as described in ref. 45. These include the matrix elements of the asymptotic Hamiltonian, the centrifugal kinetic energy, and the atom–molecule interaction potential from eqn (2)

All of these matrix elements can be evaluated as described in the seminal papers by Volpi and Bohn⁴⁴ and Krems and Dalgarno.⁴⁵ For example, to compute the matrix elements of the atom–molecule interaction potential, we expand it in Legendre polynomials

where V_λ are the R-dependent expansion coefficients, and use eqn (12) of ref. 45.

As a final step, to obtain the adiabatic radial solutions F_ii0(R) in eqn (5) we plug eqn (6) into the time-independent Schrödinger equation, which yields the coupled-channel equations in the adiabatic basis (the ACC equations)^34,54,57,83

where integration in the radial derivative coupling matrix elements 〈…〉 is performed over r. The boundary conditions on the radial solutions F(R) arewhere S_ii0 is the scattering matrix element and l_i and are the orbital angular momentum quantum number and wavevector of the ith adiabatic channel, respectively (l₀ and k_i0 for the initial channel). In the limit R → ∞ the adiabatic potentials tend to the eigenvalues of Ĥ_as, which define the internal (rovibrational and Zeeman) energy levels of the isolated diatomic molecule ε_i = ε_i(R → ∞). We will assume that the initial state is fixed throughout, and omit the subscript i₀ from now on.

B. The diabatic-by-sector formulation

The ACC eqn (12) include the matrix elements of the radial derivative operators and , which couple different adiabatic states. These non-adiabatic couplings are challenging to handle numerically because they diverge at avoided crossings,³⁴ which are common in atom–molecule systems (see below).

To avoid dealing with the non-adiabatic couplings, we use the diabatic-by-sector (DBS) approach commonly employed in quantum reactive scattering calculations,^55,56 in which the entire range of R is split into small sectors. Within each sector, the adiabatic states are independent of R, so non-adiabatic couplings can be neglected and eqn (12) becomes

The log-derivative matrix is propagated within each sector using eqn (14) and transformed from the adiabatic basis of the current sector to that of the next by a unitary transformation. At the final atom–molecule separation (R = R_final), the log-derivative is transformed to the eigenbasis of the asymptotic Hamiltonian Ĥ_as and scattering boundary conditions are applied to obtain the S-matrix. The details of this procedure are described below.

As shown in Fig. 1, the nth sector is defined by three points: the initial point R = a_n, the midpoint R = c_n, and the endpoint R = b_n. Consecutive sectors share a border, such that the initial point of the (n + 1)th sector is the endpoint of the nth sector (a_n+1 = b_n). The adiabatic eigenvalue problem is solved at the midpoint of each sector by expanding the adiabatic states in the primitive diabatic basis [eqn (7) and (9)]. The matrix elements of the adiabatic Hamiltonian in this basis are evaluated analytically as described in, e.g., ref. 45. Numerical diagonalization of the resulting matrix provides the adiabatic energies and eigenvectors as a function of R. The former are used in eqn (14) whereas the latter are employed to construct the sector-to-sector overlap matrix as described below.

Fig. 1 Schematic of the diabatic-by-sector procedure for propagating the log-derivative matrix in the adiabatic basis. The nth sector is defined on an interval R ∈ [a_n,b_n] by the starting point R = a_n, the midpoint R = c_n, and the endpoint R = b_n. The endpoint of the nth sector is the starting point of the following sector (a_n = b_n−1). The blue arrows show the propagation and sector-to-sector transformation of the log-derivative matrix.

The log-derivative matrix in the adiabatic basis of the nth sector defined at an R-value within the sector (R = R_n) is

Y(R_n) = F′(R_n)[F(R_n)]⁻¹, (15)

The subscript n denotes the adiabatic basis in which the log-derivative is represented. To propagate the log-derivative matrix within a sector from R = a_n to R = b_n, we use the methods proposed by Johnson and Manolopoulos.^84,85 Following ref. 85, we define the reference potential of the nth sector in the adiabatic basis defined at R = c_n as

W_ref(c_n) = 2μT^T_nH_ad(c_n)T_n − 2μE1, (16)

where 1 is the identity matrix, T_n is the adiabatic eigenvector matrix of sector n, and T^T_nH_ad(c_n)T_n is the diagonal matrix of adiabatic eigenvalues ε_i(c_n).

Note that H_ad is the adiabatic Hamiltonian written in the primitive basis with matrix elements [H_ad]_ij = 〈χ_i|Ĥ_ad|χ_j〉. The residual coupling matrix at any value of R in the nth sector is

U(R) = 2μT^T_n[H_ad(R) − H_ad(c_n)]T_n. (17)

The residual coupling matrix approaches zero in the limit of vanishing sector width ΔR, and hence has a negligible effect on the propagation of the log-derivative matrix, provided that ΔR is made small enough. We find that for ultracold Mg + NH collisions studied in this work, it is more computationally efficient to use a smaller sector width to avoid the expensive transformation required to compute the residual coupling matrices. Thus, we set U(R) → 0 in all calculations presented in this paper.

The residual coupling matrix vanishes at the sector midpoint (R_n = c_n), so the quadrature matrix in the improved log-derivative method⁸⁵

also vanishes, eliminating one inversion that is required when using a diabatic basis. In eqn (18), h = ΔR/2 is the steps size, equal to the half-width of the sector. The trade-off for eliminating the inversion is the diagonalization that needs to be performed to solve the adiabatic eigenvalue problem (eqn (6)). We note that this inversion, like the diagonalization required for the adiabatic treatment, is independent of the total energy, and therefore only needs to be calculated when the adiabatic Hamiltonian changes, like when changing the magnetic field strength or scaling the interaction potential.

Once the log-derivative matrix is propagated across the sector, Y(a_n) → Y(b_n), it is transformed to the adiabatic basis of the next sector as⁵⁵

Y(a_n+1) = O^T_n→n+1Y(b_n)O_n→n+1, (19)

where the overlap matrix is computed from the adiabatic eigenvector matrices from subsequent sectors:

O_n→n+1 = T^T_nT_n+1. (20)

The sector-to-sector similarity transformation (19) can be regarded as consisting of two steps, the first of which transforms the log-derivative matrix from the adiabatic basis of the nth sector to the primitive basis, and the second from the primitive basis to the adiabatic basis of the (n + 1)th sector.

On reaching the last propagation sector centered at R_n = R_Ns in the asymptotic region of large R, the log-derivative matrix is transformed back to the diabatic basis, and then to a basis, which diagonalizes the asymptotic Hamiltonian (3).⁴⁵ The asymptotic transformation matrix is a product of two orthogonal matrices

where C is the matrix of eigenvectors of the asymptotic Hamiltonian. As this transformation brings the log-derivative matrix to the proper diabatic basis of eigenstates of the asymptotic Hamiltonian, one can now apply the standard K-matrix boundary conditions as described by Krems and Dalgarno⁴⁵ to obtain the K and S-matrices using the standard expressions, e.g. ref. 84. The key observables of interest—the integral cross sections and transition probabilities between the different internal molecular levels induced by atom–molecule collisions—are obtained from the S-matrix using well-documented expressions.^34,45

C. Truncating the adiabatic basis during propagation

As noted above, rigorous CC calculations of non-reactive molecular collisions in external fields typically use the same R-independent (diabatic) basis throughout the entire range of R. A more computationally efficient procedure would be to reduce the size of the basis (diabatic or adiabatic) during propagation. Stechel, Walker, and Light⁶¹ implemented such a procedure by sampling the elements of the sector-to-sector overlap matrix in the adiabatic basis (20) but provided no details on either the implementation or performance of the method. In addition, their approach cannot be implemented in the diabatic basis, for which the sector-to-sector overlap matrix reduces to the identity matrix. In this section, we propose and describe an alternative method based on sampling the log-derivative matrix, which works in both the diabatic and adiabatic formulations.

In the remainder of this paper, we will refer to the general class of methods that reduce the size of the basis during propagation as R-dependent basis truncation (RBT). We refer to our version as log-derivative-based RBT and the version presented by Stetchel, Walker, and Light as overlap-based RBT. Below, we focus on log-derivative-based RBT in the adiabatic basis applied to ultracold atom–molecule collisions. Additional details regarding the performance of RBT are relegated to Appendix B. The RBT procedure in the diabatic basis, along with the original overlap-based RBT,⁶¹ is described in Appendix B.

As noted above, in log-derivative-based RBT we sample the elements of the log-derivative matrix as it is propagated from small to large R. In each propagation sector centered at R = c_n, the propagation basis size M_n defines the dimension of the log-derivative matrix (dim(Y) = M_n × M_n). Note that the propagation basis size is less than or equal to the primitive basis size used to solve the adiabatic eigenvalue problem defined by eqn (6). Our aim is to reduce the propagation basis size to lower the computational cost of matrix operations required to propagate the log-derivative matrix across a sector and to transform it to the next sector.

In a given sector, the ith channel is locally open if ε_i(R) < E and locally closed otherwise, where E is the total energy. Let M^(LO)_n and M^(LC)_n be the numbers of locally open and locally closed channels in the nth sector, respectively (see Fig. 2). A weakly closed channel is any closed channel that becomes locally open in any sector, i.e., satisfies ε_i(c_n) < E for at least one value of c_n. We start RBT at the R-value corresponding to the minimum of the isotropic part of the interaction potential, V₀(R). After that condition is met, the following algorithm describes how the log-derivative matrix is sampled to determine how many channels can be removed in the nth sector.

The outermost for i loop runs over the locally closed channels, starting with the channel with the greatest adiabatic energy in the nth sector (we assume that the adiabatic eigenvalues are sorted in the order of increasing energy). The next for j loop runs over the locally open channels and computes the maximum absolute value (max_{_}value) of the log-derivative matrix element between the ith closed channel and the locally open channels.

The following if … then statement checks whether max_value is less than or equal to a predetermined RBT threshold τ_RBT. If so, the closed channel is removed, and the propagation basis size of the following sector M_n+1 is reduced by one. Otherwise, the closed channel is kept, and the log-derivative is transformed to the next sector as usual.

M_n+1 = M_n

for i from M_n to M_n − M^(LC)_n + 1

 max_value = 0

 for j from 1 to M^(LO)_n

  Find the maximal matrix element [Y(b_n)]_ji → max_value

 end if

 if max_value ≤ τ_RBT then

  remove the ith channel

  M_n+1 = M_n+1 − 1

 else

  exit loop over i

 end if

end for

Once the RBT algorithm determines which channels are to be removed, the corresponding rows and columns of the log-derivative matrix are deleted in the sector-to-sector transformation step, in which the log-derivative matrix is transformed from the basis of one sector to that of the following sector [eqn (19)].

For example, suppose we find that the basis can be reduced from M_n to M_n+1 < M_n after propagating the log-derivative through the nth sector. The M_n × M_n+1 overlap matrix O_n→n+1 = T^T_nT_n+1 is then constructed from the first M_n columns of T_n and the first M_n+1 columns of T_n+1. Using this overlap matrix to transform the log-derivative matrix via eqn (19) reduces the dimension of the log-derivative matrix from M_n × M_n to M_n+1 × M_n+1.

Fig. 2 Schematic representation of the log-derivative matrix (Y) in the adiabatic basis as a square matrix of dimension M_n × M_n. The log-derivative matrix is divided into four submatrices whose dimensions are based on the number of locally open channels (M^(LO)_n) and the number of locally closed channels (M^(LC)_n). The matrix is symmetric, so the green sections are transposes of one another.

The above truncation procedure progressively reduces the size of the log-derivative matrix as it is propagated from small to large R. This, in turn, reduces the cost of the two inversions required to propagate the log-derivative matrix across a sector from to O(M_n³).

III. Results

A. Adiabatic coupled-channel calculations for ultracold Mg + NH collisions in a magnetic field

In this section, we apply the rigorous adiabatic approach developed in Section II to calculate the cross sections for ultracold Mg + NH collisions in a magnetic field. The Mg–NH interaction is moderately anisotropic and the NH molecule is well-represented as a rigid rotor, enabling fully converged CC calculations using standard diabatic basis sets.^42,43,86,87 This makes ultracold Mg + NH collisions a convenient test system, allowing us to benchmark the adiabatic approach against the standard diabatic CC calculations.

Fig. 3 shows the magnetic field dependence of the lowest rotational-Zeeman levels of an isolated NH(³Σ) molecule. The rotational levels are arranged in manifolds, with the rotational spacing between the manifolds being much larger than the intramolecular spin–spin interaction.^45,88 We can therefore label each manifold by its rotational quantum number N. The dashed line shows the low-field seeking Zeeman energy level of the ground rotational state, which is the initial state used in all calculations presented here. Inelastic transitions from this state (also known as Zeeman or spin relaxation) cause loss of NH molecules from a magnetic trap.^45,89–91

Fig. 3 The N = 0 (a), N = 1 (b), and N = 2 (c) energy levels of NH (³Σ⁻) plotted against the magnitude of the applied magnetic field. The initial state for scattering calculations is shown by the dashed line.

Fig. 4(a) and (b) show the adiabatic potentials for the Mg + NH complex. Because the Mg–NH interaction is dominated by the isotropic V₀(R) term in eqn (11), as shown in Fig. 4(c), the deepest adiabatic potentials closely follow the shape of V₀(R). The avoided crossings at short-range are mediated by the anisotropic part of the interaction potential. At long-range, the avoided crossings occur due to the different values of L in each adiabatic channel, resulting in long-range centrifugal barriers being shifted with respect to each other. Significantly, the presence of avoided crossings at both short- and long-range suggests that the non-adiabatic coupling terms in eqn (12) cannot be neglected at any value of R. The non-adiabatic couplings must be accounted for using, e.g., the diabatic-by-sector (DBS) method as described below.

Fig. 4 (a) Lowest adiabatic potentials of the Mg–NH collision complex. The inset zooms into the minimum of the six deepest adiabatic potentials to highlight the Zeeman splitting. (b) Medium- to long-range adiabatic potentials of the Mg–NH collision complex. All potentials have a well-defined value of M_tot = 1 and are calculated at a magnetic field of B = 1000 G. The zero of energy corresponds to the asymptotic value of the Mg–NH interaction potential (R → ∞). (c) Legendre expansion coefficients of the Mg–NH interaction potential.

Fig. 5 shows the results of our CC calculations for Mg + NH collisions. The cross sections computed using the adiabatic method are in excellent agreement with reference results obtained from the diabatic CC treatment,^42,43,86,87 even near the scattering resonance shown in Fig. 5(a). This shows that the fully converged adiabatic CC formulation is as accurate as the rigorous diabatic CC approach^44,45 when applied to cold and ultracold atom–molecule collisions in a magnetic field.

Fig. 5 Elastic (dashed lines) and total inelastic (solid lines) cross sections as a function of collision energy computed for ultracold Mg + NH collisions calculated via the adiabatic formulation using R-dependent truncation with a threshold of 0.001a₀⁻¹. Open shapes are cross sections calculated using the standard diabatic treatment with the uncoupled SF representation at magnetic field strengths 10 G (brown diamonds), 100 G (navy circles), 1000 G (teal squares). The upper panel (a) shows cross sections near a resonance at B = 100 G, and the lower panel (b) shows cross sections at ultracold energies.

We now turn to comparing the computational performance of the adiabatic basis to that of the standard diabatic basis. Because the computational time scales as the size of the basis cubed, O(M³), we compare the two approaches by reducing the number of channels used to propagate the log-derivative matrix and quantifying the error in the cross sections. To do this, we solve the adiabatic eigenvalue problem, eqn (6), at the midpoint of every sector (R = c_n) using the uncoupled-SF representation (eqn (9)) with basis functions including all primitive states with N ≤ 6, L ≤ 8, and M_tot = 1. We then reduce the number of adiabatic channels used to propagate the log-derivative matrix to a constant, R-independent value. For the diabatic basis, we simply reduce the number of diabatic channels, keeping it the same at all R, starting from the highest-energy, largest-N basis states.

Fig. 6 shows the detailed effect of reducing the adiabatic basis size in an R-independent way. To this end, we plot the ratios of the calculated inelastic cross sections to the benchmark values as a function of the basis size. The top panel of Fig. 6 shows results for a non-resonant collision energy E_col = 10⁻³ cm⁻¹, i.e., the energy away from the scattering resonances shown in Fig. 5(a). The bottom panel of Fig. 6 shows the results for the collision energy at the peak of the leftmost resonance feature shown in Fig. 5(a), which occurs at E_col = 0.033 cm⁻¹. We observe that at every basis size, the adiabatic basis produces cross sections with less error and provides high-quality approximate results with much fewer basis functions. Table 1 lists the number of basis functions required to compute the inelastic cross sections below a given percent error. The computational advantage of propagation in the adiabatic basis is approximated by cubing the ratio of the size of the diabatic basis (M^di) to the size of the adiabatic basis (M^ad). Table 1 shows that using the adiabatic basis results in the most significant advantage if larger percent errors can be tolerated.

Fig. 6 Ratios of the total inelastic cross sections calculated via a constant propagation basis size for the full propagation to the benchmark values computed with 954 basis states (N_max = 6, L_max = 8). (a) Off-resonant collision energy 10⁻³ cm⁻¹. (b) Resonant collision energy 0.033 cm⁻¹. Blue circles—adiabatic treatment, brown crosses—diabatic treatment. The vertical dashed lines represent cumulatively the number of states in each rotational manifold.

Table 1 The propagation basis size required to compute inelastic cross sections at different levels of accuracy in the adiabatic (M^ad) and diabatic (M^di) bases. The basis size was kept constant for the propagation

% Error	Ecol = 10^−3 cm^−1
	M^ad	M^di
0.1	540	720	2.4
1	380	540	2.9
5	220	480	10.4
10	200	380	6.9
50	80	200	15.6

---

As shown in Fig. 6(a), the diabatic basis provides results within ≃20% of the fully converged benchmark values with only three rotational manifolds in the basis (N_max = 2), and that the adiabatic basis can really only push about one rotational manifold lower to get similar accuracy.

Because the Mg–NH interaction is moderately anisotropic, one may wonder whether the superior performance of the adiabatic basis will hold for deeper, more anisotropic potentials that require far more rotational states for convergence.^22,23,36 To test this hypothesis, we performed the same convergence test as that shown in Fig. 6 replacing the original Mg–NH PES with a much deeper and more anisotropic PES generated by multiplying all Legendre expansion coefficients in eqn (11) by a factor to 10. To scale up the relative anisotropy, we scaled the anisotropic Legendre expansion coefficients (λ > 0 in eqn (11)) by a additional factor of two. The scaled potential V′(R) can then be written as

To accommodate such a deep and anisotropic PES, we solved the adiabatic eigenvalue problem in the total rotational angular momentum (TRAM) representation as our primitive basis^19,32,36,39 (see Appendix C and ref. 32 for computational details and matrix elements in the TRAM basis). Note that the uncoupled basis cannot provide converged results for such highly anisotropic PESs.^22,37 Approximately twice as many rotational manifolds are required to get converged cross sections using the scaled up PES (N_max = 12) compared to using the original Mg + NH PES (N_max = 6). Fig. 7 shows that the adiabatic basis requires 2–3 fewer rotational manifolds compared to the diabatic basis to compute cross sections that are within 20% of the benchmark cross sections, resulting in a substantial reduction in the number of adiabatic channels even for the much deeper and more anisotropic PES. This is in contrast to the 1–2 rotational manifolds that can be removed from the adiabatic basis when using the unscaled potential.

Fig. 7 Total inelastic (a) and elastic (b) cross sections for cold Mg + NH collisions calculated with a constant propagation basis size M and normalized to fully converged CC values for the scaled-up PES [eqn (22)]. The calculations were performed using the TRAM basis. The fully converged benchmark values were obtained with 1002 basis states (N_max = 15, J_r,max = 4). Panel (c) shows the elastic and total inelastic cross sections as a function of collision energy for the scaled-up PES.

Fig. 8 shows the Mg + NH cross sections computed using the adiabatic basis of a fixed size M. As for the data shown in Fig. 6, the fully converged primitive basis was used to solve the adiabatic eigenvalue problem, and the M channels with the lowest adiabatic eigenvalues were propagated. At collision energies away from the resonance [see Fig. 8(a)], the benchmark cross sections are well reproduced with just 200 adiabatic basis functions. Reducing the number of adiabatic potentials to 128 channels (i.e., the 24 open channels and 104 weakly closed channels) does not strongly reduce the accuracy of the elastic cross sections, which are remarkably well reproduced with only the 24 open channels (see Fig. 8) away from resonances. However, the inelastic cross sections with just the open channels are underestimated by a factor of ≃30.

Fig. 8 Cross sections for Mg + NH collisions vs. collision energy for different R-independent propagation basis sizes. (a) Elastic and total inelastic cross sections at ultracold temperatures. (b) Elastic cross sections near resonances. (c) Total inelastic cross sections near resonances. All calculations were done at a magnetic field strength of B = 100 G.

A total of 240 adiabatic basis functions are required to fully resolve the resonances with only a slight shift. With 200 basis functions, there is a noticeable shift in the resonance positions and several spurious narrow resonances appear. With 128 adiabatic basis functions, the main resonance features are replaced by narrow resonances. This shows the importance of the strongly closed adiabatic channels in the vicinity of scattering resonances. Still, the Wigner s-wave limiting values of elastic (inelastic) cross sections can be reproduced with only 24 (128) adiabatic channels, a factor of 39.7 (7.5) smaller than the full number of adiabatic states .

B. R-Dependent basis truncation (RBT)

Thus far we have used a fixed-size adiabatic basis at all R. We now invoke the RBT procedure described in Section IIC and analyze its performance. We set an RBT threshold that results in a small (<5%) error in the cross sections. Another key RBT parameter is the initial propagation basis size used before RBT is turned on. We choose the initial basis size that produces nearly-exact results near the resonance (Fig. 9), approximate results near the resonance (Fig. 10), and nearly-exact results at ultracold collision energies (Fig. 11).

Fig. 9 (a) Elastic and (c) total inelastic cross sections for cold Mg + NH collisions in a magnetic field of strength B = 100 G calculated using RBT with the R-dependent propagation basis size M_n shown in panel (b). The RBT thresholds used for the adiabatic and diabatic treatments are 10⁻²a₀⁻¹ and 10⁻⁴a₀⁻¹, respectively. The initial propagation basis size before RBT is turned is M₀ = 420 for the adiabatic treatment and M₀ = 532 for the diabatic treatment. (d) Computational gain of the adiabatic and diabatic RBT over the standard R-independent propagation as a function of R, where is the number of basis functions in the standard (fully uncoupled) basis. The horizontal dashed lines in (b) represent cumulatively the number of states in each rotational manifold. The gray circles are benchmark results obtained via the full uncoupled space-fixed basis in the standard diabatic treatment.

Fig. 10 Same as Fig. 9 with the same RBT thresholds but with initial basis sizes of M₀ = 240 and M₀ = 480 for the adiabatic and diabatic treatments, respectively.

Fig. 11 Same as Fig. 9 but for ultracold collision energies and with RBT thresholds of 10⁻¹a₀⁻¹ and 10⁻³a₀⁻¹ and initial basis sizes of M₀ = 200 and M₀ = 532 for the adiabatic and diabatic treatments, respectively.

The lower right panel of each figure shows the estimated computational gain offered by RBT over the standard CC calculation with basis functions, given by . The computational gan increases gradually with R, as the number of propagated channels drops. The most significant gains of 4–5 orders of magnitude occur at R ≥ 15a₀ in the adiabatic basis.

These results show that the adiabatic RBT procedure is more efficient than the diabatic one at all R. This is especially true at long range, where truncating to the open channels in adiabatic RBT is nearly 5 orders of magnitude more computationally efficient than standard CC, and 2–3 orders of magnitude more efficient than diabatic RBT. Because basis set truncation mostly occurs in the long-range region, it is the short-range region that dominates the computational time in RBT.

To estimate the overall computational time for the propagation of the log-derivative matrix, we calculated the quantity

where N_S is the total number of sectors and M_n is the number of basis functions used for the nth sector. We then take ratios of these sums to estimate the computational gains of adiabatic RBT (γ_aRBT) and diabatic RBT (γ_dRBT) over the benchmark fully converged CC calculation with 954 basis states (γ₉₅₄ = 954N_S). The ratio γ₉₅₄/γ_aRBT is 63.2, 306, and 539 for the basis size profiles from Fig. 9(b), 10(b), and 11(b), respectively. Comparing adiabatic RBT to diabatic RBT, the ratios γ_dRBT/γ_aRBT from the same three plots are 5.11, 18.3, and 27.4, respectively.

This suggests that adiabatic RBT offers the largest computational advantage when approximate cross sections (deviating slightly from the accurate resonance positions) can be tolerated, in keeping with the results from Table 1. We note that these estimates are only for the propagation part of the adiabatic CC calculation, and they do not account for the computational overhead associated with solving the adiabatic eigenvalue problem or the additional floating-point operations required for propagation in the adiabatic basis.

To address the latter, we measure the CPU time required to propagate the log-derivative for ten collision energies using three different approaches: the diabatic basis without RBT, the diabatic basis with RBT, and the adiabatic basis with RBT. For the calculations with RBT, we use the same RBT parameters used to produce the results in Fig. 10. The results are listed in Table 2. Adiabatic RBT is faster that the standard diabatic treatment by a factor of 56, and outperforms diabatic RBT by a factor of 3. However, it should be noted that our measurements of CPU time should still be regarded as estimates because of differences in implementation of propagation in the adiabatic vs. diabatic bases.

Table 2 Measurements of the CPU time needed to propagate the log-derivative matrix for different approaches. The RBT parameters, i.e., the initial basis size M₀ and the RBT threshold τ_RBT, are the same as those used to produce the results for Fig. 10. The values of the step size ΔR₁ were chosen from Fig. 12 to reduce the error from the R-grid to be less than 1%

Method	ΔR1 (a0)	M0	τRBT (a0^−1)	CPU time (s)
Diabatic	0.1	954	n/a	2140
Diabatic	0.1	480	10^−4	116
Adiabatic	0.02	240	10^−2	38.1

---

Unlike adiabatic RBT, diabatic RBT can be directly compared to the benchmark CC calculation with 954 diabatic basis functions. The ratio γ₉₅₄/γ_dRBT ranges from 12 to 20 [see Fig. 9–11], showing that diabatic RBT offers up to a 20-fold reduction in computational time when approximate cross sections can be tolerated. Table 2 shows that diabatic RBT is a factor of 18 times faster than the standard diabatic approach. The agreement between the CPU time measurements and our estimates using eqn (23) in the diabatic basis highlights the efficiency of diabatic RBT and validates our approach for estimating computational savings.

IV. Summary and conclusions

We have developed a rigorous adiabatic approach to ultracold atom–molecule collisions in an external magnetic field. The approach is based on propagating the log-derivative matrix across a grid of sectors. At each sector, one (1) solves the adiabatic eigenvalue problem at the midpoint of the sector, (2) propagates the log-derivative matrix across the sector, and (3) transforms the log-derivative matrix to the adiabatic basis of the next sector. Finally, at the end of the last sector, scattering boundary conditions are applied to the log-derivative matrix in the asymptotic basis. We applied the approach to ultracold Mg + NH collisions in a magnetic field, and found excellent agreement with standard CC calculations with a factor of ≃2 reduction in the number of scattering channels. We also developed an efficient log-derivative-based RBT method to reduce the size of the basis during propagation of the log-derivative matrix.

As shown in Section III for cold Mg + NH collisions, the adiabatic basis can be truncated significantly more aggressively than the diabatic basis, regardless of whether RBT is used, and especially when an appreciable error (<50%) in the calculated integral cross sections can be tolerated. That is, quality scattering observables can generally be obtained with far fewer basis functions in the adiabatic basis than the diabatic basis. In particular, by reducing the size of the basis without RBT, we found that the adiabatic approach can resolve the resonances near a collision energy of 0.04 cm⁻¹ with only 200 adiabatic basis functions, half of the basis functions required in the diabatic formulation. Away from resonances, the adiabatic approach needs only the 128 open and weakly-closed channels to get accurate inelastic cross sections and only the 24 open channels for the elastic cross sections (Fig. 8). By contrast, the diabatic approach needs at least 400 basis functions to obtain accurate cross sections either close to or away from resonances. This gain is expected because, unlike their diabatic counterparts, the adiabatic basis functions contain essential information about the scattering dynamics at fixed values of R. Although we showed that this result holds for a deeper and more anisotropic version of the Mg + NH PES, it would be interesting to explore this for a realistic collision, such as Rb + SrF.³⁶

Comparing the computational performance of the adiabatic vs. diabatic basis sets is challenging because the two approaches differ in the number of floating-point operations required per sector at the same basis size. We estimate that the adiabatic approach as presented here requires approximately 15 times more floating-point operations for the propagation part of the calculation than the standard diabatic method at the same basis size (see Appendix A for details). Considering this extra computational cost and using eqn (23) and the results in Fig. 10 and 11, we estimate the adiabatic basis is 15–30 times more computationally efficient than the standard diabatic basis with 954 basis functions, provided small shifts in resonance positions can be tolerated. Using the same metrics, the computational performance of the adiabatic RBT is on par with that of the diabatic RBT. This conclusion is supported by the measurements of CPU time shown in Table 2, which show a factor of 56 speed-up when using adiabatic RBT versus the standard diabatic approach.

Our results show that the adiabatic approach provides the largest computational advantage at long range because the adiabatic basis can be truncated down to just the open channels already at R ≥ 15a₀. Recent quantum reactive scattering calculations in the absence of external fields also reported a notable reduction in the number of adiabatic states (860 → 500) when projecting the wavefunction of the LiNa₂ reaction complex from adiabatically adjusting principal axis hyperspherical (APH) to Delves coordinates.⁹² The reduction observed in ref. 92 is, however, much less dramatic than that reported here (420 → 24), likely due to a large number of open channels in the ultracold Li + NaLi → Na + Li₂ reaction.

In notable contrast to the adiabatic basis, the diabatic basis does not produce accurate cross sections if truncated below about 90 basis functions due to the coupling between rotational states with different N in the asymptotic Hamiltonian (see Appendix B). The ability to use only the open channels to propagate the log-derivative matrix at long range is an attractive feature of the adiabatic approach because the number of open channels is generally much smaller than the total number of channels, and does not depend on the atom–molecule interaction anisotropy (or N_max). This also suggests that the adiabatic approach could be most advantageous for ultracold molecular collisions mediated by long-range interactions, such as NH + NH,^33,93 where the log-derivative matrix must be propagated out to very large atom–molecule separations (R_max = 500a₀³³).

The concept of universality in ultracold collision physics states that all of (or most of) the physics can be encapsulated in a small number of short-range parameters^42,94–96 that are independent of either the collision energy or magnetic field. These universal parameters are typically described in the framework of multichannel quantum defect theory, which leverages the separation of energy scales in ultracold collisions to arrive at a simplified description of two-body collision physics⁴² in terms of a few short-range parameters. However, this requires one to assume that some degrees of freedom (such as nuclear spins) play a spectator role. The recent experimental observation of hyperfine-to-rotational energy transfer in ultracold Rb + KRb collisions suggests that short-range couplings between the spin and rotational degrees of freedom can be non-negligible, posing a challenge for MQDT-FT.⁴²

Unlike MQDT-FT,⁴² our RBT approach does not rely on approximations. By systematically truncating the log-derivative matrix as it is propagated from small to large R and then stopping the propagation at an intermediate value of R_open, where the log-derivative matrix is maximally truncated to the open-channel basis, we obtain the log-derivative matrix Y(R_open), which contains all information about scattering observables. This provides an alternative, numerically exact way of condensing the complex physics of multichannel atom–molecule collisions into a small number of short-range parameters, the matrix elements of Y(R_open). For these parameters to be maximally useful, they should be independent of the collision energy and external fields. Whether or not this is the case remains to be explored. It would also be interesting to extend the efficient basis sets for solving the adiabatic eigenvalue problem (such as the TRAM representation³²) to include additional degrees of freedom, such as molecular vibrations and multiple potential energy surfaces.

Conflicts of interest

There are no conflicts to declare.

Data availability

Data for this article are available within the body of the article and in the form of Fig. 1–15 and Tables 1 and 2.

Appendices

Appendix A: Computational details

In all calculations, we propagate the log-derivative matrix from R = 4a₀ to R = 25a₀ with a step size of ΔR₁ = 0.01a₀. We then propagate to R = 100a₀, where we apply S-matrix boundary conditions. The step size for propagating the log-derivative matrix for R > 25a₀ is ΔR₂ = 0.1a₀. The values of parameters used in our CC calculations are listed in Table 3. We employ the same 954-channel uncoupled SF basis given by eqn (9) in benchmark diabatic CC calculations, as well as to solve the adiabatic eigenvalue problem (6). This basis includes all states with N ≤ 6, L ≤ 8, and M_tot = 1 and gives fully converged results. We choose the upper Zeeman state of the N = 0 rotational manifold as our initial state, and perform calculations for the M_tot = M_N + M_L + M_S = 1 block of states.

Table 3 Values of parameters used in our Mg + NH scattering calculations

Parameter	Value
Rotational constant of NH, Be	16.32176 cm^−1
Spin–rotation constant of NH, γsr	−0.05467 cm^−1
Spin–spin constant of NH, λSS	0.9197 cm^−1
Reduced mass of Mg–NH, μMg–NH	9.23267993 amu

---

The matrix operations involved in the propagation part of the calculation scale as the size of the basis cubed. There are two matrix inversions in the diabatic approach, whereas the adiabatic basis requires an orthogonal transformation in addition to the two inversions. Symmetric matrix inversion scales as ≃M³ and orthogonal transformations as ≃3M³.⁹⁷ From these rough estimates and a simple speed test, we find that it takes roughly two times more floating-point operations to propagate the log-derivative across a sector in the adiabatic basis than in the diabatic basis. In addition, we found that the diabatic-by-sector approach requires about 5–10 times more sectors than the diabatic method (Fig. 12). This amounts to an overall factor of ≃10–20 increase in the computational time required to propagate the log-derivative matrix in the adiabatic basis compared to the diabatic basis. We take the midpoint of this range, and conclude that the adiabatic approach requires ≃15 times more floating-point operations to propagate the log-derivative through a range of R-values for the same basis size compared to the diabatic approach. Fortunately, this increase is more than offset by the more efficient performance of RBT in the adiabatic basis, as described below.

Fig. 12 Percent error in the inelastic cross section plotted as a function of the step size ΔR₁ with ΔR₂ kept fixed (see text). The total inelastic cross section with a grid spacing of ΔR₁ = 0.001a₀ is used as the benchmark value in computing the error values.

Appendix B: Truncating the basis during propagation

Fig. 13(a) shows the number of basis states used for propagation as a function of R, for log-derivative-based RBT in the adiabatic basis. These RBT “trajectories,” or basis size profiles, are plotted for different RBT thresholds that produce cross sections of varying accuracy. The corresponding errors introduced by RBT are shown in Fig. 13(b) and (c) for resonant and non-resonant collision energies, respectively. Analyzing the RBT errors gives insight into which adiabatic states are important in different regions of R. Compare, for example, the basis size profile for τ_RBT = 10⁻²a₀⁻¹ (green curve), which introduces a very small (<0.1%) error with that for τ_RBT = 0.63a₀⁻¹ (orange curve), for which the error is large (≃100%). These two basis size profiles mainly differ in that the green curve retains ≈400 adiabatic channels until R = 10a₀ and the adiabatic channels in the N = 1 manifold until R = 15a₀. The orange curve truncates the adiabatic basis to just the channels in the N = 0 and N = 1 manifolds by R ≈ 8a₀, and just the open channels by R ≈ 13a₀. The significant difference in accuracy of cross sections resulting from these two basis size profiles indicates the importance of the additional adiabatic channels retained in the basis for τ_RBT = 10⁻²a₀⁻¹ in the two regions of R discussed above.

Fig. 13 Log-derivative-based RBT in the adiabatic basis. (a) Basis size M_n against atom–molecule separation R for three different RBT thresholds. The horizontal dashed lines indicate the cumulative number of channels in each rotational manifold up to N = 6. (b) Percent error in the elastic (filled red circles) and total inelastic (open blue circles) cross sections against RBT threshold at E_col = 0.033 cm⁻¹ and B = 100 G. The vertical dashed lines show the RBT thresholds for which we plot the basis size profiles in (a). The horizontal dashed line marks 1% error as a guide to the eye. (c) Same as (b) but for E_col = 10⁻³ cm⁻¹.

We also tested the performance of overlap-based RBT, in which the matrix elements of the overlap matrix (eqn (20)) are sampled to determine the channel to be removed. We found the best implementation to be close to that proposed by Stetchel, Walker, and Light.⁶¹ They suggest removing the αth channel if the quantity

is less than the user-defined RBT threshold τ_RBT, where α is the index of the channel being tested and O_jα are matrix elements of the overlap matrix (eqn (20)). As in log-derivative-based RBT in the adiabatic basis, we start the truncation from the adiabatic channels with the highest energy, and continue in the order of decreasing energy. We stop truncating once we find the first state that does not satisfy β_α < τ_RBT. Fig. 14 shows the basis size profiles and errors in inelastic cross sections with respect to τ_RBT for overlap-based RBT. An interesting aspect of overlap-based RBT is that the error in the cross sections changes rapidly from ≈10⁻²% to ≈100% when τ_RBT is varied slightly from τ_RBT = 0.004a₀⁻¹ to 0.005a₀⁻¹. Comparing the basis size profiles in Fig. 13(a) and 14(a) for similar levels of accuracy, we find that log-derivative-based RBT produces basis size profiles that deviate only slightly from those of overlap-based RBT. We note that basis set profiles in overlap-based RBT can be obtained without propagating the log-derivative matrix.

Fig. 14 Same as Fig. 13 but for overlap-based RBT, which is necessarily in the adiabatic basis.

The solid lines in Fig. 13(a) differ from the dashed lines in that the basis size profiles were produced by sampling the log-derivative matrix at a different collision energy. Comparing these shows that the basis size profiles are largely independent of collision energy. This makes implementation easier because the log-derivative matrix needs only to be sampled at a single collision energy, and the resulting basis size profiles can be used for other collision energies (this was only tested for collision energies below 0.1 cm⁻¹). The overlap matrix is independent of the collision energy, and so are overlap-based RBT basis size profiles.

Diabatic RBT. A major advantage of log-derivative-based RBT over overlap-based RBT is that the former can be applied in the diabatic basis. Below we describe the main features of RBT in the diabatic basis (diabatic RBT), including the procedures for (i) how to define the locally open and locally closed channels, (ii) how to search for the states to be removed, and (iii) when to stop the search.

First, in diabatic RBT, we label a basis state as locally open or locally closed based on the diagonal matrix element of the reference potential, defined by eqn (16) with the adiabatic eigenvector matrix T_n replaced by the identity matrix. The ith primitive basis state is locally open if [W_ref]_ii < 0 and locally closed if [W_ref]_ii > 0.

The performance of diabatic RBT is illustrated in Fig. 15. As expected, the error in the cross sections steadily decreases as the RBT threshold is decreased. For a given accuracy, the basis size obtained with diabatic RBT is larger than that with adiabatic RBT at almost all R values. This is consistent with the results shown in Fig. 9–11 in the main text. We can also quantify the computational gain provided by diabatic RBT by using eqn (23) and calculating the ratio γ₉₅₄/γ_dRBT. From these ratios, the green and blue basis size profiles in Fig. 15 (both of which provide accurate cross sections) correspond to 4- and 5.5-fold computational gains over the standard diabatic CC basis consisting of 954 functions.

Fig. 15 Same as Fig. 13 but for diabatic RBT.

As noted above, in adiabatic RBT we start channel elimination from the highest-energy adiabatic state in each sector and continue until we encounter a channel that does not meet the cutoff condition. The entire procedure is then stopped, a procedure we refer to as the “top-down approach”. In diabatic RBT, an alternative “any-state” approach is possible, where any locally closed channel that satisfies the cutoff condition can be removed, regardless whether or not all the higher-energy channels have already been eliminated. The any-state approach is thus expected to trim the diabatic basis more aggressively as a function of R.

Fig. 16(a) shows the basis size profile for any-state diabatic RBT. Also plotted are markers (orange points) showing the indices of diabatic basis functions retained in the basis. The inset highlights how any-state diabatic RBT is able to remove all of the N = 1 basis functions while keeping the N = 2 functions, resulting in just 96 basis functions. This is a factor of 2 fewer channels than with the top-down approach, which truncates to ≃190 basis functions, retaining almost all states up to and including N = 2 functions.

Fig. 16 (a) Basis size profile (black line) for diabatic RBT at E_col = 10⁻³ cm⁻¹, B = 100 G, and an RBT threshold of 10⁻³a₀⁻¹. At each R-value, the brown points mark the indices of the diabatic basis functions retained in the basis. The inset shows the region R ∈ [18,23]a₀ to highlight the truncation of the diabatic basis functions with N = 1 while retaining the functions with N = 2. (b) The basis size profiles for adiabatic (blue) and diabatic (brown) RBT at E_col = 10⁻³ cm⁻¹ and an RBT threshold of 10⁻³a₀⁻¹. The dashed lines correspond to calculations with the spin–spin interaction in the asymptotic Hamiltonian turned off and the solid lines correspond to full CC results with the spin–spin interaction included.

The diabatic basis must retain the N = 2 primitive basis functions at long range because of the spin–spin interaction in the asymptotic Hamiltonian (4), which couples the N = 0 states that dominate the open channels to the N = 2 primitive states.⁴⁵ Fig. 16(b) compares the basis size profiles with the spin–spin interaction included (thin solid curves) vs. omitted (thick dashed curves). A significant difference occurs at R = 21a₀, where the basis size profile with the spin–spin interaction omitted truncates to just the 24 open channels, whereas the full calculation does not truncate below 96 diabatic channels. The basis size profile for adiabatic RBT is not significantly affected by removing the spin–spin interaction because this coupling is already incorporated in the adiabatic basis functions.

RBT parameters. To determine the RBT parameters—the initial basis size M₀ and the RBT threshold τ_RBT—it is instructive to highlight two features of the RBT procedure.

First, RBT is highly robust with respect to the choice of threshold. Fig. 13–15 show that a broad range of RBT thresholds yields both accurate cross sections and substantial basis truncation. As a result, the threshold can be chosen conservatively without sacrificing accuracy. For example, exploratory scans over collision energy or other parameters can be performed using a relatively low threshold, such as τ_RBT = 0.0016a₀⁻¹ for overlap-based RBT, for which one can safely expect accurate cross sections in addition to significant computational savings.

Second, the resulting basis size profiles are largely independent of the initial basis size. Fig. 17 shows the basis size profiles obtained using the same RBT threshold but different initial basis sizes; once the profiles overlap, they are nearly identical. Because the RBT threshold is easy to choose and the two parameters are effectively independent, we found it optimal to use the initial basis size as the primary parameter that controls the accuracy of RBT calculations. This parameter can be determined using a standard basis convergence test as described above [see Fig. 8].

Fig. 17 Basis size profiles computed using log-derivative-based RBT in the adiabatic basis for different initial basis sizes. All other parameters are fixed.

Appendix C: Matrix elements of the spin–spin interaction in the TRAM basis

Here, we provide the expressions for the matrix elements of the adiabatic Hamiltonian in the TRAM basis. The TRAM is a vector sum of angular momenta for mechanical rotation, Ĵ_r = N̂ + L̂. The TRAM basis functions used as our primitive basis are eigenstates of Ĵ_r², Ĵ_z, N̂², and L̂², and can be written in terms on the uncoupled states |NM_N〉 |LM_L〉 asThe matrix elements for all interactions in the adiabatic Hamiltonian [eqn (2)] for collisions of ¹S-atoms with ³Σ-molecules can be found in ref. 32 except for the spin–spin interaction given by⁴⁵where λ_SS is the spin–spin interaction constant, Y_lm(R̂) are spherical harmonics, and [Ŝ ⊗ Ŝ]_q⁽²⁾ is a tensor product. The matrix elements in the TRAM basis areThe matrix elements of the spherical harmonics in eqn (C4) in the TRAM basis are³²where the reduced matrix element is given byand the terms in the parenthesis and curly brackets are 3j and 6j symbols, respectively. The matrix elements of the tensor product of the electronic spin with itself are⁴⁵The reduced matrix element isFor S = 1, the 6j-symbol is equal to 1/6, and the reduced matrix element is .⁴⁵ Combining eqn (C5)–(C9) with eqn (C4), we arrive at the final expression for the matrix elements of the spin–spin interaction in the TRAM basis:

Acknowledgements

We are grateful to Chris Greene and Jose P. D'Incao for stimulating discussions. This work was supported by the NSF CAREER award no. PHY-2045681.

References

1. S. F. Yelin, K. Kirby and R. Côté, Schemes for robust quantum computation with polar molecules, Phys. Rev. A:At., Mol., Opt. Phys., 2006, 74, 050301.
2. A. V. Gorshkov, S. R. Manmana, G. Chen, E. Demler, M. D. Lukin and A. M. Rey, Quantum magnetism with polar alkali-metal dimers, Phys. Rev. A:At., Mol., Opt. Phys., 2011, 84, 033619.
3. V. V. Albert, J. P. Covey and J. Preskill, Robust Encoding of a Qubit in a Molecule, Phys. Rev. X, 2020, 10, 031050.
4. S. L. Cornish, M. R. Tarbutt and K. R. A. Hazzard, Quantum computation and quantum simulation with ultracold molecules, Nat. Phys., 2024, 20, 730–740.
5. R. V. Krems, Cold Controlled Chemistry, Phys. Chem. Chem. Phys., 2008, 10, 4079–4092.
6. N. Balakrishnan, Perspective: Ultracold molecules and the dawn of cold controlled chemistry, J. Chem. Phys., 2016, 145, 150901.
7. R. V. Krems, Molecules in Electromagnetic Fields, Wiley-VCH, 2019.
8. Y. Liu and K. K. Ni, Bimolecular Chemistry in the Ultracold Regime, Annu. Rev. Phys. Chem., 2022, 73, 73–96.
9. T. Karman, M. Tomza and J. Pérez-Ríos, Ultracold chemistry as a testbed for few-body physics, Nat. Phys., 2024, 20, 722–729.
10. L. D. Carr, D. DeMille, R. V. Krems and J. Ye, Cold and ultracold molecules: science, technology and applications, New J. Phys., 2009, 11, 055049.
11. J. L. Bohn, A. M. Rey and J. Ye, Cold molecules: Progress in quantum engineering of chemistry and quantum matter, Science, 2017, 357, 1002–1010.
12. D. DeMille, N. R. Hutzler, A. M. Rey and T. Zelevinsky, Quantum sensing and metrology for fundamental physics with molecules, Nat. Phys., 2024, 20, 741–749.
13. T. P. Softley, Cold and ultracold molecules in the twenties, Proc. R. Soc. A, 2023, 479, 20220806.
14. T. Langen, G. Valtolina, D. Wang and J. Ye, Quantum state manipulation and cooling of ultracold molecules, Nat. Phys., 2024, 20, 702–712.
15. Y. Liu, M. G. Hu, M. A. Nichols, D. Yang, D. Xie, H. Guo and K.-K. Ni, Precision test of statistical dynamics with state-to-state ultracold chemistry, Nature, 2021, 593, 379–384.
16. J. J. Park, Y. K. Lu, A. O. Jamison, T. V. Tscherbul and W. Ketterle, A Feshbach resonance in collisions between triplet ground-state molecules, Nature, 2023, 614, 54–58.
17. H. Yang, J. Cao, Z. Su, J. Rui, B. Zhao and J. W. Pan, Creation of an ultracold gas of triatomic molecules from an atom-diatomic molecule mixture, Science, 2022, 378, 1009–1013.
18. S. Haze, J. L. Li, D. Dorer, J. P. D’Incao, P. S. Julienne, E. Tiemann, M. Deiß and J. Hecker Denschlag, Controlling few-body reaction pathways using a Feshbach resonance, Nat. Phys., 2025, 21, 228–232.
19. Y. X. Liu, L. Zhu, J. Luke, M. C. Babin, M. Gronowski, H. Ladjimi, M. Tomza, J. L. Bohn, T. V. Tscherbul and K.-K. Ni, Hyperfine-to-rotational energy transfer in ultracold atom–molecule collisions of Rb and KRb, Nat. Chem., 2025, 17, 688–694.
20. N. Bigagli, W. Yuan, S. Zhang, B. Bulatovic, T. Karman, I. Stevenson and S. Will, Observation of Bose-Einstein condensation of dipolar molecules, Nature, 2024, 631, 289–293.
21. R. Bause, A. Christianen, A. Schindewolf, I. Bloch and X. Y. Luo, Ultracold Sticky Collisions: Theoretical and Experimental Status, J. Phys. Chem. A, 2023, 127, 729.
22. T. V. Tscherbul, J. Kłos and A. A. Buchachenko, Ultracold spin-polarized mixtures of ²Σ molecules with S-state atoms: Collisional stability and implications for sympathetic cooling, Phys. Rev. A:At., Mol., Opt. Phys., 2011, 84, 040701.
23. M. Morita, M. B. Kosicki, P. S. Żuchowski and T. V. Tscherbul, Atom-molecule collisions, spin relaxation, and sympathetic cooling in an ultracold spin-polarized Rb(²S)–SrF(²Σ⁺) mixture, Phys. Rev. A:At., Mol., Opt. Phys., 2018, 98, 042702.
24. M. Morita, J. Kłos, A. A. Buchachenko and T. V. Tscherbul, Cold collisions of heavy ²Σ molecules with alkali-metal atoms in a magnetic field: Ab initio analysis and prospects for sympathetic cooling of SrOH(²Σ⁺) by Li(²S), Phys. Rev. A, 2017, 95, 063421.
25. C. Chin, R. Grimm, P. Julienne and E. Tiesinga, Feshbach resonances in ultracold gases, Rev. Mod. Phys., 2010, 82, 1225–1286.
26. H. Yang, D. C. Zhang, L. Liu, Y. X. Liu, J. Nan, B. Zhao and J.-W. Pan, Observation of magnetically tunable Feshbach resonances in ultracold ²³Na⁴⁰K + ⁴⁰K collisions, Science, 2019, 363, 261.
27. X. Y. Wang, M. D. Frye, Z. Su, J. Cao, L. Liu, D. C. Zhang, H. Yang, J. M. Hutson, B. Zhao, C.-L. Bai and J.-W. Pan, Magnetic Feshbach resonances in collisions of ²³Na⁴⁰K with ⁴⁰K, New J. Phys., 2021, 23, 115010.
28. J. J. Park, H. Son, Y. K. Lu, T. Karman, M. Gronowski, M. Tomza, A. O. Jamison and W. Ketterle, Spectrum of Feshbach Resonances in NaLi + Na Collisions, Phys. Rev. X, 2023, 13, 031018.
29. H. Son, J. J. Park, Y. K. Lu, A. O. Jamison, T. Karman and W. Ketterle, Control of reactive collisions by quantum interference, Science, 2022, 375, 1006–1010.
30. M. Morita, R. V. Krems and T. V. Tscherbul, Universal Probability Distributions of Scattering Observables in Ultracold Molecular Collisions, Phys. Rev. Lett., 2019, 123, 013401.
31. T. Karman, M. Gronowski, M. Tomza, J. J. Park, H. Son, Y. K. Lu, A. O. Jamison and W. Ketterle, Ab initio calculation of the spectrum of Feshbach resonances in NaLi + Na collisions, Phys. Rev. A, 2023, 108, 023309.
32. T. V. Tscherbul and J. P. D'Incao, Ultracold molecular collisions in magnetic fields: Efficient incorporation of hyperfine structure in the total rotational angular momentum representation, Phys. Rev. A, 2023, 108, 053317.
33. Y. V. Suleimanov, T. V. Tscherbul and R. V. Krems, Efficient method for quantum calculations of molecule–molecule scattering properties in a magnetic field, J. Chem. Phys., 2012, 137, 024103.
34. M. S. Child, Molecular Collision Theory, Wiley, 1974.
35. S. C. Althorpe and D. C. Clary, Quantum Scattering Calculations on Chemical Reactions, Annu. Rev. Phys. Chem., 2003, 54, 493–529.
36. M. Morita, M. B. Kosicki, P. S. Żuchowski, P. Brumer and T. V. Tscherbul, Magnetic Feshbach resonances in ultracold atom-molecule collisions, Phys. Rev. A, 2024, 110, L021301.
37. T. V. Tscherbul and A. Dalgarno, Quantum theory of molecular collisions in a magnetic field: Efficient calculations based on the total angular momentum representation, J. Chem. Phys., 2010, 133, 184104.
38. T. V. Tscherbul, Total-angular-momentum representation for atom-molecule collisions in electric fields, Phys. Rev. A:At., Mol., Opt. Phys., 2012, 85, 052710.
39. A. Simoni and J. M. Launay, Ultracold atom-molecule collisions with hyperfine coupling, Laser Phys., 2006, 16, 707–712.
40. R. Chapurin, X. Xie, M. J. Van de Graaff, J. S. Popowski, J. P. D'Incao, P. S. Julienne, J. Ye and E. A. Cornell, Precision Test of the Limits to Universality in Few-Body Physics, Phys. Rev. Lett., 2019, 123, 233402.
41. X. Xie, M. J. Van de Graaff, R. Chapurin, M. D. Frye, J. M. Hutson, J. P. D'Incao, P. S. Julienne, J. Ye and E. A. Cornell, Observation of Efimov Universality across a Nonuniversal Feshbach Resonance in ³⁹K, Phys. Rev. Lett., 2020, 125, 243401.
42. M. Morita, P. Brumer and T. V. Tscherbul, Multichannel Quantum Defect Theory with a Frame Transformation for Ultracold Atom-Molecule Collisions in Magnetic Fields, Phys. Rev. Lett., 2024, 133, 093402.
43. D. Vieira, R. V. Krems and T. V. Tscherbul, Molecular collisions and reactive scattering in external fields: Are field-induced couplings important at short range?, J. Chem. Phys., 2017, 146, 024102.
44. A. Volpi and J. L. Bohn, Magnetic-field effects in ultracold molecular collisions, Phys. Rev. A:At., Mol., Opt. Phys., 2002, 65, 052712.
45. R. V. Krems and A. Dalgarno, Quantum-mechanical theory of atom-molecule and molecular collisions in a magnetic field: Spin depolarization, J. Chem. Phys., 2004, 120, 2296–2307.
46. T. V. Tscherbul and R. V. Krems, Controlling Electronic Spin Relaxation of Cold Molecules with Electric Fields, Phys. Rev. Lett., 2006, 97, 083201.
47. M. L. González-Martínez and J. M. Hutson, Ultracold atom-molecule collisions and bound states in magnetic fields: Tuning zero-energy Feshbach resonances in He–NH (³Σ⁻), Phys. Rev. A:At., Mol., Opt. Phys., 2007, 75, 022702.
48. T. V. Tscherbul, G. C. Groenenboom, R. V. Krems and A. Dalgarno, Dynamics of OH(²Π)–He collisions in combined electric and magnetic fields, Faraday Discuss., 2009, 142, 127.
49. T. V. Tscherbul, Y. V. Suleimanov, V. Aquilanti and R. V. Krems, Magnetic field modification of ultracold molecule–molecule collisions, New J. Phys., 2009, 11, 055021.
50. Z. Pavlovic, T. V. Tscherbul, H. R. Sadeghpour, G. C. Groenenboom and A. Dalgarno, Cold Collisions of OH(²Π) Molecules with He Atoms in External Fields, J. Phys. Chem. A, 2009, 113, 14670–14680.
51. M. T. Hummon, T. V. Tscherbul, J. Kłos, H. I. Lu, E. Tsikata, W. C. Campbell, A. Dalgarno and J. M. Doyle, Cold N + NH Collisions in a Magnetic Trap, Phys. Rev. Lett., 2011, 106, 053201.
52. T. V. Tscherbul and A. Dalgarno, Quantum theory of molecular collisions in a magnetic field: Efficient calculations based on the total angular momentum representation, J. Chem. Phys., 2010, 133, 184104.
53. S. Koyu, R. Hermsmeier and T. V. Tscherbul, Total angular momentum representation for state-to-state quantum scattering of cold molecules in a magnetic field, J. Chem. Phys., 2022, 156, 034112.
54. C. H. Greene, P. Giannakeas and J. Pérez-Ríos, Universal few-body physics and cluster formation, Rev. Mod. Phys., 2017, 89, 035006.
55. R. T. Pack and G. A. Parker, Quantum reactive scattering in three dimensions using hyperspherical (APH) coordinates. Theory, J. Chem. Phys., 1987, 87, 3888–3921.
56. D. Skouteris, J. F. Castillo and D. E. Manolopoulos, ABC: a quantum reactive scattering program, Comput. Phys. Commun., 2000, 133, 128–135.
57. J. von Stecher, J. P. D’Incao and C. H. Greene, Signatures of universal four-body phenomena and their relation to the Efimov effect, Nat. Phys., 2009, 122, 417–421.
58. M. D. Higgins and C. H. Greene, Five-body recombination of identical bosons, Proc. Natl. Acad. Sci. U. S. A., 2025, 122, e2503390122.
59. H. J. da Silva, B. K. Kendrick, H. Li, S. Kotochigova and N. Balakrishnan, Nonadiabatically Driven Quantum Interference Effects in the Ultracold K + KRb → Rb + K2 Chemical Reaction, J. Phys. Chem. Lett., 2025, 16, 6171–6177.
60. J. C. Light and R. B. Walker, An R matrix approach to the solution of coupled equations for atom–molecule reactive scattering, J. Chem. Phys., 1976, 65, 4272–4282.
61. E. B. Stechel, R. B. Walker and J. C. Light, R-matrix solution of coupled equations for inelastic scattering, J. Chem. Phys., 1978, 69, 3518–3531.
62. J. M. Hutson and C. R. Le Sueur, MOLSCAT: a program for non-reactive quantum scattering calculations on atomic and molecular collisions, Comput. Phys. Commun., 2019, 241, 9–18.
63. R. Hermsmeier, J. Kłos, S. Kotochigova and T. V. Tscherbul, Quantum Spin State Selectivity and Magnetic Tuning of Ultracold Chemical Reactions of Triplet Alkali-Metal Dimers with Alkali-Metal Atoms, Phys. Rev. Lett., 2021, 127, 103402.
64. L. Wiesenfeld, Quantum nature of molecular vibrational quenching: Water-molecular hydrogen collisions, J. Chem. Phys., 2021, 155, 071104.
65. D. C. Clary, Fast Chemical Reactions: Theory Challenges Experiment, Annu. Rev. Phys. Chem., 1990, 41, 61–90.
66. D. C. Clary, Effect of rotational excitation on chemical reaction cross sections, Mol. Phys., 1983, 48, 619–629.
67. E. J. Rackham, T. Gonzalez-Lezana and D. E. Manolopoulos, A rigorous test of the statistical model for atom-diatom insertion reactions, J. Chem. Phys., 2003, 119, 12895–12907.
68. M. H. Alexander, E. J. Rackham and D. E. Manolopoulos, Product multiplet branching in the O(¹D) + H₂ → OH(²Π) + H reaction, J. Chem. Phys., 2004, 121, 5221–5235.
69. T. V. Tscherbul and A. A. Buchachenko, Adiabatic channel capture theory applied to cold atom-molecule reactions: Li + CaH → LiH + Ca at 1 K, New J. Phys., 2015, 17, 035010.
70. Y. Scribano, G. Parlant and B. Poirier, Communication: Adiabatic quantum trajectory capture for cold and ultra-cold chemical reactions, J. Chem. Phys., 2018, 149, 021101.
71. M. Pawlak, Y. Shagam, E. Narevicius and N. Moiseyev, Adiabatic theory for anisotropic cold molecule collisions, J. Chem. Phys., 2015, 143, 074114.
72. M. Pawlak, Y. Shagam, A. Klein, E. Narevicius and N. Moiseyev, Adiabatic Variational Theory for Cold Atom–Molecule Collisions: Application to a Metastable Helium Atom Colliding with ortho- and para-Hydrogen Molecules, J. Phys. Chem. A, 2017, 121, 2194–2198.
73. J. Troe, Statistical adiabatic channel model of ion-neutral dipole capture rate constants, Chem. Phys. Lett., 1985, 122, 425–430.
74. M. Quack and J. Troe, Specific rate constants of unimolecular processes II. Adiabatic channel model, Ber. Bunsen-Ges. Phys. Chem., 1974, 78, 240–252.
75. M. Quack and J. Troe, Complex formation in reactive and inelastic scattering: Statistical adiabatic channel model of unimolecular processes III, Ber. Bunsen-Ges. Phys. Chem., 1975, 79, 170–183.
76. K. Matsuda, L. D. Marco, J. R. Li, W. G. Tobias, G. Valtolina, G. Quéméner and J. Ye, Resonant collisional shielding of reactive molecules using electric fields, Science, 2020, 370, 1324–1327.
77. L. Anderegg, S. Burchesky, Y. Bao, S. S. Yu, T. Karman, E. Chae, K.-K. Ni, W. Ketterle and J. M. Doyle, Observation of microwave shielding of ultracold molecules, Science, 2021, 373, 779–782.
78. A. V. Gorshkov, P. Rabl, G. Pupillo, A. Micheli, P. Zoller, M. D. Lukin and H. P. Büchler, Suppression of Inelastic Collisions Between Polar Molecules With a Repulsive Shield, Phys. Rev. Lett., 2008, 101, 073201.
79. L. Lassablière and G. Quéméner, Controlling the Scattering Length of Ultracold Dipolar Molecules, Phys. Rev. Lett., 2018, 121, 163402.
80. T. Karman and J. M. Hutson, Microwave Shielding of Ultracold Polar Molecules, Phys. Rev. Lett., 2018, 121, 163401.
81. E. J. Rackham, F. Huarte-Larranaga and D. E. Manolopoulos, Coupled-channel statistical theory of the N(²D) + H₂ and O(¹D) + H₂ insertion reactions, Chem. Phys. Lett., 2001, 343, 356–364.
82. T. V. Tscherbul and R. V. Krems, Rigorous quantum calculations for atom–molecule chemical reactions in electric fields: From single to multiple partial wave regimes, J. Chem. Phys., 2025, 163, 184101.
83. J. Wang, J. P. D'Incao and C. H. Greene, Numerical study of three-body recombination for systems with many bound states, Phys. Rev. A:At., Mol., Opt. Phys., 2011, 84, 052721.
84. B. R. Johnson, The multichannel Log-Derivative Method for Scattering Calculations, J. Comput. Phys., 1973, 13, 445.
85. D. E. Manolopoulos, An improved log derivative method for inelastic scattering, J. Chem. Phys., 1986, 85, 6425.
86. A. O. G. Wallis and J. M. Hutson, Production of Ultracold NH Molecules by Sympathetic Cooling with Mg, Phys. Rev. Lett., 2009, 103, 183201.
87. M. L. González-Martínez and J. M. Hutson, Effect of hyperfine interactions on ultracold molecular collisions: NH(³Σ⁻) with Mg(¹S) in magnetic fields, Phys. Rev. A:At., Mol., Opt. Phys., 2011, 84, 052706.
88. J. Brown and A. Carrington, Rotational Spectroscopy of Diatomic Molecules, Cambridge University Press, 2003.
89. R. V. Krems, H. R. Sadeghpour, A. Dalgarno, D. Zgid, J. Kłos and G. Chałasiński, Low-temperature collisions of NH(X³Σ⁻) molecules with He atoms in a magnetic field: An ab initio study, Phys. Rev. A:At., Mol., Opt. Phys., 2003, 68, 051401.
90. W. C. Campbell, E. Tsikata, H. I. Lu, L. D. van Buuren and J. M. Doyle, Magnetic Trapping and Zeeman Relaxation of NH (X³Σ⁻), Phys. Rev. Lett., 2007, 98, 213001.
91. W. C. Campbell, T. V. Tscherbul, H. I. Lu, E. Tsikata, R. V. Krems and J. M. Doyle, Mechanism of Collisional Spin Relaxation in ³Σ Molecules, Phys. Rev. Lett., 2009, 102, 013003.
92. B. K. Kendrick, H. Li, M. Li, S. Kotochigova, J. F. E. Croft and N. Balakrishnan, Non-adiabatic quantum interference in the ultracold Li + LiNa → Li₂ + Na reaction, Phys. Chem. Chem. Phys., 2021, 23, 5096–5112.
93. Y. V. Suleimanov and T. V. Tscherbul, Cold NH-NH collisions in a magnetic field: Basis set convergence versus sensitivity to the interaction potential, J. Phys. B:At., Mol. Opt. Phys., 2016, 49, 204002.
94. J. P. Burke, C. H. Greene and J. L. Bohn, Multichannel Cold Collisions: Simple Dependences on Energy and Magnetic Field, Phys. Rev. Lett., 1998, 81, 3355–3358.
95. B. Gao, E. Tiesinga, C. J. Williams and P. S. Julienne, Multichannel quantum-defect theory for slow atomic collisions, Phys. Rev. A:At., Mol., Opt. Phys., 2005, 72, 042719.
96. J. F. E. Croft, A. O. G. Wallis, J. M. Hutson and P. S. Julienne, Multichannel Quantum Defect Theory for cold molecular collisions, Phys. Rev. A:At., Mol., Opt. Phys., 2011, 84, 042703.
97. L. N. Trefethen and D. Bau, Numerical linear algebra, Society for Industrial and Applied Mathematics, Philadelphia, 1997.

---