Chemical Reactions Using a Non-Equilibrium Wigner Function Approach

A three-dimensional model of binary chemical reactions is studied. We consider an ab initio quantum two-particle system subjected to an attractive interaction potential and to a heat bath at thermal equilibrium at absolute temperature T > 0. Under the sole action of the attraction potential, the two particles can either be bound or unbound to each other. While at T = 0, there is no transition between both states, such a transition is possible when T > 0 (due to the heat bath) and plays a key role as kBT approaches the magnitude of the attractive potential. We focus on a quantum regime, typical of chemical reactions, such that: (a) the thermal wavelength is shorter than the range of the attractive potential (lower limit on T) and (b) (3/2)kBT does not exceed the magnitude of the attractive potential (upper limit on T). In this regime, we extend several methods previously applied to analyze the time duration of DNA thermal denaturation. The two-particle system is then described by a non-equilibrium Wigner function. Under Assumptions (a) and (b), and for sufficiently long times, defined by a characteristic time scale D that is subsequently estimated, the general dissipationless non-equilibrium equation for the Wigner function is approximated by a Smoluchowski-like equation displaying dissipation and quantum effects. A comparison with the standard chemical kinetic equations is made. The time τ required for the two particles to transition from the bound state to unbound configurations is studied by means of the mean first passage time formalism. An approximate formula for τ, in terms of D and exhibiting the Arrhenius exponential factor, is obtained. Recombination processes are also briefly studied within our framework and compared with previous well-known methods.


Introduction
While the foundations of equilibrium statistical mechanics are well established [1][2][3][4][5], those of non-equilibrium statistical mechanics continue to be an open and active important research area [1,[6][7][8].The monographs by [9][10][11] provide a wide and accessible perspective on classical equilibrium and irreversible thermodynamics.Within the realm of quantum thermodynamics (see [12] and the references therein), where the laws of thermodynamics from quantum mechanics and the theory of open quantum systems remain a mainstay [13], particular attention has been put to the Markovian assumptions and the dissipative dynamics proposed by Lindblad and Gorini-Kossakowski-Sudarshan [14,15].One crucial hallmark is that the off-equilibrium evolution of statistical systems displays stochasticity: see [16][17][18] in the classical regime and [13,[19][20][21] in the quantum one.
It is the goal of the present study to focus on non-equilibrium statistical mechanics in the quantum regime, which faces many open fundamental problems.A key question is how to formulate suitably approximate dynamical equations, displaying stochasticity, by starting out from a quantum mechanical framework and carrying out adequate long time approximations.Such a task turns out to be very difficult for isolated and large nonequilibrium quantum systems.
Fortunately, certain rewarding results have been obtained for a closed composite system consisting of the system S int of actual interest and a larger "heat bath" (HB) in quantum mechanics, under various conditions.Those results required considering the quantum states of S int and of the HB altogether and various conditions on Hamiltonians and that the states of S int be spread over many different energies.In short, even if the composite is initially in a pure nonequilibrium state, the reduced state of S int will tend to be canonical.Some of those results (which are related to and support our present work) are: (i) derivation of the canonical density operator for S int for large times [22] and for a suitable overwhelming majority of wave functions (typicality) [23]; (ii) S int will approach an equilibrium state and remain close to it for almost all times, independently both of the HB and, under suitable conditions, also of the initial state of S int [24].For equilibration of isolated systems, see [25].For further strengthening results, which include a number of extensions, small S int and equilibration for finite times, under wider conditions on energy spreads of states, as well as the interplay between quantum information and thermodynamics, see [26][27][28][29] and the references therein.
Chemical reactions constitute a very important and wide class of phenomena in which quantum mechanics and (equilibrium and non-equilibrium) statistical mechanics play essential roles.See [4,5] and the references therein for studies devoted to the equilibrium chemical constant.Much research has also been devoted to non-equilibrium chemical reactions: see, in particular, [9,11,17] for stochastic and thermodynamical approaches.The non-equilibrium rate constant is treated in [17] by means of a thorough analysis of the Kramers equation.Pulsed lasers allow one to probe, with very short temporal resolution, off-equilibrium dynamics in chemical reactions (see [30,31] and the references therein).
In biological systems, one prominent example is the thermal denaturation of DNA: a physical non-equilibrium process through which the two homologous strands, initially bound to each other, split off into two separate single ones, at a certain (melting) temperature [32,33].In a previous work [34], and motivated by non-equilibrium DNA thermal denaturation, the present authors treated the dynamics of two interacting classical three-dimensional macromolecular chains.Each chain was assumed to be Gaussian (except for possible additional weak interactions along it).These chains were immersed in a medium at thermal equilibrium and initially bound to each other by potentials, which were attractive for medium and large distances (say, Morse-like).As temperature increased from room values up to about the melting one, the two strands became dissociated, and the mean first passage time of this process was calculated.Here, we shall generalize nontrivially the methods used in [34] to a different process, namely a non-equilibrium chemical reaction of two quantum particles (atoms and small molecules) in a medium, by resorting to a somewhat simplified model.
The strategy in our present study will parallel that of [34] very closely.Then, for a better understanding of the present work, it is adequate to start by summarizing below the main successive conceptual steps followed in [34]: 1.The two macromolecular chains were described, from the outset, by a classical (non-negative) probability distribution, depending on all spatial positions of all atoms and on time and evolving through an irreversible Kramers-like master equation (stochasticity being thereby taken care of).Neither Planck's constant h nor thermal wavelengths were included.However, the formulation did include an important length scale, namely the bond length yielding the average distance between two successive atoms along each chain.2. The overall center of mass (CM) motion was factored out off-equilibrium from the remaining spatial variables (the set [z] of all positions of the atoms along the two chains and the relative position y of the two centers of mass of the chains).The interest focused then on the non-equilibrium reduced probability distribution depending on [z], y and time t.The inclusion of attractive potentials between the two chains in the corresponding Kramers-like master equation seemed to be, to the best of our knowledge, a novel feature.3. Upon integrating over all [z] (leaving y unintegrated), the classical reduced Boltzmann equilibrium distribution for the two interacting chains was used as the generating function for an infinite family of orthogonal polynomials in [z] (depending on y, parametrically).That generating function was not a Gaussian one, due to the various potentials.4. The orthogonal polynomials were used to define non-equilibrium moments (y-and t-dependent) for the non-equilibrium reduced probability distribution.5.The Kramers-like equation yielded a hierarchy for the non-equilibrium moments.6.At temperatures below, but close to, the melting one (i.e., DNA denaturation), the long-time dynamics was approximated by a Smoluchowski-like equation (containing a y-dependent effective potential) for the lowest non-equilibrium moment.7. The application of the mean first passage time formalism [16,17,35] to the Smoluchowski-like equation in step 6 enabled us to study approximately the time duration of non-equilibrium thermal denaturation of the two strands, initially bound to each other, towards configurations of two separate single strands.
We shall firstly outline the contents of the present work by emphasizing the correspondence with the above steps 1-7 in [34] and, secondly, underline the differences with [34].Section 2 deals with the non-equilibrium evolution of two quantum particles in three spatial dimensions, subject to an attractive potential between them and in contact with a large HB at thermal equilibrium.Throughout our study, we assume that the state of the system, for sufficiently long times, approaches towards its own canonical density operator at thermal equilibrium with the HB, in agreement with the results in [22][23][24][25][26][27][28][29].We shall simplify the treatment by omitting the detailed analysis of the states of the HB (just retaining the fact that it establishes the temperature), and in an effective way, we will focus exclusively on the quantum states of the two particles.The quantum Hamiltonian, nonequilibrium Wigner function and evolution equation (t-reversible and without ab initio dissipation) for the system, the separation (at both equilibrium and off-equilibrium) of the overall center of mass from the relative motion variables (position x and momentum q) of the two particles and certain spectral properties are treated in Sections 2.1-2.3 (counterparts of steps 1 and 2 in [34]).The construction of an infinite family of orthogonal polynomials in the relative momentum q, generated by the equilibrium relative Wigner function as the weight function, is carried out in Sections 2.4 and 2.5 (counterpart of step 3 in [34]).The moments (x-and t-dependent) of the non-equilibrium relative Wigner function are constructed in Section 2.6 (counterpart of step 4 in [34]).Section 3 provides the dynamical equations for the non-equilibrium moments implied by those for the non-equilibrium relative Wigner function.The infinite hierarchy is presented and analyzed in Section 3.1 (counterpart of step 5) in [34]).The restrictions to the lowest non-equilibrium moment and to an irreversible Smoluchowski-like equation for it, under various approximations, are considered and justified in Sections 3.2 and 3.3 (counterpart of step 6) in [34]).Section 4 deals with approximate kinetic equations.Section 5 compares approximately thermal and chemical equilibria.Sections 6 and 7 are devoted to the application of the mean first passage time formalism to the Smoluchowski equation for the lowest non-equilibrium moment (counterpart of step 7) in [34]).Section 8 compares briefly our approach with other well-known methods and, in so doing, treats chemical recombination.Section 9 summarizes our conclusions.
Let us quote some key differences between the classical system considered in [34] and the quantum one to be treated here.In the quantum case: (a) even if, in general, Wigner functions could be negative in some domain [36,37], we shall argue that the actual equilibrium Wigner function could still be used as a weight function to generate orthogonal polynomials, by invoking a suitable mathematical framework (Section 2.4); (b) upon constructing the orthogonal polynomials, one integrates over the relative momentum q (while in [34], one integrated over all positions [z] of the atoms along the two chains); (c) even if the system approaches thermal equilibrium for sufficiently long times [22][23][24][25][26][27][28][29], the nonequilibrium evolution equation for the Wigner function (our starting point in Section 2) does not display either irreversibility or a thermalizing evolution in general (while the Kramers-like master equation in [34] did display them from the outset).Then, approximations (for short thermal wavelengths and long times) are necessary here in order to arrive at the actual irreversible Smoluchowski-like equation for the lowest non-equilibrium moment in Section 3.3.In turn, for a simpler presentation, those approximations underlying the contents of Section 3.3 are carried out successively and separately, for a one-dimensional case, in Appendices C-F.One underlying basic question is why the starting point of our analysis is the non-equilibrium Wigner equation without explicit dissipation, instead of a master equation displaying the latter from the very outset: Appendix G will be devoted to discuss this issue, a posteriori.

