next up previous
Next: Final analysis Up: Theory Previous: Coordinates and Basis set

Hamiltonian and Propagation

The triatomic Hamiltonian in Jacobi coordinates for the H + O2 system in the BF frame is given by[33,49,50,61,62,63,67,68]

\hat{H} = \frac{\hbar^2}{2\mu_R}\frac{\partial^2}{\partial R...
...\mu_R R^2} + \frac{\hat{j^2}}{2\mu_r r^2} + V
\end{displaymath} (7)

where $\mu_R$ and $\mu_r$ are the reduced mass of the H + O2 collision complex and the reduced mass of the O2 molecule, respectively. We propagate the wave function in time using the time-dependent Schrödinger equation

\begin{displaymath}i\hbar \frac{\partial \Psi (R,r,\vartheta;t)}{\partial t} = \hat{H}\Psi (R,r,\vartheta;t).
\end{displaymath} (8)

Using the expansion of the wave function, given in Eq. (2), and taking the inner product with the basis functions we get a set of equations of motion for the coefficients $C_{\lambda \nu j \Omega} (t)$, which are simplified using the orthogonality properties of the Wigner D-matrices. Furthermore, consistent with the DVR approach the integrals in R and r involving multiplicative operators are evaluated within the quadrature approximation, resulting in diagonal matrix representations for these operators. This results in the following set of equations of motion of the coefficient $C_{\lambda \nu j \Omega} (t)$.

$\displaystyle i\hbar \frac{\partial C_{\lambda \nu j \Omega} (t)}{\partial t} =$ - $\displaystyle \sum_{\lambda'}^{N_R} T^{R}_{\lambda \lambda'} C_{\lambda' \nu j
\Omega} (t) -\sum_{\nu'}^{N_r} T^{r}_{\nu \nu'} C_{\lambda \nu' j
\Omega} (t)$  
  + $\displaystyle C_{\lambda \nu j \Omega} (t) \frac{J(J+1) + j(j+1) -2\Omega^2}{2
\mu_R R_\lambda^2}$  
  + $\displaystyle C_{\lambda \nu j \Omega} (t) \frac{j(j+1)}{2 \mu_r r_\nu^2} +
\su...}} \overline{V}^\Omega_{j j'}
(R_\lambda,r_\nu) C_{\lambda \nu j' \Omega} (t)$  
  - $\displaystyle C_{\lambda \nu j \Omega+1} (t) \frac{c^{+}_{j\Omega}
..._{J\Omega} \left( 1+\delta_{\Omega,1} \right)^\frac{1}{2}}{2\mu_R
R_\lambda^2}.$ (9)

$T^{R}_{\lambda \lambda'} = \hbar^2/2\mu_R
\left\langle\phi^{-}_{\lambda} (R) \left\vert \partial^2/ \partial R^2
\right\vert \phi^{-}_{\lambda'} (R)\right\rangle$ and $T^{r}_{\nu \nu'} =
\hbar^2/2\mu_r \left\langle\phi^{-}_{\nu} (r) \left\vert \partial^2/
\partial r^2 \right\vert \phi^{-}_{\nu'} (r)\right\rangle$ can be evaluated analytically[55,56,57,73]. $\overline{V}^\Omega_{j j'} (R_\lambda,r_\nu)$ is the reduced potential in DVR-point $(R_\lambda,r_\nu)$ and is defined as

\overline{V}^\Omega_{j j'} (R_\lambda,r_\nu) = \int_{-1}^{1}...
\Theta^{\Omega}_{j'} (\vartheta) \,d\cos\vartheta.
\end{displaymath} (10)

This expression is evaluated using a Gauss-Legendre quadrature.

The last two terms in the equations of motion are the Coriolis terms and represent the Coriolis coupling in the system. From its form the sparsity of the Coriolis coupling matrix is immediately clear. There is only coupling between the substates $\Omega $ and $\Omega-1$, and the substates $\Omega $ and $\Omega+1$. The coefficients of the form $c^{\pm}_{ab}$ are given as

c^\pm_{ab} = \left[ a(a+1) - b(b\pm1) \right]^\frac{1}{2}.
\end{displaymath} (11)

The time-evolution of the wave function is accomplished by using a symplectic integrator of the (m=6,n=4) type, introduced by Gray and Manolopoulos in 1996 (see Ref. [74]). This means that we take short and discrete time steps $\Delta t$, each of which requires 6 evaluations of $\hat{H}\Psi$.

At the end of each time step we absorb the wave function at the edges of the grid according to $C_{\lambda \nu j \Omega} (t) \rightarrow A^{R}_\lambda
A^{r}_\nu C_{\lambda \nu j \Omega} (t)$ with $A^{R}_\lambda = 1$ for $R\leq
R_{abs}$ and $A^{R}_\lambda = \exp \left[ -B_R
\left(R_\lambda-R_{abs}\right)^2 \right]$ for $R \geq
R_{abs}$[61]. A similar definition holds for $A^{r}_\nu$. The damping coefficients BR and Br are determined by trial and error, where they must be sufficiently strong to dampen the wave function as it approaches the edge of the grid, but not so strong as to cause artificial reflections. This procedure is equivalent to using an imaginary potential to absorb the wave function[75,76]

The propagation of the wave function for J>0 is facilitated by the Coriolis-coupled parallel method, introduced recently by Goldfield and Gray[49]. This method is most efficient on parallel computers with fast connections between the processors, e.g., the IBM SP2 or the Cray T3E. In this method we distribute the portions of the wave function corresponding to different $\Omega $ substates over different processors. Communication between processors is needed only when the Coriolis terms in Eq. (7) are calculated. The overhead generated is negligible compared to the computational work associated with the kinetic/potential terms in Eq. (7). This method requires the same amount of CPU-time per processor as a J=0 calculation, allowing the calculation to finish in approximately the same amount of wall-time. We use the Message-Passing Interface (MPI) library[77,78,79] to perform the communication between the processors. We wish to emphasize here that implementing this method is extremely simple. About 15 lines of Fortran code were needed to convert a serial program into a parallel program. For more details we refer to Ref. [49].

next up previous
Next: Final analysis Up: Theory Previous: Coordinates and Basis set
Anthony J. H. M. Meijer