Entanglement Effect and Angular Momentum Conservation in a Non-separable Tunneling Treatment

The famous transition state (TS) theory high pressure rate constant expression reads¹ (we use atomic units, $\hbar=1$)

$$k=\frac{k_{B}T}{2\pi}\frac{Z^{\#}}{Z_{R}}e^{-V^{\#}/k_{B}T},$$

(1.1)

where $T$ and $k_{B}$ are temperature and the Boltzmann constant, $Z_{R}$ is the partition function of the reactant(s), and $V^{\#}$ is the barrier height. The partition function for the TS dividing surface $Z^{\#}$ can be factorized into rotational and vibrational contributions in the framework of the rigid-rotor harmonic-oscillator (RRHO) approximation at the saddle point,

	$\displaystyle Z^{\#}$	$\displaystyle=$	$\displaystyle Z^{\#}_{rot}Z^{\#}_{vib},$		(1.2)
	$\displaystyle Z^{\#}_{rot}$	$\displaystyle=$	$\displaystyle 2\sqrt{2\pi}\sqrt{I_{x}I_{y}I_{z}}\beta^{-3/2},$		(1.3)
	$\displaystyle Z^{\#}_{vib}$	$\displaystyle=$	$\displaystyle\prod_{i}\frac{1}{2\sinh(\omega_{i}\beta/2)},$		(1.4)

where $\beta=1/k_{B}T$, $I_{i}$ are the principal inertia moments at the saddle point, and $\omega_{i}$ are the frequencies for motions transverse to the reaction coordinate.

Ad hoc approaches for accounting for the effects of tunneling often assume a separability between the effects of tunneling and the contribution of the transverse modes to the transition state partition function. Commonly, this assumption is implemented in terms of a separability of the reaction coordinate and the transverse modes. Then an additional tunneling factor enters into Eq. 1.2,

	$\displaystyle Z^{\#}$	$\displaystyle=$	$\displaystyle Z^{\#}_{tun}Z^{\#}_{rot}Z^{\#}_{vib},$		(1.5)
	$\displaystyle Z^{\#}_{tun}$	$\displaystyle=$	$\displaystyle\beta\int dEP(E)e^{-\beta E},$		(1.6)

where $P(E)$ is the one-dimensional barrier transmission coefficient at energy $E$. In the spirit of the RRHO approximation, a parabolic barrier approximation can be used, for which the tunneling transmission coefficent reads,

$$P(E)=\frac{1}{1+e^{-2\pi E/\omega_{b}}},$$

(1.7)

where $\omega_{b}$ is the imaginary frequency at the saddle point and the saddle point is taken as the zero of energy. Substituting Eq. 1.7 into Eq. 1.6 one obtains the classic Wigner tunneling expression,²

$$Z^{\#}_{tun}=\frac{\omega_{b}\beta/2}{\sin(\omega_{b}\beta/2)}.$$

(1.8)

While Eq. 1.8 takes into account local tunneling effects, it fails at low temperatures because the integral in Eq. 1.6 diverges at low energies for the parabolic barrier. This divergence occurs at temperatures lower than the crossover temperature, $T_{c}$,

$$k_{B}T_{c}=\frac{\omega_{b}}{2\pi}.$$

(1.9)

To avoid this divergence one should use a more realistic model for the transmission coefficient. In the one-dimensional WKB approximation, the transmission coefficient reads,

$$P(E)=e^{-\tilde{S}(-E)},\ E<0,$$

(1.10)

where $\tilde{S}(E)$ is the abbreviated action over the periodic trajectory along the reaction coordinate $r$ with the energy $E$ in the inverted barrier potential $-V(r)$,

$$\tilde{S}(E)=2\int_{r_{-}}^{r_{+}}dr\sqrt{2m_{r}(E+V(r))},$$

(1.11)

where $r_{\pm}$ are the left and right turning points and $m_{r}$ is the mass associated with the reaction coordinate. Evaluating the integral in Eq. 1.6 in the stationary phase approximation (SPA) one obtains,

$$Z^{\#}_{tun}=\sqrt{2\pi}\beta\left(\frac{d^{2}\tilde{S}(E)}{dE^{2}}\right)^{-1/2}e^{-\tilde{S}(E)+\beta E},\ \frac{d\tilde{S}}{dE}=\beta,$$

(1.12)

which is valid in the ”deep tunneling” regime where the Wigner expression fails.

From Eqs. 1.11 and 1.12 one can see that, in contrast with the localized TS dividing surface of relevance at higher temperatures, in the deep tunneling regime the contribution from tunneling through the TS is related to the dynamics over a relatively wide range of coordinates. As a result, in the deep tunneling regime, the separability assumption between the reaction coordinate and the transverse modes becomes questionable and it is desirable to construct a theory that avoids this assumption. Over the years, tremendous effort has been devoted to constructing such a theory. For conciseness, we review only a few results of particular relevance to our further discussion.

Long ago Langer³ suggested that the tunneling rate constant can be calculated via analytical continuation of the system partition function $Z=Z_{R}+\frac{i}{2}Z_{I}$, considered as a path integral in the configurational space of the closed paths. The imaginary part in this integral comes from the vicinity of the trajectory that corresponds to the saddle point in that space and the path integral is considered in the quadratic approximation. The TS partition function $Z^{\#}$ in Eq. 1.1 then reads,

$$Z^{\#}=2\pi Z_{I},$$

(1.13)

Affleck⁴ considered the transition from low to high temperatures for one-dimensional systems and has confirmed that Eq. 1.13 coincides with Eq. 1.12 at temperatures below the crossover.

Many prior studies focused on the evaluation of the microcanonical rate constant rather than the canonical one. Miller⁵, who was one of the pioneers in the study of muti-dimensional tunneling effects in chemical reactions, applied his correlation function formalism⁶ and Gutzwiller’s expression for the Green’s function⁷ to obtain the expression for the “cumulative reaction probability”, which we prefer to call the microcanonical number of states, $N(E)$, in the spirit of RRKM theory.⁸ Notably, it was later realized that the expression obtained for $N(E)$ did not agree with the exact expression in the separable multi-dimensional case and a correction was suggested.⁹ Curiosly enough, Miller in his original paper⁵ also suggested an expression for the thermal rate constant that is correct in the separable case but his derivation was inaccurate.¹⁰ Furthermore, as we see, this expression appears to be missing key terms of relevance to non-separable systems. Alternative, non-separable tunneling path approximations have been developed for the ground vibrational state in the framework of the vibrational adiabaticity assumption, starting with the seminal paper of Marcus and Coltrin.¹¹ Truhlar and co-wokers, who have originally developed many of these approaches,^{12, 13} extensively used them in thermal rate constant calculations.^{14, 15, 16, 17, 18, 19}

Benderskii et al.^{20, 21} further discussed the derivation of the rate constant within Langer’s approach. By neglecting the interaction between longitudinal and transverse fluctuations, they obtained an explicit expression for the rate constant that coincided with Miller’s expression for the thermal rate constant.⁵ It includes the product $1/\prod_{i}2\sinh(u_{i}/2)$ as a prefactor, which looks formally similar to the vibrational partition function, cf. Eq. 1.4. Furthermore, the stability parameters $u_{i}$, which were first introduced by Gutzwiller,^{5, 7} approach their classical values $u_{i}=\beta\omega_{i}$ when the temperature approaches $T_{c}$, making the analogy with the vibrational partition function complete. Jackels et al.²² contributed to the practical aspects of the tunneling path calculations through a potential expansion in the vicinity of the minimum energy path (MEP). Richardson^{23, 24} has used an improved Green’s function approximation within Miller’s correlation function expressions.²⁵ Many papers have successfully used the so called ring polymer molecular dynamics approach for the tunneling rate constant calculation.^{26, 27} This approach implements Langer’s theory directly by approximating the partition function path integral with a large number of beads.

To summarize, two generic classes of results are available in the literature. These two classes differ in the presumed starting point for the derivation. The papers based on Miller’s correlation function formalism express the rate constant in terms of a Boltzmann weighted integral over energy. The integrand is expressed via Green’s functions, for which different approximations are used. The advantage of this approach is that it is better justified and allows for a smooth transition from low to high temperatures. Its drawback is that current implementations do not reproduce the exact separable result, cf. Eq.1.5, because of the difficultities in evaluating the Green’s function.

On the other hand, the papers based on Langer’s approach typically either use different kinds of decoupling approximations to obtain the final rate constant or calculate the fluctuation part of the path integral directly. We did not find in the literature a closed rate expression for the non-separable reactive system. So the main goal of this paper is to provide such an expression that can be substituted for Eq. 1.5 in the deep tunneling regime. Also, all analytical papers we have seen, with one exception, do not take into account the conservation of angular momentum. The one exception is a paper by Miller⁵, which, as we will discuss later, is flawed. Therefore, the second goal of this paper is to fill this gap by providing an expression for the analog of the rotational partition function in the deep tunneling regime, and which goes smoothly to the classical rotational partition function at higher temperatures. We will also work through various ramifications of our results, focusing on perturbative expansions.

The layout of the paper is as follows. In Sec. 2 we will discuss the general features of the problem and consider the one-dimensional rate expression. In Sec. 3 we will consider the coupled reaction coordinate + vibrational modes problem when the energy is the only integral of motion. The conservation of angular momentum is considered in Sec. 4. The paper concludes with a lengthy general discussion, Sec. 5 and some brief summarizing remarks, Sec. 6.

We will use Langer’s approach³ to calculate the TS partition function at low temperatures. Let us consider the Hamiltonian system with $N$ degrees of freedom with unit masses in the potential $\bar{V}(q)$,

$$\hat{H}=-\frac{1}{2}\sum_{i}\frac{\partial^{2}}{\partial q_{i}^{2}}+\bar{V}(q).$$

(2.1)

In what follows we will often use vector notation. For example, $\dot{q}^{2}$ will mean $\sum_{i}\dot{q}_{i}^{2}$.

The density matrix $\hat{\rho}=e^{-\beta\hat{H}}$ can be written in the form of the path integral,²⁸

	$\displaystyle\rho(q_{1},q_{2};\beta)$	$\displaystyle=$	$\displaystyle\overbrace{\underbrace{\int\cdots\int}}\mathcal{D}q\exp(-S[q]),$		(2.2)
	$\displaystyle S[q]$	$\displaystyle=$	$\displaystyle\int_{0}^{T}dt(\dot{q}^{2}/2-V(q)),$		(2.3)
	$\displaystyle T$	$\displaystyle=$	$\displaystyle\beta,\ q(0)=q_{1},\ q(T)=q_{2}.$

The function $S[q]$ formally coincides with the action of the classical particle in the inverted potential $V(q)=-\bar{V}(q)$ and $T$ denotes the trajectory period in this section (not to be confused with the temperature).

The partition function of the system can be written as $Z=Tr[\hat{\rho}]=\int d^{N}\!q\rho(q,q;T)$. The tunneling trajectory, also known as the instanton trajectory, $q_{0}(t)$ corresponds to the saddle point of the action $S$ in the configurational space of the closed paths. The instanton trajectory is the key object in Langer’s theory. Considering small variations around the instanton path, one comes to the conclusion that it satisfies the classical equations of motion,

$$\ddot{q}_{0}=-\frac{\partial V}{\partial q}_{0},\ S(q_{1},q_{2},t)=S[q_{0}].$$

(2.4)

Then, taking into account that

$$\dot{q}_{2}=\frac{\partial S(q_{1},q_{2},t)}{\partial q_{2}},\ \ \dot{q}_{1}=-\frac{\partial S(q_{1},q_{2},t)}{\partial q_{1}},$$

(2.5)

and that the first derivative of $S$ over $q_{1}=q_{0}(0)=q_{2}=q_{0}(T)$ should be zero, one finds that the initial and final momenta coincide, $\dot{q}_{0}(0)=\dot{q}_{0}(T)$, that is the instanton trajectory is periodic. In the case of a separable reaction coordinate the instanton trajectory coincides with that described in the introduction for the WKB approximation.

The tunneling action in the quadratic approximation about the instanton trajectory reads,

$$S^{\prime}[\delta q]=\frac{1}{2}\int_{0}^{T}dt(\delta\dot{q}^{2}-\sum_{i,j}\frac{\partial^{2}V(q_{0}(t))}{\partial q_{i}\partial q_{j}}\delta q_{i}\delta q_{j}).$$

(2.6)

In the quadratic approximation the partition function $Z_{I}$ can be written as

	$\displaystyle iZ_{I}=\exp(-S(T))\int d^{N}\!q^{\prime}\overbrace{\underbrace{\int\cdots\int}}\mathcal{D}q\exp(-S^{\prime}[q]),$		(2.7)
	$\displaystyle q(0)=q(T)=q^{\prime},$

where $S(T)=S[q_{0}]$. The expression under the integral can be viewed as the infinite dimensional gaussian integral,

$$Z_{I}\sim\int dq_{t}dq_{2}...exp(-\lambda_{1}q_{1}^{2}-\lambda_{2}q_{2}^{2}-...).$$

(2.8)

One $\lambda_{i}$ in this equation is negative. It corresponds to the unstable mode at the saddle point. After analytical continuation it should give a purely imaginary value. We will see that it happens automatically when we calculate the path integral in Eq. 2.7.

There is another $\lambda_{i}$ that is identically zero. This feature is related to the fact that the Lagrangian equations, Eq. 2.4, are invariant under a time shift and therefore the instanton trajectory is not a point in the configurational space of the closed paths but rather a circle, which is obtained by an arbitrary time shift transformation to the individual instanton trajectory. The linearized trajectory, which corresponds to the infinitesimal time-shift and which has zero $\lambda_{i}$, is the derivative of the original instanton $q_{0}(t)$ over the shift $\tau$, $\dot{q}_{0}(t)$. The way to handle this divergence is to leave the corresponding gaussian integral,

$$\int\exp(-S^{\prime}[\dot{q}_{0}]x^{2})dx,$$

(2.9)

in Eq. 2.7 intact and substitute it with the time shift transformation period, which is $T=\beta$.

Actually, any continuously parametrized transformation that preserves the Lagrangian in Eq. 2.6 will generate an additional dimension for the instanton manifold and a corrseponding $\lambda_{i}$ with zero value. So, for the angular momentum conservation case one will have four zero $\lambda_{i}$’s. To handle the three additional divergent gaussian integrals in Eq. 2.7 one again leaves these integrals, which correspond to three infinitesimal orthogonal rotations, intact and substitutes them with the total solid volume of the three dimensional rotational group, which equals to $8\pi^{2}$.

To proceed we use Feynman’s result about path integrals with quadratic action²⁸. He showed that the path integral $K(q_{2},q_{1},t)$ with the action $S^{\prime}[q]$ given by Eq. 2.6 for trajectories, which start at time $t=0$ with $q=q_{1}$ and finish at time $t$ with $q=q_{2}$ can be written as

	$\displaystyle K(q_{2},q_{1},t)=U(t)e^{-S^{\prime}(q_{1},q_{2},t)},$		(2.10)
	$\displaystyle S^{\prime}(q_{1},q_{2},t)=S^{\prime}[\bar{q}],$

where $\bar{q}$ is the trajectory that satisfies the linearized equations of motion,

$$\ddot{\bar{q}}_{i}=-\sum_{j}\frac{\partial^{2}V(q_{0}(t))}{\partial q_{i}\partial q_{j}}\bar{q}_{j}.$$

(2.11)

with $\bar{q}(0)=q_{1}$ and $\bar{q}(t)=q_{2}$. Further we will consider only the trajectories which satisfy Eq. 2.11 and will omit the bar over $q$ notation. The function $U(t)$ does not depend on $q_{1}$ and $q_{2}$. In Appendix A we show that $U(t)$ is expressed in terms of the final coordinates of $N$ trajectories $q^{(1)}(t),q^{(2)}(t),\ ...\ q^{(N)}(t)$, which start with the initial velocities $\dot{q}^{(1)}(0),\dot{q}^{(2)}(0),\ ...\ \dot{q}^{(N)}(0)$, and zero initial coordinates as (see also Ref.²⁹),

$$U(t)=(2\pi)^{-N/2}\sqrt{\frac{|\dot{q}^{(1)}(0),\dot{q}^{(2)}(0),...,\dot{q}^{(N)}(0)|}{|q^{(1)}(t),q^{(2)}(t),...,q^{(N)}(t)|}}.$$

(2.12)

Substituting Eqs. 2.10 and 2.12 into Eq. 2.7 one obtains the final expression for $Z_{I}$,

$$iZ_{I}=e^{-S(T)}U(T)\int^{\prime}d^{N}\!q^{\prime}e^{-S^{\prime}(q^{\prime})},$$

(2.13)

where $S^{\prime}(q^{\prime})=S^{\prime}(q^{\prime},q^{\prime},T)$. The prime notation on the integral means that the divergent terms must be substituted with the appropriate values as described above. It is easy to see that the following expression for the tunneling action $S^{\prime}(q^{\prime})$ holds,

$$S^{\prime}(q^{\prime})=\frac{1}{2}q^{\prime}(\dot{q}(T)-\dot{q}(0)).$$

(2.14)

Thus, all necessary quantities, cf. Eqs. 2.12 and 2.14, are expressed via the initial and final values of the coordinates and momenta. In what follows we will refer to them at $t=0$ unless it is explicitly stated otherwise.

Using the obtained relations let us consider the one-dimensional case. We will choose the starting point of the instanton trajectory away from the turning point, $\dot{q}_{0}\neq 0$. To calculate $S^{\prime}(q)$, Eq. 2.14, we need to find the trajectory for which $q(T)=q$. But we know such a trajectory. It is $\dot{q}_{0}(t)$. We note that up to the $\dot{q}_{0}$ factor the corresponding integral coincides with Eq 2.9.

To caculate $U(T)$ we need to find another trajectory that satisfies Eq. 2.11. There are different ways to do so. We choose a way that proves useful later for the multi-dimensional case. Let us consider $q_{0}$ as a function of energy. Then the trajectory

$$q_{E}(t)=\frac{\partial q_{0}(t,E)}{\partial E}$$

(2.15)

satisfies the linearized equations of motions, Eq. 2.11, because $q_{0}(t,E)$ satisfies Eq. 2.4 and it is linearly independent of $\dot{q}_{0}$. The latter can be seen from the fact that the linearized energy, which is an integral of motion for the linearized equations of motion, Eq. 2.11, is zero for the $\dot{q}_{0}$ trajectory and unity for the $q_{E}$ one. Without loss of generality one can assume that $q_{E}(0)=0$. Then, from the energy expression $E=\dot{q}^{2}/2+V(q)$ one obtains

$$\dot{q}^{E}=1/\dot{q}_{0}.$$

(2.16)

To obtain the value of $q_{E}(t)$ at $t=T$ one notes that $q_{0}(t,E+\delta E)$ is also a periodic trajectory with a slightly different period $T(E+\delta E)=T+\frac{dT}{dE}\delta E$. Then it is easy to see that

$$q_{E}(T)=-\frac{dT}{dE}\dot{q}_{0}.$$

(2.17)

Note the sign in this equation. This means that the tunneling trajectory $q_{0}(t)$ has passed through the focal point and the determinant ratio in Eq. 2.12 has a negative sign. Substituting Eqs. 2.16 and 2.17 into Eq. 2.12 one obtains

$$U(T)=\sqrt{-1}\left(\dot{q}_{0}\sqrt{2\pi\frac{dT}{dE}}\right)^{-1}.$$

(2.18)

Substituting Eq. 2.18 and $\beta$ for $\int\exp[-S^{\prime}(\dot{q}_{0})x^{2}]dx$ into Eq. 2.7 and then into Eq. 1.13 one obtains,

$$Z_{tun}=\beta\sqrt{2\pi}\left(\frac{dT}{dE}\right)^{-1/2}e^{-S(T)}.$$

(2.19)

To match with Eq. 1.12, one notes that the full action $S(q^{\prime},q^{\prime\prime},t)$, Eq. 2.3, considered as a function of time, is related with the abbreviated action $\tilde{S}(q^{\prime},q^{\prime\prime},E)$,

$$\tilde{S}(q^{\prime},q^{\prime\prime},E)=\int\dot{q}dq,$$

(2.20)

considered as a function of energy, by Legendre’s transform,

$$\tilde{S}=S+Et,\ \frac{\partial S}{\partial t}=-E,\ \frac{\partial\tilde{S}}{\partial E}=t.$$

(2.21)

Then, taking into account Eq. 2.5 and that the tunneling trajectory is periodic, one obtains that $T=\frac{d\tilde{S}}{dE}$ and Eq. 1.12 is recovered.

In this section we will denote $\dot{q}$ as $p$ and the state vector $(q,p)$ as $v$. Let us consider the linear $2N$-dimensional transformation $\hat{M}$, known as the monodromy matrix,⁷ from $v(0)$ to $v(T)$

$$v(T)=\hat{M}v(0).$$

(3.1)

For Hamiltonian systems the monodromy transformation is symplectic, as it preserves the natural symplectic form, $\hat{\Omega}=q\wedge p$,³⁰

$$\hat{\Omega}(v_{1},v_{2})=q_{1}p_{2}-q_{2}p_{1}.$$

(3.2)

We will call the operation $\hat{\Omega}(v_{1},v_{2})$ the symplectic product of vectors $v_{1}$ and $v_{2}$ and denote it as $v_{1}\odot v_{2}$.

The eigenvalues of a symplectic transformation come in pairs, $\lambda_{+}\lambda_{-}=1$. ³⁰ To each $\lambda_{i}\neq 1$ there corresponds a unique eigenvector $v_{i},\ \hat{M}v_{i}=\lambda_{i}v_{i}$. For the unit eigenvalue the situation is more complicated and the corresponding eigenspace is transformed through itself.

Now we turn to the analysis of the monodromy transformation for our problem. The eigenvectors with $\lambda\neq 1$ correspond to the vibrational transverse modes. With energy being the only integral of motion, there are $2(N-1)$ of them. We will denote them as $v_{i,\pm}=(q_{i,\pm},p_{i,\pm}),\ i=1,..,N-1$. We will also denote $\lambda_{i,+}$ as $\lambda_{i}$. To proceed we will assume that the tunneling trajectory has a turning point where the system comes to complete rest, $\dot{q}_{0}=0$, and, for a moment, that the $t=0$ moment corresponds to it. Then, due to time-reversal symmetry, the vectors $v_{i,+}$ and $v_{i,-}$ correspond to the trajectories that are time-reversed relative to each other, so that $q_{i,+}=q_{i,-}=q_{i}$ and $p_{i,+}=-p_{i,-}=p_{i}$. Considering the symplectic product $v_{i,s}\odot v_{j,r}$ before and after the monodromy transformation for each pair of $v_{i,s}$ and $v_{j,r}$ one finds that

$$p_{i}q_{j}=\delta_{i,j}$$

(3.3)

with a proper normalization. If $p_{i}q_{i}<0$ one can always exchange $v_{i,+}$ and $v_{i,-}$.

The two lineary independent vectors which are left are associated with the reaction coordinate. They correspond to the $\dot{q}_{0}(t)$ and $q_{E}(t)=\frac{\partial q_{0}(t,E)}{\partial E}$ linearized trajectories and constitute the eigenspace for $\lambda=1$. Using the same arguments that have been used to derive Eq. 2.17, one obtains

$$\hat{M}v_{E}=v_{E}-\frac{dT}{dE}v_{0}.$$

(3.4)

Using Eq. 3.2 and the linearized energy form,

$$H(q,p)=g_{0}q,\ \ g_{0}=\frac{\partial V}{\partial q}_{0},$$

(3.5)

($g_{0}$ is the potential gradient at the turning point) the following additional orthogonality relations can be obtained,

	$\displaystyle v_{0}\odot v_{i}$	$\displaystyle=$	$\displaystyle g_{0}q_{i}=0,$		(3.6)
	$\displaystyle v_{E}\odot v_{i}$	$\displaystyle=$	$\displaystyle q_{E}p_{i}=0,$		(3.7)
	$\displaystyle q_{E}g_{0}$	$\displaystyle=$	$\displaystyle 1.$		(3.8)

We will proceed with calculation of the prefactor $U(T)$. One can use $p_{i}$ as starting momenta in Eq. 2.12 for the unstable trajectories. Then the corresponding coordinates, obtained after $t=T$, are given by

$$q^{f}_{i}=\frac{1}{2}(\lambda_{i}-\lambda_{i}^{-1})q_{i}.$$

(3.9)

Unfortunately, one cannot use the turning point $q_{\pm}$ as the starting point for one of the reaction coordinate trajectories $\dot{q}_{0}(t)$ and $q_{E}(t)$, because $q_{\pm}$ is a focal point for the periodic trajectory that has been started from it, cf. Eq. 2.18. Therefore we will make a small time-shift $\tau$ from the turning point and use perturbation theory to calculate the determinants ratio in Eq. 2.12. The calculation is rather cumbersome and therefore has been moved to Appendix B. The result is shown below,

	$\displaystyle U(T)$	$\displaystyle=$	$\displaystyle(2\pi)^{-N/2}\tau^{-1}\prod_{i}\left[\frac{1}{2}(\lambda_{i}-\lambda_{i}^{-1})\right]^{-1/2}\sqrt{-\frac{|g_{0},p_{i}|}{|g_{0},q_{i}|}}$
		$\displaystyle\times$	$\displaystyle\left(\frac{dT}{dE}g_{0}^{2}+\frac{2}{g_{0}^{2}}\sum_{i}(p_{i}g_{0})^{2}\frac{\lambda_{i}+\lambda_{i}^{-1}-2}{\lambda_{i}-\lambda_{i}^{-1}}\right)^{-1/2}.$

Here and in what follows we use the shorthand notation for the group of vectors under the determinant, for example, $|g_{0},p_{i}|=|g_{0},p_{1},p_{2},...|$. It is worthwhile noting again the minus sign under the square root in this equation.

Now we turn to calculating the action integral in Eq. 2.13. For the action one can use Eq. 2.14. Due to Eq. 3.3 the action is diagonalized if one uses $q_{i}$ as the new basis. Expressing the linearized trajectory $v(t)$, for which $q(0)=q(T)=q_{i}$, in terms of the linear combination of the eigenvectors, $v(t)=C_{i,+}v_{i,+}(t)+C_{i,-}v_{i,-}(t)$, at $t=0$ and $t=T$, one obtains $C_{i,+}+C_{i,-}=\lambda_{i}C_{i,+}+\lambda_{i}^{-1}C_{i,-}=1$, which gives $C_{i,+}=1/(1+\lambda_{i}),\;C_{i,-}=\lambda_{i}/(1+\lambda_{i})$. Substituting $C_{i,\pm}$ into Eq. 2.14, $\dot{q}(T)-\dot{q}(0)=(C_{i,+}\lambda_{i}-C_{i,-}/\lambda_{i}-C_{i,+}+C_{i,-})p_{i}$, the expression for the action coefficient of the diagonal term reads $(\lambda_{i}-1)/(\lambda_{i}+1)$, where we again have used Eq. 3.3.

The zero-mode term in the integral in Eq. 2.13 is reduced to the standard form, Eq. 2.9, if one uses $\dot{q}_{0}=-g_{0}\tau$ as the basis vector in it, which corresponds to the $\dot{q}_{0}(t)$ trajectory. Substituting the zero-mode standard integral with $\beta$, and using Eq. 3. one obtains the final expression for $Z^{\#}$,

	$\displaystyle Z^{\#}$	$\displaystyle=$	$\displaystyle\beta\sqrt{2\pi}F_{tun}\prod_{i}[2\sinh(u_{i}/2)]^{-1}e^{-S(\beta)},$		(3.11)
	$\displaystyle F_{tun}$	$\displaystyle=$	$\displaystyle\left(\frac{d^{2}\tilde{S}}{dE^{2}}+\frac{2}{g_{0}^{4}}\sum_{i}(p_{i}g_{0})^{2}\tanh(u_{i}/2)\right)^{-1/2},$		(3.12)
	$\displaystyle S(\beta)$	$\displaystyle=$	$\displaystyle\tilde{S}(E)-E\beta,\ \frac{d\tilde{S}}{dE}=\beta,$		(3.13)

where we have used the conventional stability parameters $u_{i}$, $\lambda_{i}=e^{u_{i}}$. One can see that the tunnneling prefactor term given by Eq. 3.12 differs from the corresponding one-dimensional one $d^{2}\tilde{S}/dE^{2}$, cf Eq. 2.19, by the presence of an additional term.

In the case when the potential is invariant under spatial rotations there are three additional integrals of motion for Eq. 2.4, which are the angular momentum $J=q\times\dot{q}$ components. We will use the vector notation $a\times q$, which could mean a vector product $a\times q_{i}$ for each $i$-th atom or the sum $\sum_{i=1}^{N_{a}}a_{i}\times q_{i}$ depending on the nature of $a$.

We will assume that $t=0$ corresponds to the turning point. Then, in the linearized form, the angular momentum integral of motion is given by

$\displaystyle J(q,p)=q_{0}\times p.$

(4.1)

Together with Eq. 3.5 they can be used to classify linearly independent trajectories.

There are $2(N-4)$ unstable modes ($N=3N_{a}-3$ is the number of vibrations and rotations), which correspond to the transverse vibrations. As in the previous section, $v_{i,\pm}=(q_{i},\pm p_{i})$ with the eigenvalues ($\hat{M}v_{i,\pm}=\lambda_{i,\pm}v_{i,\pm}$) $\lambda_{i,\pm}=\lambda_{i}^{\pm 1}$. The vectors $p_{i}$ and $q_{i}$ satisfy the orthogonality relations, Eq. 3.3.

The two vectors associated with the reaction coordinate are $v_{0}=(0,-g_{0})$ and $v_{E}=(q_{E},0)$. Together with $v_{i}$ they satisfy Eqs. 3.6-3.8. While $v_{0}$ is the eigenvector of the monodromy transformation $\hat{M}v_{0}=v_{0}$ with unit eigenvalue, $v_{E}$ satisfies Eq. 3.4. It is worth noting that, due to rotational invariance of the potential, the turning point actually is not a point anymore but rather a three-dimensional sphere in the $q$-subspace. As a result, the vector $q_{E}$ is defined up to an arbitrary linearized rotation $n\times q_{0}$. We will choose it so that

$$q_{E}\times q_{0}=0.$$

(4.2)

There are six additional vectors associated with the angular momentum. Three of them are obtained by infinitesimal rotations of the instanton trajectory,

$$v_{\xi}=(q_{\xi},0),\ q_{\xi}=n_{\xi}\times q_{0},\ \xi=x,y,z,$$

(4.3)

whery $n_{\xi}$ is the unit vector in the direction of $\xi$. Because the potential is invariant relative to rotations, vectors $q_{\xi}$ satisfy the relations,

$$q_{\xi}g_{0}=0.$$

(4.4)

Vectors $v_{\xi}$ are the eigenvectors of the monodromy transformation $\hat{M}v_{\xi}=v_{\xi}$ with unit eigenvalue. Applying the monodromy transformation to $v_{\xi}\odot v_{i,\pm}$ one obtains,

$$v_{\xi}\odot v_{i,\pm}=q_{\xi}p_{i}=0,\ q_{0}\times p_{i}=0.$$

(4.5)

Three other vectors $v^{J}_{\xi},\ \xi=x,y,z$ correspond to non-zero values of the angular momentum, Eq. 4.1. We will choose them in the following way. In the $q$-space the vectors $q_{E}$, $q_{\xi}$, and $q_{i}$ constitute the complete basis. The vectors $g_{0}$ and $p_{i}$ are complemetary to $q_{E}$, $q_{\xi}$, and $q_{i}$, cf. Eqs. 3.3, 3.6-3.8, 4.4, and 4.5. But they are incomplete. We will add to them the vectors $p_{\xi},\ \xi=x,y,z$ so that they will form a complete basis that is complementary to the basis $q_{E}$, $q_{\xi}$, and $q_{i}$,

	$\displaystyle p_{\xi}q_{\eta}$	$\displaystyle=$	$\displaystyle\delta_{\xi,\eta},$		(4.6)
	$\displaystyle p_{\xi}q_{i}$	$\displaystyle=$	$\displaystyle 0,$		(4.7)
	$\displaystyle p_{\xi}q_{E}$	$\displaystyle=$	$\displaystyle 0.$		(4.8)

Then the vectors $v^{J}_{\xi}$ can be defined as,

$$v_{\xi}^{J}=(0,p_{\xi}),\ \xi=x,y,z.$$

(4.9)

Substituting $v_{\xi}$ into Eq. 4.1 one finds that $Jv_{\xi}=n_{\xi}$. The symplectic product $v_{\xi}\odot v^{J}_{\eta}$ is given by,

$$v_{\xi}\odot v^{J}_{\eta}=\delta_{\xi,\eta}.$$

(4.10)

and the following products are zeros,

$$v_{E}\odot v^{J}_{\xi}=0,\ v_{i}\odot v^{J}_{\xi}=0,$$

(4.11)

cf. Eqs. 4.6-4.8.

The result of the application of the monodromy matrix to $v^{J}_{\xi}$ can be written in its most general form as,

$$\hat{M}v_{\xi}^{J}=v_{\xi}^{J}+C_{0}^{\xi}v_{0}+\sum_{i,s}C_{i,s}^{\xi}v_{i,s}+\sum_{\eta}C_{\xi,\eta}v_{\eta},$$

(4.12)

where we have used the integrals of motion, Eqs. 3.5 and 4.1, to exclude $v_{E}$ and $v_{\eta}^{J},\ \eta\neq\xi$ terms. Considering the symplectic products $v_{E}\odot v^{J}_{\xi}$ before and after the monodromy transformation one concludes that $C_{0}^{\xi}=0$. Similarly, considering the symplectic products $v_{i,\pm}\odot v^{J}_{\xi}$, one finds that $C_{i,\pm}^{\xi}=0$. As a result, Eq. 4.12 is considerably simplified,

$$\hat{M}v_{\xi}^{J}=v_{\xi}^{J}+\sum_{\eta}C_{\xi,\eta}v_{\eta}.$$

(4.13)

Considering the symplectic products $v_{\xi}^{J}\odot{v}^{J}_{\eta}$ before and after the monodromy transformation and using Eq. 4.10 one finds that the matrix $C_{\xi,\eta}$ is symmetric,

$$C_{\xi,\eta}=C_{\eta,\xi}.$$

(4.14)

From Eq. 4.13 it follows that if the matrix $C_{\xi,\eta}$ is not degenerate, the six-dimensional unit eigenspace related to the angular momentum is factorized into three two-dimensional subspaces $\tilde{v}_{\xi}^{J},\ \sum_{\eta}C_{\xi,\eta}v_{\eta},\ \xi=x,y,z$.

Now we turn to the calculation of the prefactor $U(T)$, Eq. 2.12. For the same reason as in the previous section we apply to the system a small time-shift and consider the corresponding trajectories using perturbation theory. The corresponding calculation is pretty much the same as in Appendix B for the system without angular momentum conservation. The differences are traced in Appendix C. Thus, the expression for $U(T)$ is given by, cf. Eq. C.,

	$\displaystyle U(T)$	$\displaystyle=$	$\displaystyle(2\pi)^{-N/2}\tau^{-1}\sqrt{-\frac{|g_{0},p_{i},p_{\xi}|}{|g_{0},q_{i},\sum_{\eta}C_{\xi,\eta}q_{\eta}|}}\prod_{i}\left[\frac{1}{2}(\lambda_{i}-\lambda_{i}^{-1})\right]^{-1/2}$
		$\displaystyle\times$	$\displaystyle\left(\frac{dT}{dE}g_{0}^{2}+\frac{2}{g_{0}^{2}}\sum_{i}(p_{i}g_{0})^{2}\frac{\lambda_{i}+\lambda_{i}^{-1}-2}{\lambda_{i}-\lambda_{i}^{-1}}\right)^{-1/2}$

where $C_{\xi,\eta}$ is the matrix from Eq. 4.13. One can simplify this expression by multiplying the numerator and denominator under the square root with $|g_{0}q_{i}q_{\xi}|$ and then consider the matrix product $(g_{0}q_{i}q_{\xi})^{T}\times(g_{0}p_{i}p_{\xi})$, whose determinant gives $g_{0}^{2}$,

	$\displaystyle U(T)$	$\displaystyle=$	$\displaystyle\sqrt{-1}(2\pi)^{-N/2}\tau^{-1}|g_{0}q_{i}q_{\xi}|^{-1}|C|^{-1/2}\prod_{i}\left[\frac{1}{2}(\lambda_{i}-\lambda_{i}^{-1})\right]^{-1/2}.$
		$\displaystyle\times$	$\displaystyle\left(\frac{dT}{dE}+\frac{2}{g_{0}^{4}}\sum_{i}(p_{i}g_{0})^{2}\frac{\lambda_{i}+\lambda_{i}^{-1}-2}{\lambda_{i}-\lambda_{i}^{-1}}\right)^{-1/2}$

The integration of the exponent in Eq. 2.13 is performed in a similar manner as in the previous section. However, to handle the divergence associated with the rotational invariance of the tunneling trajectory, one has to recover the original integrals over the trajectories corresponding to the infinitesimal orthogonal rotations and then substitute them with the total solid angle, as described in Sec. 2. This is achieved by using $q_{\xi},\ \xi=x,y,z$ as the basis vectors in the integral in Eq. 2.13. Collecting all the terms, one comes to the following expression,

$$\int^{\prime}d^{N}\!q^{\prime}e^{-S^{\prime}(q^{\prime})}=2(2\pi)^{-N/2}\beta\tau|g_{0}q_{i}q_{\xi}|\prod_{i}\sqrt{\frac{1}{2}\frac{\lambda_{i}+1}{\lambda_{i}-1}},$$

(4.17)

where the prime index on the integral indicates that the divergences have been replaced by the proper values.

Substituting Eqs. 4. and 4.17 into Eq. 2.13 and then into Eq. 1.13, one obtains for $Z^{\#}$,

$$Z^{\#}=4\pi{\beta}F_{tun}\prod_{i}[2\sinh(u_{i}/2)]^{-1}|C|^{-1/2}e^{-S(\beta)},$$

(4.18)

where the tunneling prefactor $F_{tun}$ is given by Eq. 3.12 and the matrix $C_{\xi,\eta}$ is given by Eq. 4.13.