General Features
We shall consider two non-relativistic spinless quantum particles of masses m j (with j = 1, 2) in three spatial dimensions, inside a large finite volume Ω (outside which, an infinitely repulsive potential is assumed).Those individual particles can be either atoms or small molecules.We emphasize that Ω is large only at the microscopic scale.In Cartesian coordinates, the position vector and the momentum operator of the j-th particle are x j = (x j,1 , x j,2 , x j,3 ) and p j = (p j,1 , p j,2 , p j,3 ), where p j,α = −ih(∂/∂x j,α ), α = 1, 2, 3 and h is Planck's constant.We shall denote x = (x 1 , x 2 ).
There are no external forces, and all interactions "seen by the system" are described, in principle, by a time-independent and velocity-independent two-body instantaneous potential V between the particles.The total (Hermitian) quantum Hamiltonian is: Let P = ∑ 2 j=1 p j be the total momentum.Then, one has: [H, P] = 0 ([A, B] = AB − BA being the commutator of the operators A and B).
We shall separate off the standard CM: its position vector being X = M −1 [∑ 2 j=1 m j x j ] with M = ∑ 2 j=1 m j .Let us introduce the relative position vector x = x 2 − x 1 and the reduced mass Then, one has: Here, H CM is the free CM Hamiltonian, and H re is the Hamiltonian of the reduced-mass particle.We now express x in terms of the standard radial (r = |x|) and angular coordinates (θ, ϕ).Then: where l is the orbital angular momentum operator.Spatial integrations over X and x will be carried out inside large finite volumes Ω.The eigenfunctions corresponding to both H CM and H re will vanish, by assumption, at the surfaces enclosing Ω.

Non-Equilibrium Statistical Formulation and Equilibrium Distributions
Let the two-particle system be immersed in an HB at thermal equilibrium, which is at absolute temperature T. The temporal evolution for t > 0 is given by the density operator ρ = ρ(t) (a statistical mixture of quantum states).It fulfills the (t-reversible) operator equation ∂ρ/∂t = (ih) −1 [H, ρ] with initial condition ρ(t = 0) = ρ in .The Hermitian and positive-definite linear operators ρ(t) and ρ in act on the Hilbert space spanned by the set of all eigenfunctions of H. Unless otherwise stated, we shall not impose that ρ(t) be normalized.
We now resort to the non-equilibrium Wigner function [1,8,[38][39][40][41][42][43], which reads as: where x = (x 1 , x 2 ), q = (q 1 , q 2 ) and dx = ∏ 2 i=1 ∏ 3 α=1 dx i,α .The time evolution of the Wigner function is furnished by: with: where dq = ∏ 2 i=1 ∏ 3 α=1 dq i,α .As Ω is large, we shall approximate spatial integrals by those for an infinite volume when such approximations are harmless, unless some specific discussion is required.Strictly speaking, as Ω is not infinite, the q are discretized momenta, and dq in Equation ( 7) should be interpreted as a six-fold series.However, as Ω is large, we shall disregard the small spacings in dq and approximate it as a six-fold integral (thus, varying continuously).A similar remark and interpretation will apply and be understood whenever integrations over momenta occur (namely, in Equations ( 10), ( 16), (18), (19), and so on).We shall accept that all integrals (or all series) over momenta converge for large values of the latter: explicit expressions for W CM,eq below, together with (11), (12) and the computations in Appendix A will support this assumption.The initial condition W in (x, q) is determined by ρ in through Equation (5).The equilibrium Wigner function W eq (x, q) is given by Equation (5) for the actual ρ eq .Let Q = q 1 + q 2 and q = (m 1 q 2 − m 2 q 1 ) /M be the CM and relative momentum vectors, respectively.The actual W eq (x, q) factorizes as: with W re,eq = W re,eq (x, q) being the relative equilibrium Wigner function.W CM,eq is proportional to exp −β(Q 2 /2M) .Correspondingly, there is also factorization off-equilibrium: W(x, q; t) = W CM (X, Q; t)W re (x, q; t), with ∂W CM /∂t = −(Q/M)∂W CM /∂X, and the (t-reversible) Wigner function W re of the relative particle fulfills: We shall introduce the fixed, physically relevant and x-independent momentum q eq = (2m/β) 1/2 .Recall that, in the classical case, the equilibrium (or Boltzmann's) canonical distribution describing thermal equilibrium of two classical particles with an HB is (proportional to a) Gaussian in q: that is, exp −β ∑ 2 j=1 ∑ 3 α=1 q 2 j,α /(2m j ) + V , with q j,α being now the classical momenta.In the high temperature or small β (quasiclassical) regime, Wigner [38] obtained successive approximations for the equilibrium quantum distribution.Based on [38], we shall give the leading terms in the small β expansion for W re,eq (x = (x 1 , x 2 , x 3 ), q = (q 1 , q 2 , q 3 )): H 1 (q α /q eq )H 1 (q α /q eq ) , (11) where: and: with c 0 = (2πh) −3 .H 2 and H 1 denote standard Hermite polynomials.We shall refer to the contributions associated with c 1 and c 2 as quantum contributions in the quasiclassical regime.
It is a reasonable assumption that the initial conditions ρ in and W in (x, q) factorize into CM and relative ones.The initial condition for W re is W re, in .
a is understood to be the range of V. 3. V(r) and all d n V(r)/dr n , for n = 1, 2, 3, . .., are continuous for all r > 0. 4. V(r) does give rise to only one bound state (bound spectrum).Thus, the relevance of the region where V < 0 is larger than that of the hard core.
As we have mentioned, each individual particle could be either an atom or a small molecule.Hence, on physical grounds, we assume that V(r) is an effective potential between them (which, in particular, includes and averages over Coulomb interactions).Let ϕ j = ϕ j (x) denote a suitably-normalized eigenfunction of H re with corresponding eigenvalue E j and j denoting a set of labels.H re has both a discrete spectrum (due to the above assumptions, just one bound state) and a denumerably-infinite number of discrete states.For the discrete spectrum E j = E d < 0 and ϕ j = ϕ d (x).The denumerably-infinite discrete spectrum has a small spacing, and it becomes a continuous one (sweeping the continuous positive real axis) as Ω −1 → 0. We shall denote it by CS, even if the small Ω −1 remains positive.The eigenfunctions corresponding to the continuous spectrum of H re are ϕ j = ϕ k (x), with j ≡ k being an almost continuous wavevector, and eigenvalues E j = E k ≥ 0. The CS eigenfunctions are normalized through: Hence, ϕ d , and all CS ϕ k span two separate Hilbert subspaces H d and H CS .Accordingly, ∑ j will include the contribution of both the single discrete eigenfunction plus that of a three-fold infinite summation over the whole CS ones.For the CS only: While at T = 0, there is no transition between the bound state and the CS ones, such a transition is indeed possible for T > 0 (due to the HB) and plays a key role as k B T approaches |E d |, as will be discussed below.Thus far, we have implemented the counterparts of steps 1 and 2 in [34].We emphasize that the contributions of the bound states have disappeared in the high temperature (quasi-classical) regime corresponding to Equations ( 11)- (14), which should not be employed in regions where bound states be relevant.However, there seems to be no compelling reason for not using Equations ( 11)-( 14) for rough or zeroth order approximations in regions where the contributions due to the bound states are negligible.

W re,eq as a Quasi-Definite Functional of Momenta
As W CM,eq is Gaussian in Q, it gives rise, as a weight function, to an infinite number of orthogonal polynomials in Q: the standard Hermite ones [44].Furthermore, notice that neither W re,eq , nor W re can be warranted to be nonnegative in general [38,39,42].A necessary and sufficient condition for the Wigner function associated with a wave function to be nonnegative is that the latter be a Gaussian distribution [36].However, the domain in which W re,eq < 0 may occur cannot be large and has to be consistent with the fact that the marginals of both W re,eq and W re are nonnegative.
Let us now invoke, in outline, a mathematical framework, based on the theory of orthogonal polynomials [45], the physical interest of which will be appreciated below.For simplicity, we shall focus on the one-dimensional case, leaving out the direct extension to three dimensions.Let us consider a kernel K = K(y) (which could be negative), a set of functions f = f (y) and the following functionals We will assume that all integrals over y are convergent.Let µ n = L K [y n ], n = 0, 1, 2, 3, . .., and consider the set of all (S + 1) × (S + 1) matrices M S (where S = 0, 1, 2, 3, . ..) with the (i, j)-th element equal to µ i+j (where i, j = 0, ..., S) and their respective determinants: det(M S ).By definition, the functional L K is quasi-definite if det(M S ) = 0 for any S = 0, 1, 2, 3, . . .[45].If L K is a quasi-definite functional, we resort to a theorem in [45], which implies the following: (i) the existence of an infinite family of orthogonal polynomials, denoted here as H K,n = H K,n (y), for n = 0, 1, 2, . . ., with weight function K (even if K < 0 in some domain of y); and (ii) the H K,n 's satisfy a basic three-term linear recurrence relation, omitted here.
Let y = q/q eq .The argument in Appendix A supports (albeit, it does not prove) that W re,eq determines a quasi-definite functional of y for any x.Based on that argument, we shall assume henceforth that, regarding their q-dependences, the Wigner function W re,eq determines a quasi-definite functional L W re,eq (for any x and t).The interest of this assumption is obvious: if it holds (as we shall suppose), it implies the existence of an infinite family of orthogonal polynomials H [n] (y) generated by W re,eq (Section 2.5).It is possible, but unnecessary, to extend the assumption to W re .

Orthogonal Polynomials H [n] Generated by W re,eq
Let us define [n] = (n 1 , n 2 , n 3 ), with the n α being non-negative integers (α = 1, 2, 3).Use will also be made of the additional shorthand notations: [0] = (0, 0, 0), [1 α ], which represents a vector with a one in the α position with the rest of the components being zero, and [2 α ] similar to [1 α ], but with a two in the α position.We introduce (unnormalized) orthogonal polynomials H [n] in y determined by W re,eq , which acts as a (in general, non-Gaussian) weight function.We choose H [0] (y) = 1.We shall also employ additional useful notations, like: . In general, the H [n] 's are constructed recurrently as follows.We impose for [n] = [n ], and any x (left unintegrated), that: where: The orthonormalized polynomials are , with a (x-dependent) normalization factor h [n] defined through: Expressing the H [n] 's as sums of products of standard Hermite polynomials is also possible, but less convenient here.The orthogonal polynomials determined by the classical relative Boltzmann distribution f re,eq in ( 12) are proportional to products of the standard Hermite polynomials H n (y α ) [44].Then, for the latter, the coefficients in the (classical) counterpart of ( 17) can be directly obtained and are x-independent.Section 2.4 and the present one have implemented the counterpart of step 3 in [34].

Moments: Off-Equilibrium and at Equilibrium
We now turn to the non-equilibrium moments.They are defined through: When W re = W re,eq , all equilibrium moments with [n] = [0] vanish, except: W eq,CS,[0] is the CS contribution.This justifies, a posteriori, the interest of having introduced the H [n] 's and the W eq,[n] 's.In the t-evolution, not far from thermal equilibrium, and recalling [22][23][24][25][26][27][28][29], one expects that W [0] would be dominant and that any W [n] , for any [n] = [0], be small (the more negligible the larger t and ∑ 3 α=1 n α are).The initial condition W in,[n] for W [n] is obtained by replacing W re by W re,in in Equation (19).W re and any W [n] have the same dimension as h−3 .This section has implemented the counterpart of step 4 in [34].

Non-Equilibrium Hierarchy
The transformation of Equations ( 9), ( 10) and ( 19) into an infinite linear hierarchy for the nonequilibrium moments W [n] can be carried out through direct computations and cancellations, increasingly cumbersome as ∑ 3 α=1 n α grows.The equations for the lowest order moments are: where: The general (t-reversible) hierarchy implied by Equations ( 9) and ( 10) for any [n] is: with: Having in mind that ,n 3 and W n 1 ,n 2 ,n 3 +1 , and analogously for the corresponding interpretation of the sum ∑ (24) (which we do not make explicit to avoid unnecessary notational complications).Here, M [n],[(n+1) γ ] , and any M [n],[n−n ] come from the first and the second terms on the right-hand side of Equation ( 9), respectively.We shall omit, for simplicity, the M [n],[n−n ] 's.It is simpler to give and discuss the counterparts of the M [n],[n−n ] 's in the one-dimensional case (see Appendix C for details).The structure of the hierarchy ( 24) is a genuine consequence of quantum mechanics.All coefficients in (24) are expressed in terms of V and of quantities computed out of the equilibrium solution W re,eq .

Order of Magnitude Estimates
From now on, we shall consider values typical of microscopic scales and phenomena.The range a of V is a few Å (for instance, 3 up to 12 Å).We shall write the reduced mass as m = n × m ne , with m ne and n being, respectively, the neutron mass and a positive (approximate) integer (not to be confused with the n's appearing as subscripts in polynomials and moments).If the individual quantum particles are atoms, then n is less than about one half of the largest atomic mass number (say, n is less than about 125).If the individual quantum particles are molecules, then n can be 125.We shall introduce the length scale δx = a/l, where l > 2, characterizing approximately the smallest scale of appreciable variations for V.For instance, δx = (1/5)a, so that 0.6 Å< δx < 2.5 Å.Moreover, the magnitude V 0 of |V|, averaged over the domain in which V < 0, lies between one and 10 electron volts.Notice that V 0 has a similar order of magnitude as the energy |E d | for the bound state (see Section 2.3).δV(< V 0 ) will denote the average variation of V within its range (where V = 0) in a scale δx.Similarly, δ n V, n = 2, 3, ..., will denote the average variation of δ n−1 V in a scale δx (δ 1 V ≡ δV).For estimates, one could regard that |δV| be about one order of magnitude smaller than V 0 and that |δ n V| be about one order of magnitude smaller than |δ n−1 V|, n = 2, 3, ....
Let λ th = h[β/2m] 1/2 = h/q eq be some suitable thermal wavelength.We shall assume that n and T are such that the following conditions are fulfilled: (a) λ th is smaller than δx, which sets a lower limit on nT; and (b) (3/2)k B T is smaller than or at most V 0 , which sets an upper limit on T. Both (a) and (b) correspond to the regime of typical chemical reactions (that is, a quantum regime, but with small λ th ).Thus, even if the relative particle is in the quantum regime (but not in the classical high-temperature one), λ th is, on average, smaller than δx and a.
By using computations in Appendix B, and as zeroth order approximations, we shall take −(1/2) for |x| > a + r 0 (dominated by the CS and, in turn, approximated through the leading term c 0 f re,eq in ( 11)) and −(λ 2 th /(δx) 2 ) for r 0 < |x| < a + r 0 (dominated by the bound state), both of them independent of the index α and of the angular coordinates.We shall disregard [2 α ],[0] and approximations for it within the hard core domain |x| < r 0 (and, more generally, contributions from the latter domain), because they will be neither required nor relevant.
All operators to be considered (namely, ) below in this section, which act on the W [n] 's in Equation (24), have dimension (time) −1 .Their orders of magnitude (denoted as (τ * ) −1 ) will now be estimated.Their inverses (namely, the τ * 's) can be regarded as effective evolution times, and henceforth, we will refer to them as such.
We now turn to the orders of magnitude involved upon applying (24).The orders of magnitude of M [n],[n−1] can be expected to be respectively similar to the various contributions of M [1 α ],[0] .On a similar basis, they can be justified by considering explicit expressions in Equations (C5)-(C9) for the counterparts of 's in the one-dimensional case (Appendix C).On the other hand, the specific quantum contributions M [n],[n−n ] , n > 1 have structures proportional to (h 2n /q 2n +1 eq ) times spatial partial derivatives of V of order 2n + 1 (for various n , 2n + 1 growing as 2 −1 ∑ 3 α=1 n α − 1), as exemplified through Equations (C2) and (C9).Such structures yield contributions smaller than that of (1/q eq )(∂V/∂x α ) by factors (λ th /δx) 2n × (δ 2n V/δV), by virtue of Assumption (a) above.Moreover, one could expect |(δ 2n V/δV)| to be smaller than unity and to decrease as n grows.We anticipate that the quantum contributions M [n],[n−n ] , n > 1, will be neglected from now on (see the analysis in Appendices C and D, in one spatial dimension).The main conclusion is that, under the above Assumption (a), the overall order of magnitude involved upon applying 24) is about the same as for M [1 α ],[0] .Namely, the sum of (i) plus (ii) plus (iii).
There are two other effective evolution times respectively associated with the operator D[[1 α ]; ] (Section 3.3 and Appendices D and E) and to the mean first passage time (Section 6).After having introduced the latter two, the estimates for all those effective evolution times will be compared and discussed a posteriori in Section 7.
The quantum contributions due to ( 13) and ( 14) in the quasiclassical regime, for |x| outside the region where the bound state is concentrated, can be regarded as negligible compared to unity in Equation (11), as V decreases very fast there.

Approximations: Small Thermal Wavelength and Long Time
By using the assumption and the order of magnitude estimates in Section 3.2, one analyzes the non-equilibrium hierarchy (24) for the W [n] 's and proceeds to the small thermal wavelength quantum regime (STWQR).In the latter, x varies by units of order δx > λ th .Any 24) becomes the approximate (t-reversible) three-term hierarchy: See Appendix D for explanations of the argument yielding (27) in the one-dimensional case.We emphasize that ( 27) is still a quantum hierarchy, since it contains all M [n],[n−1] , which are quantum-mechanical.For simplicity (and without essential loss of generality), we shall assume the initial condition . The solution of the resulting approximate non-equilibrium hierarchy ( 27) is given, through a Laplace transform from t to the variable s, in terms of products of certain s-dependent generalized operator continued fractions D[[n]; s], which satisfy the three-dimensional generalization of the one-dimensional Equation (D3).Next, we proceed to the long-time approximation (thus introducing irreversibility) for t longer than a certain largest effective evolution time τ * max , to be estimated in Section 7, and also by following the one-dimensional case of Appendices D and E (in particular, Equation (E1)).One gets (for fixed and small s = > 0) the following irreversible approximate quantum (Smoluchowski-like) equation for the lowest non-equilibrium moment W [0] (x; t): with D[[1 α ]; ] being a quantum operator and assuming approximately the initial condition ] ], also displays quantum effects (compare with Equation (E1)).Our Equation ( 28) embodies stochasticity and displays a structure typical of diffusion-convection-reaction equations: it is linear; convection is determined by the V-dependent term; and external sources are absent in it.Section 3.2 and the present section have implemented the counterpart of step ( 6) in [34].

Towards Kinetic Equations
In what follows, we will consider and interpret Equation ( 28) by replacing approximately the operator D[[1 α ]; ] by a positive and α-independent constant D, by assuming that [2 α ],0 = [2],0 (r) does not depend on α nor on the angular variables and by supposing that (the angular variables are absent).Such Equation (28), with the latter approximation, will be equivalent to Equation (46) in Section 6.Let r 1 < Ω 1/3 be some radial coordinate just slightly larger than r 0 + a, so that V is negligible at r 1 and, of course, for r > r 1 .Let S(r 1 ) be a spherical surface with center at r = 0 and radius r 1 .Let n = (n 1 , n 2 , n 3 ) be the outer normal unit vector at a generic point on S(r 1 ), and let Ω(r 1 ) be the sphere enclosed by S(r 1 ).First, we shall integrate Equation (28) over Ω(r 1 ) (at fixed t) and apply Green's theorem so as to transform the volume integral on the resulting right-hand side into a surface integral over S(r 1 ) (with differential surface element dS(r 1 )).The t-dependence will not be written explicitly.We have: where the subscript r 1 denotes the derivative at r = r 1 , and we have neglected the contribution of V.In the rightmost side of Equation ( 30), we have replaced with small ∆r (such that r 1 − 2 −1 ∆r < r 0 + a).Second, we shall integrate Equation ( 28) over the volume (denoted as Ω − Ω(r 1 )) corresponding to r 1 ≤ r ≤ [(3Ω)/(4π)] 1/3 .The latter volume is enclosed by S(r 1 ) and by the spherical surface (denoted as S(Ω)) with center at r = 0 and radius [(3Ω)/(4π)] 1/3 (also at fixed t).We shall also apply Green's theorem so as to transform the latter volume integral into two surface integrals over S(r 1 ) and S(Ω).Furthermore, let n(Ω) = (n 1 (Ω), n 2 (Ω), n 3 (Ω)) be the outer normal unit vector at a generic point of S(Ω).Then, by proceeding as in the derivation of Equation ( 30), we get: For large Ω, by neglecting I(Ω) and approximating like in Equation (30), Equation (31) becomes: Notice that, consistently: which supports a posteriori disregarding I(Ω).For adequately large t, equilibrium sets in and the right-hand sides of the approximate Equations ( 33) and ( 30) vanish.We shall analyze the resulting equilibrium equation in Section 5. Next, assume for the time being that, for adequately large t and as rough approximations: (i) for x inside Ω(r 1 ), spatial variations of W [0] are not significant, and the latter is dominated by the contribution of the bound state; so that Ω(r 1 ) are not important, and the latter is dominated by the contribution of the CS, so that Ω−Ω(r 1 ) d 3 xW [0] (Ω − Ω(r 1 ))W [0] (r 1 + 2 −1 ∆r).The volume of the sphere Ω(r 1 ) has been denoted as Ω(r 1 ), and so on, for Ω − Ω(r 1 ).Then, Equations ( 33) and ( 30) become two approximate linearly-coupled kinetic equations for and: For illustrative purposes, let us compare the above expressions with standard chemical kinetics [11].Let us suppose that particles 1 and 2 are similar (or identical, ignoring possible effects due to quantum-mechanical symmetrization).For a chemical reaction involving two chemical species A 1 , A 1 and A d :

and its inverse reaction
the phenomenological kinetic equations read as [11]: where C 1 and C d are the t-dependent concentrations of A 1 and A d , with k + and k − being suitable forward and backward rate constants.We stress that ( 37)-( 38) are non-linearly-coupled equations for C 1 and C d .In principle, the W [0] 's in different regions arising in the approximate linear Equations ( 35) and (36) (with the CM evolution factored out) do not appear to coincide with the concentrations in (37) and (38), so that, strictly speaking, there should be no essential contradiction or problem of principle in the simultaneous approximate validity of both pairs of equations.In spite of that, one could still ask whether there could be some analogy pointing out towards the approximate compatibility between the linear (35) and (36) and the non-linear (37) and (38).We shall outline an argument indicating that such a consistency appears to hold, at least in a semi-quantitative sense and adequately close to equilibrium.Recall that the W [0] 's have the same dimension as h−3 .We shall introduce the t-dependent quantities a d and a CS through Here, k is an arbitrary constant, to be fixed from the outset, while k is another constant to be determined below.The underlying physical purpose is to interpret a d and a CS as proportional to concentrations inside Ω(r 1 ) and Ω − Ω(r 1 ), respectively.For large t, we argue that a CS equals the constant equilibrium value a CS,eq plus a smaller (decreasing) t-dependent contribution, so that: a CS (da CS /dt) a CS,eq (da CS /dt).Then, Equations ( 35) and (36) become: Notice that k has factored out.Next, we compare the phenomenological Equations ( 37)-(38) to Equations ( 39) and ( 40), interpreting a d as C d and a CS as C 1 .The two sets are compatible with each other if k a CS,eq fulfills: Ω(r 1 ) = (Ω − Ω(r 1 ))2k a CS,eq .The latter equation and the restriction of CS at equilibrium determine k and a CS,eq (once k has been fixed).Then, for large t, the non-linear Equations ( 37) and ( 38) and the linear Equations ( 39) and ( 40) would be compatible.

Thermal and Chemical Equilibria
Let us now return to Equations ( 15) and (20).Upon integration over the large finite volume Ω: Then, the right-hand side of Equation ( 41) contains the contribution of the discrete spectrum plus the one due to the whole, almost continuous spectrum, which is finite (as long as Ω is finite, but diverging like Ω if Ω 1/3 → ∞).If, for a finite (adequately large) Ω, T decreases within the range specified in Section 3.2, then the discrete contribution to W re,eq dominates: the two particles remain bound to each other.As T increases (for fixed large Ω), the almost continuous spectrum contribution to W re,eq dominates.Those two extreme situations and the intermediate ones could correspond to thermal equilibrium, but there may be, in general, an imbalance between the discrete and the CS contributions to (41), unless some further condition be fulfilled.To grasp the latter, notice that there exist key ranges of intermediate T and Ω 1/3 , such that the orders of magnitude of the discrete and the CS contributions are not significantly different from each other.From Equation ( 20), the discrete contribution to (41) is: The order of magnitude of the CS contribution in Equation ( 41) can be estimated by disregarding the bound state term and approximating W re,eq W eq,CS,[0] (recall (20)) through ( 11) and ( 12) (disregarding c 1 and c 2 ).Therefore, with Ω − Ω(r 1 ) Ω, a zeroth-order approximation to the CS contribution in (41) is: There would be a comparative similarity of the orders of magnitude of ( 42) and ( 43) provided that exp (−βE d ) be about Ω/(8π 3/2 λ 3 th ), which holds for suitable Ω, n and T. We shall now discuss the size of the microscopically large Ω.One could argue physically that our two-particle system is contained in some macroscopic gaseous system at T, containing a very large number N of similarly interacting pairs.The gaseous system occupies a macroscopically large volume Ω Ω.Then, one could choose Ω Ω /N.For typical diatomic gases at adequately low pressures and T's not smaller than room temperature, one has that Ω is about 10 4 -10 5 Å 3 per pair.Such an order of magnitude for Ω could allow for the compatibility between ( 42) and (43).
Based on the above comparison, we now turn to the possibility of another kind of equilibrium for large t, namely chemical equilibrium.In fact, at equilibrium, the vanishing of the right-hand sides of the approximate Equations ( 33) and ( 30) can be recast as: At equilibrium, we approximate W [0] (r 1 − 2 −1 ∆r) by q −3 eq exp (−βE d ) Ω(r 1 ) −1 and W [0] (r 1 + 2 −1 ∆r) by q −3 eq (8π 3/2 λ 3 th ) −1 , consistently with ( 42) and ( 43), and we replace the [2 α ],[0] 's by their zeroth order approximations in Section 3.2.Then, (44) becomes: where Ω(r 1 ) is of order a 3 (see Section 3.2).Equation ( 44) is, within our rough approximations (say, about one order of magnitude more or less), compatible with the approximate similarity of the orders of magnitude of Equations ( 42) and ( 43) when the latter occurs (which, in turn, holds if Ω(r 1 )/Ω is about 2(λ th ) 2 /(δx) 2 )).Hence, there is a simultaneous coexistence of thermal and chemical equilibria in the actual simple framework.Moreover, it appears adequate to regard Equation ( 44) (in which no chemical potentials have been introduced) as a counterpart of the law of mass action for equilibrium chemical reactions.Notice also that ) by a factor of about Ω(r 1 )/Ω: as a side remark, this justifies the neglect of the CS contribution compared to the bound state one, when estimating [2 α ],0 for r 0 < |x| < r 0 + a in Appendix B.

Transitions among a Bound and Continuous Spectrum: Mean First Passage Time
In a similar fashion as in Section 4, let us replace approximately the quantum operator D[[1 α ]; ] by a positive and both α-and x-independent constant D, and suppose that [2 α ],[0] depends solely on the radial coordinate r and that it is negative.Let us restrict to W [0] = W [0] (r, t), independent by assumption on the angular coordinates.Implicitly, we are accepting that W [0] (r, t) is an adequate (not to say, dominant) approximation to the full W [0] .Then, Equation ( 28) becomes: It will be useful to replace W [0] (r, t) by a new distribution f = f (r, t) defined through W [0] (r, t) = r −2 f (r, t).Then, Equation (46) reduces to: We are interested in finding approximately the time required for the relative particle to proceed from the domain in which the bound state is concentrated to a continuous spectrum configuration (that is, large r).For that purpose, we shall apply the mean first passage time formalism (MFPT) [16,17,35].We shall extend the treatment in [34].An MFPT τ(r) (which can be regarded as another effective evolution time) is the solution of the so-called adjoint equation associated with Equation ( 47), namely of: ∂V ∂r 2q eq [2],[0] mr ∂τ(r) ∂r Given that r 0 is an adequately small radius, such that V > 0 for 0 ≤ r < r 0 (i.e., the hard core), and that for the appreciably larger interval r 0 < r < a + r 0 , V < 0 and |V| takes on its largest values, let r 2 be considerably larger than r 0 + a (and so, than the r 1 in Sections 4 and 5), so that V is negligible at r 2 .Since r 2 is smaller than [(3Ω)/(4π)] 1/3 , we shall look for the solution τ(r) of Equation (48), such that ∂τ(r)/∂r = 0 at r = r 0 ("reflecting boundary condition") and τ(r) = 0 at r = r 2 ("absorbing boundary condition").One finds readily for r 0 < r < r 2 that: In Section 7 and Appendix F, we shall be able to approximate the operator D[[1 α ]; ] by a positive and α-independent function D(r), instead of by a constant D directly.One would like to be able to control the further approximation of the function D(r) by a suitable constant D upon evaluating the MFPT, so as to perform estimates with Equation (49).For that purpose, the following remark regarding Equation ( 49) will be useful.Let us come back to Equation ( 28), replace in it the operator D[[1 α ]; ] by a function D(r) and rewrite the resulting Equation (28) as Here, L S is the operator defined by the right-hand side of Equation ( 28) when the operator D[ [1 α ]; ] is replaced by the function D(r).Let L + S be the adjoint of L S .An MFPT τ L (r) for L + S can be introduced, by extending [16,17,35], as the solution of the differential equation L + S τ L (r) = −1 in r 0 < r < r 2 , with similar boundary conditions ∂τ L (r)/∂r = 0 at r = r 0 ("reflecting boundary condition") and τ L (r) = 0 at r = r 2 ("absorbing boundary condition").A calculation similar to that yielding Equation ( 49) leads to: with the same J(r ) as in Equation ( 49).The interest of Equation ( 51) is two-fold: either D(r ) inside the integral could now be approximated by some constant D (some average of D(r ) in r < r < r 2 ) and, so, allow for the use of τ(r) and Equation ( 49) (as done below in this section) or variations of D(r ) inside the integral could be considered in certain cases (see Section 7).We rewrite Equation (50) as: Recall that it was supposed that V(r) has no discontinuities and that it becomes negligible as r > r 0 + a.As stated in Section 3.2, V 0 lies between 1 and 10 eV, and T satisfies Assumptions (a) and (b) there.Let r be larger than r 0 + a (with r < r 2 ).Then, it is permissible to replace [2 α ],[0] by its zeroth order approximation, given in Section 3.2.Notice that [2 α ],[0] for r < r 0 is not needed.After that approximation, one performs the three integrations over r in (52).We take V(r ) 0 for r > r 0 + a, V(r 0 + a) −V 0 and V(r ) 0 for r > r > r 0 + a.Hence, Equation (52) simplifies to: Then, approximating D(r ) by some (average) constant D and using τ(r) together with Equation ( 49), we get: The last contribution in the rightmost side of (54) has been neglected by invoking that exp (βV 0 ) is adequately larger than unity (unless k B T is about V 0 ) and the approximate relationships in Section 5.
As a rough approximation, we replace exp (V(r 0 + a) − V(r )) for any r 0 < r < r 0 + a inside the integral in (the β-independent) J 1 by some (r and T independent) κ.The magnitudes of J 1 , the constant κ and related quantities will be discussed briefly in Section 8.Then, Equation (54) reduces to: Consistently with the assumptions in Section 2.3 and to fix the ideas, the shape of V resembles that of the Morse potential.Moreover, τ(r) −1 could be interpreted as some sort of rate constant for the transition between the bound state and the continuous spectrum states.At this point, we recall the exponential-like features of the rate constant in the classical Arrhenius formula and in the transition state theory [35,46].We see that our quantum-mechanical approach has led approximately to the exponential Arrhenius factor, namely exp (βV 0 ).Notice that Equation (51), by extending the approximations for τ(r), also leads to the exponential Arrhenius factor for τ L (r).However, our approach fails so far to predict τ(r) (or τ L (r)) strictly, as long as it involves an unknown prefactor, namely either the function D(r) or the average constant D (approximating roughly the unknown operator D[[1 α ]; ]).Through further approximations, estimates of D and of τ(r) will be provided in Section 7.

Estimating Effective Evolution Times
We disregard the fact that D[[1 α ]; ] is an operator and approximate it by the function D(r).The latter can be estimated by extending directly to three dimensions the one-dimensional analysis in Appendix F. This yields the following rough estimate for D(r) (in seconds): In turn, by recalling Section 3.2, we shall interpret −M as the product of the constant h/(λ th mδx) (associated with M [0],[(1 α ] , that is with (q eq /m)(∂/∂x α )) times the overall order of magnitude estimate for M [1 α ],[0] (namely, the sum of the estimates (i) plus (ii) in Section 3.2, for r 0 < |x| < r 0 + a and |x| > r 0 + a, respectively), which depends on r = |x|.Then, Equation (57) does depend on r.For r 0 < r < r 0 + a, the dominant contributions give (δV in eV, λ th , n and δx being discussed in Section 3.2): which is T-independent.
For r > r 0 + a, the overall order of magnitude estimate for M [1 α ],[0] is about (1/2)h/(λ th mδx).Then: which depends on T like T −1/2 .Therefore, for adequately small λ th (fulfilling Assumptions (a) and (b) in Section 3.2, the estimate in (59) would decrease and become closer to that in (58).Then, provided that the orders of magnitude of the latter two are not too different from each other, the constant D considered in Section 6 and, in particular, in Equations ( 49), ( 54) and ( 56) is some average of (58) and (59).Recall that the average time for ironing out local thermal inhomogeneities over one nanometer in a gas under normal conditions is about 10 −11 s [3].
Then, returning to Section 3.3 and as in Appendix F, it is reasonable to interpret τ * max D/2 and to expect that the long-time approximation be adequate for t > D/2.
Turning to the estimate for τ(r) in the rightmost side of (56), we get (r 0 + a < r < 10a): in which D is an average of ( 58) and (59) when they do not differ much.Notice that if (58) is appreciably smaller than (59), then the use of Equation ( 51) would be more adequate than that of Equation (49).For r 0 < r < r 0 + a, the integral in (51) could then be restricted to r 0 +a r dr (q eq /m) 2 D(r )r 2 J(r ), where D(r ) should be replaced by (58).For r > r 0 + a, D(r ) should be replaced by (59) in the integral in (51).Both resulting integrals could be estimated as in Section 6, so as to yield further values of the MFPT.This section and the preceding one have implemented the counterparts of step 7 in [34].

Comparison with Other Approaches: Recombination and Possible Extensions
Notice that, through Condition 1 in Section 2.3, we have restricted to a potential V(r), which, in particular, is attractive (< 0) in the interval r 0 < r < (3Ω/4π) 1/3 (say, with a global structure resembling that of a Morse-like one).From the conditions in Section 2.3 and through the developments and approximations up to and including Section 3, we have arrived at the quantum Equation (28), which, in particular, provided the basis for the MFPT analysis in Section 6.Such a study can be compared to that in Subsection VII.C.4 in [35] for a radial classical Smoluchowski equation with a damping or friction parameter (γ f r , given from the outset) and a similar spherically-symmetric potential (Morse-like), in connection with diffusive problems and, specifically, with a recombination chemical reaction (see, in particular, Equation (7.30) and Figure 24 in [35]).The corresponding MFPT analysis yields Equations (7.31)-(7.33) in [35].One key difference between the classical Equations (7.30)-(7.31) in [35] and our Equations ( 46)-( 48) is that the latter contain the quantum contribution [2],[0] .Another difference is that our approach allows one to estimate, after a number of approximations, the quantum counterpart of γ f r , namely the values of D(r) (or D): see Appendix F and Section 7. Notice that the values of [2],[0] in r 0 < r < a + r 0 (dominated by the bound state), estimated in Section 3.2, are responsible for the Arrhenius factor exp (βV 0 ) in (56).Our (quantum) solution ( 49) and (50) is also applicable to similar diffusive problems and is the counterpart of (with the same boundary conditions as) the classical solution in Equation (7.33) in [35].
It is worth emphasizing that the MFPT solution in Equation (7.32) in [35] corresponds to the opposite boundary conditions (absorption at r 0 and reflection far outside the range of the potential), and it has been applied to recombination processes in chemical reactions: see [35] and the references therein for generalizations.This suggests the interest of extending succinctly our study in Section 6 above so as to find the MFPT solution τ rec (r) for Equations (47) and (48) for boundary conditions ∂τ rec (r)/∂r = 0 at r = r 2 (reflection) and τ rec (r) = 0 at r = r 0 (absorption boundary) and to compare it with [35].The operator D[[1 α ]; ] is approximated by the same positive and α-independent function D(r) as in Section 6.A new computation, similar to that yielding Equations ( 50) and (51), gives readily for r 0 < r < r 2 : which are the actual quantum counterparts of Equations (7.32), (7.34) and (7.35) in [35].Here, Equations ( 61) and (62) would give approximately the (recombination) time required for the relative particle to proceed from the domain corresponding to the continuous spectrum configuration (r 0 + a < r) to the one in which the bound state is concentrated (that is, r 0 < r < r 0 + a).Equations ( 61) and ( 62) can be estimated, by approximating D(r ) by a constant (D) and extending readily the approximations leading from ( 49) and ( 50) to (54).One finds: where: with exp (−βV 0 ) being an Arrhenius-like factor with a negative exponent, decreasing as T decreases.Notice that Dτ rec,1 (r) depends on temperature only through β exp (−βV 0 ), while Dτ rec,2 (r) is temperature-independent.The possible magnitudes of the exponentials exp for various r and r , which appear in J 1 in Equation (55) (and so, the related constant κ in Equation ( 56)), in τ rec,1 (r) and in τ rec,2 (r) depend rather strongly on whether V varies rapidly or smoothly and on the sign of V(r ) − V(r ).Let V(r ) − V(r ) > 0 (as, otherwise, those exponentials could be neglected).
For suitably smooth V (say, rather small δV), those exponentials would not be much greater than unity (for instance, κ could not exceed 10 2 or 10 3 ), but they could be certainly large for rapidly-varying V.
In any case, regardless of the magnitudes of those exponentials, we stress that J 1 (and hence, κ), together with the integrals in τ rec,1 (r) and τ rec,2 (r) are temperature independent.A numerical study of those integrals lies outside our scope here.The quantum Equations ( 63)-( 65) can be compared with the classical Equations ( 74) and (75) in [35], which appear to yield a different temperature dependence.
Recall that our methods rely on Assumption (a) in Section 3.2, imposing a lower limit on temperature.The activated barrier crossing problem (say, for transitions from a bound state to continuous state configurations) has attracted enormous research attention in connection with chemical reactions [17,35,46]: in one simple version, it corresponds to a potential V(r), which presents a "hard core" for 0 ≤ r < r 0 ; it is attractive in the interval r 0 < r < r 3 , repulsive in the interval r 3 < r < r 4 and vanishes fast as r → +∞.Our Equation (28) and the arguments up to and including Section 3 would also hold if Condition 1 in Section 2.3 is replaced by: 1*.V(r) is repulsive (> 0) for 0 ≤ r < r 0 ("hard core", with adequately small r 0 ), attractive (< 0) in the interval r 0 < r < r 3 , repulsive (> 0) in the interval r 3 < r < (3Ω/4π) 1/3 and vanishes fast as r → (3Ω/4π) 1/3 .The estimate of [2],[0] in r 0 < r < a + r 0 (dominated by the bound state) in Section 3.2 should now be revised, but it would not appear to spoil the validity of the corresponding Equation (28).The MFPT analysis in Section 6 up to and including Equations ( 46)-(51) would also apply (as would that leading to Equations ( 61) and ( 62)).However, the analysis in Sections 4 and 5 and the approximate estimates in Section 6 beyond (51) should be revised.Therefore, in order to apply our methods to the activated barrier crossing problem, those parts of our analysis should have to be extended, which lie outside our scope here.

Conclusions and Discussions
In this work, we have put forward a nontrivial generalization of the methods previously used in the study of non-equilibrium DNA thermal denaturation, presented in [34], to provide a simple model for a binary chemical reaction starting from a quantum mechanical framework.The close correspondence between the methods and successive steps in [34] and the present work, as well as the main differences between them, have been highlighted.We have treated a quantum two-particle system in three spatial dimensions, subject to an attractive potential V and to a "heat bath" (HB) at thermal equilibrium at absolute temperature T > 0. We have focused on a quantum regime typical of chemical reactions, such that: (a) the thermal wavelength λ th is shorter than the range of the attractive potential; (b) the energy for the discrete bound state |E d | (about the magnitude of the attractive potential within the negative region) is, in magnitude, ≥ (3/2)k B T. Our starting point has been the non-equilibrium time-reversible Wigner equation without explicit dissipation (as justified in Appendix G).We have separated the overall center of mass (CM) of the relative motion.We have focused on the dynamics of the non-equilibrium Wigner function W re for the relative particle motion.We have considered an infinite family of orthogonal polynomials H [n] (depending on the relative momentum), generated by the relative equilibrium Wigner function W re,eq (non-Gaussian and non-positive, in general).These orthogonal polynomials have generated (by integrating over the relative momentum) non-equilibrium moments W [n] of W re .In turn, the general non-equilibrium Wigner equation for W re has implied a linear nonequilibrium hierarchy for the W [n] 's (time-reversible, as well, and depending on the relative position).Under the above Assumptions (a) and (b), and for sufficiently long times t ≥ τ * max , the linear nonequilibrium hierarchy has been approximated by the irreversible Smoluchowski-like Equation ( 28 This has led to effective linear kinetic equations allowing, in turn, for approximate comparisons (with standard phenomenological non-linear kinetic equations for large t) and between thermal versus chemical equilibria.Then, we proceed to the simpler Equation (46), which depends only on t and on the radial separation r between both particles and still preserves a quantum signature via the term [2 α ], [0] .The order of magnitude of D has been estimated.Equation (46) has allowed us to apply the mean first passage time (MFPT) formalism [16,17,35] which, in turn, has given rise to an explicit representation for the time τ = τ(r) required for the two particles to proceed from the bound state to the continuous spectrum configurations, as a function of r.An explicit formula for τ(r), in terms of D and displaying the genuine Arrhenius exponential factor exp [V 0 /(k B T)], has been obtained and discussed briefly.A consistent extension, in which the operator D[[1 α ]; ] is replaced by a suitable function of the radial distance, has also been outlined.The MFPT has also been applied briefly to chemical recombination.Our MFPT studies are compared also briefly to other MFPT studies by other authors [35].For the sake of a simpler presentation without loss of generality, the various approximations leading to the Smoluchowski-like equation for W [0] have been carried out in detail for a one-dimensional case in Appendices C-F.We point out that no use of non-equilibrium Wigner functions, moments and hierarchies was made in the works [22][23][24][25][26][27][28][29].Among several open problems (and hence, possible future uses of the methods in the present paper), we quote: (a2) extensions to include more than one bound state between the two particles, thereby allowing for transitions among those states (eventually, for instance, transitions resembling vibrational ones approximately) and (b2) generalizations to more general two-particle chemical reactions (like A + B → C + D). d 3 yW re,eq y n α d 3 yW re,eq , (B1) By recalling Equation (20), one finds formally: d 3 x d 3 yW re,eq ≥ 0. In turn, d 3 x d 3 yW re,eq appears in the quantum theory of the second virial coefficient (see, for instance, Section 14.3 of [3]): it is the sum of the discrete spectrum contribution (finite) plus that of the continuous one divergent as d 3 x).A sum of two terms with similar properties occurs for (3).
In principle, it is unclear whether the left-hand side of (B2) is ≥ 0 for any x.However, by using Equation (1) and d 3 xϕ * j (x)H re ϕ j (x) = E j d 3 xϕ * j (x)Hϕ j (x) for both the discrete and the continuous spectra, and by integrating (B2) with d 3 x, the V − E j term in (B2) equals −(h 2 /2m)(∂ϕ j /∂x α )(∂ϕ * j /∂x α ).One finds the right-hand side of (B3), as it stands.Recall that the contribution of the almost continuous spectrum is finite, as Ω is finite.Then, one can assert that, for finite Ω, the right-hand side of (B3) is well defined and ≥ 0.

Appendix C. One-Dimensional Non-Equilibrium Hierarchy (1)
We shall now study the counterpart of the model in Sections 2-6 in one spatial dimension (x), by omitting unnecessary details.There is just one particle, subject to a real potential V = V(x), with V(x) = V(−x) and V(x) < 0 in −a < x < a, and to a HB at T, in the large finite interval −L/2 ≤ x ≤ L/2 (a < L).The Hamiltonian is H = −(h 2 /2m)(∂ 2 /∂x 2 ) + V. Now, there are a discrete spectrum (d) and an almost continuous one (j = CS), with energies E j and eigenfunctions ϕ j = ϕ j (x) and spanning the Hilbert subspaces H d and H CS .Correspondingly, ∑ j denotes sums over all eigenfunctions in j = d and j = CS.For t > 0, the exact (t-reversible) dissipationless quantum master equation for the general off-equilibrium Wigner function W [38,39] is: with initial condition W in .The equilibrium density operator ρ eq = exp[−βH] determines the equilibrium Wigner function W eq (x, q): As in Section 2.2 and without further discussion, as L is large and unless otherwise stated, we shall approximate spatial integrals by those in −∞ < L < +∞ and series over momenta by integrations over them in −∞ < q < +∞.
We shall introduce the (unnormalized) polynomials in y (= q/q eq ) H Q,n = H Q,n (y) (n = 0, 1, 2, 3, . . .), orthogonalized in y (for fixed x) by using the equilibrium distribution W eq as the weight function.One has: H Q,n (y) = y n + ∑ j n,n−j y n−j , with coefficients n,n−j .The actual H Q,n (y)'s lead to the non-equilibrium moments (n = 0, 1, 2, . ..): The initial condition W in,n for W n is obtained by replacing W by W in in Equation (C4).For W = W eq (x, q), Equation (C4) yields W eq,n = 0 if n > 0, while W eq,0 = q −1 eq ∑ j exp −βE j ϕ j (x)ϕ * j (x).The transformation of the one-dimensional Equations (C1) and (C2) into a linear hierarchy for the nonequilibrium moments W n can be carried out through computations and cancellations, increasingly cumbersome as n grows.We have obtained the first five equations in that quantum hierarchy: The very fact that the full quantum equation for n = 4 does contain a term of order h2 in W 1 implies that the quantum hierarchy is not a three-term hierarchy.The general (t-reversible) hierarchy implied by Equations (C1) and (C2) for any n is: The M n,n 's for n = 0, 1, 2, 3, 4 are identified upon comparing Equation (C10) and Equations (C5)-(C9).M n,n =0 = 0 for any n, except for n = 1 and M n,n−n = 0 if n is even.In the exact nonequilibrium quantum hierarchy (C10), the contributions from W n+1 always have the same structures (−(q eq /m)∂W n+1 /∂x, with n-independent coefficients).On the other hand, the contributions from W n (0 < n ≤ n − 1) do carry n-dependent coefficients, which increase with n.In particular (and leaving aside other contributions) the equation for ∂W 5 /∂t can be shown to contain, in its right-hand side, (h 2 /q 3 eq )(∂ 3 V/∂x 3 )W 2 as the highest spatial derivative of V, while that for ∂W 6 /∂t contains (h 2 /q 3 eq )(∂ 3 V/∂x 3 )W 3 and (h 4 /q 5 eq )(∂ 5 V/∂x 5 )W 1 , and so on.In turn, the nonvanishing n,n−2 for low order n = 2, 3, 4, 5 are: y n = dyW eq (x, q)y n dyW eq (x, q) .(C16) A key feature of the non-equilibrium hierarchy (C10) is that all coefficients in it are expressed in terms of V and of quantities computed out of the equilibrium solution W eq .The complicated structure of the hierarchy (C10) is a genuine consequence of quantum mechanics.By invoking [22][23][24][25][26][27][28][29], (C10) should simplify necessarily for very long t, with W n → 0 for n > 0, while W 0 = 0 with W 0 → W eq,0 and M 1,0 W eq,0 = 0.One expects that, for adequate t and by solving the hierarchy (C10) for n > 0, all W n for n > 0 can be expressed, through suitable linear operators, in terms of W 0 and of suitable initial conditions W in,n , for n > 0. In particular, one would get W 1 as a linear functional of W 0 (and of all W in,n , n > 0): ], in short.Then, the hierarchy (C10) boils down to ∂W 0 /∂t = −(q eq /m) ∂W 1 [W 0 ; [W in,n , n > 0]]/∂x (linear in W 0 ), plus the initial condition W in,0 .These remarks could be directly generalized to the hierarchy (24).We shall not undertake the latter very complicated problem, in general.Rather, we shall concentrate on an approximate version of it in the following Appendices, under conditions relevant for chemical reactions.
where we have supposed that (−M 1,2 M 2,1 ) > 0 (as it seemingly happens, so as to be consistent with classical Brownian motion and, hopefully, in the case of interest here), taken the "+ root upon solving for D[1; , x] (so D[1; , x] > 0), assuming that < 2(−M 1,2 M 2,1 ) 1/2 , and finally, used (ii) above.The rightmost side of (F2) provides the rough estimate for the function D[1; , x] and, hence, for the order of magnitude of the contribution of the operator D[1; ].At a later stage, one could approximate the function D[1; , x] by some average constant, D.Then, the long-time approximation could be expected to hold for t > D/2.This sort of estimate is used directly in Section 7, for the three-dimensional case.
We now assume that − 2,0 ≥ 0 for any x and that the operator D[1; ] is approximated by the function D[1; , x] > 0. In order to connect with Brownian motion, we shall introduce the following non-equilibrium quantum (Helmholtz) free energy in terms of W 0 = W 0 (x; t): by integrating in the large finite interval −L/2 ≤ x ≤ L/2.N 0 is a suitable constant, while c = c(x) is x-dependent, but t-independent.Both N 0 and c will be determined below suitably.This definition of A in Equation (F6) generalizes the non-equilibrium classical (Helmholtz) free energy treated in [47] for classical Brownian motion.We consider ∂A/∂t, use Equation (E1) and perform a partial integration over x.Then: , x](q eq /m) 2 (− 2,0 ) W 0 At this stage, we shall choose c, such that: Upon comparing Equation (F8) with M 1,0 W eq,0 = 0, we find: c = −k B T ln W eq,0 , (F9) having chosen, for simplicity, a vanishing integration constant.As W 0 ≥ 0 and D[1; , x] > 0 and with the choice for c in Equations (F9), one gets the expectedly consistent behavior of A: ∂A/∂t ≤ 0, which expresses irreversibility.The constant N 0 will be determined with the additional requirement that, at thermal equilibrium, A, as defined by Equation (F6), coincides with the equilibrium quantum (Helmholtz) free energy, which reads: A eq = −k B T ln 1 q eq ∑ j exp −βE j dxϕ j ϕ j * .(F10) Therefore, Equation (F6), when W 0 = W eq,0 , should become (F10).For that to be true, it is necessary that: ln N 0 = ln dxW eq,0 dxW eq,0 . (F11)

Appendix G. Alternative Starting Point: A Master Equation with Ab Initio Dissipation
We shall outline here the reasons for having based our analysis on the non-equilibrium Wigner Equations ( 6) and (7) (or ( 9) and ( 10)) without explicit dissipation.One first reason is that the research in [22][23][24][25][26][27][28][29] yielded an approach to the equilibrium canonical quantum distribution, independently both on the HB and, essentially, also on the initial state.Therefore, our starting point has been independent of the interactions and mechanisms determined by external sources (except on temperature, imposed by the HB).Two further supporting requirements are: (a1) the existence of bound states should be allowed and their role displayed explicitly in formulations and in approximations; (b1) if one tries a master equation displaying dissipation from the very outset, the corresponding equilibrium quantum distribution should be independent of the dissipation mechanism built in that master equation.We shall explore succinctly the compatibility of several Wigner function master equations containing explicit ab initio dissipation with (a1) and (b1).For simplicity, we shall treat directly the evolution of the relative particle in the one-dimensional case.The evolution of the CM can be factored out as in Section 2.2 and will be disregarded (in particular, it is trivially independent of bound states).Accordingly, we shall modify Equation (C1) for the relative particle by adding the dissipation term DW so as to deal with the master equation: ∂W(x, q; t) ∂t = − q m ∂W(x, q; t) ∂x M Q W is the same as in (C2).DW is the dissipation term, D being a suitable operator considered below.In all possibilities, γ D will denote a positive friction constant.
Possibility (2) [41,48]: to first non-leading order in the high temperature or small β (quasiclassical) regime.The equilibrium distribution obtained by Wigner [38] in such a regime for the one-dimensional case is indeed the equilibrium distribution for (C1), and it also satisfies: DW eq (x, q) = 0 to the same order [41] (allowed by the fact that the actual D is V-dependent).Therefore, (G3) fulfills (b1) in the first non-leading order in the quasiclassical regime.The interest of (b1) was stressed, in particular, in [41].However, (G3) does not satisfy (a1): even if V is attractive, the contribution of the bound states is not explicitly displayed in the small β regime.Possibility (3): where γ D,1 is a constant.Equation (G4) with γ D,1 = 0 is the Caldeira-Leggett model for quantum Brownian motion [13,49].Notice that Equation (G4) with a suitable γ D,1 = 0 is regarded in [13] as an extension of [49], enabling the latter to belong to the dissipative dynamics theories in [14,15].Furthermore, Equation (G4) with either γ D,1 = 0 or γ D,1 = 0 fails to fulfill (a1) and (b1), for similar reasons as above.There could still exist perhaps an open possibility in the dissipative dynamics framework in [14,15], namely when the corresponding equilibrium distribution (C3) fulfills DW eq = 0 (thereby satisfying (b1).Such theories have been treated succinctly in [13], without providing any explicitly manageable expression for D, thereby precluding extending our analysis in this paper.Moreover, it is fully open whether the procedure yielding such an operator D [13] would hold for attractive potentials, so as to be consistent with (a1).
To conclude: we have been unable to find Wigner function master equations with explicit ab initio dissipation that could fulfil (a1) and (b1) in general (although we have not excluded completely that they could exist, at least formally, in the framework of [14,15]).For the above reasons, our study has relied on the non-equilibrium Wigner Equations ( 6) and (7) (or ( 9) and ( 10)) without explicit dissipation.
)with ∂ n,n−2 /∂x = d n,n−2 /dx, ∂V/∂x = dV/dx.The n,n−2 's in the first four Equations (C5)-(C8) contain quantum effects.On the other hand, Equation (C9) (for n = 4) acquires an additional term of quantum origin multiplying W 1 , and so, it differs from (C5)-(C8).The reason for that difference is that the quantum corrections in Equations (C1) and (C2) manifest themselves only at order h2 and, then, in turn, in the equations in the hierarchy at orders n ≥ 4.