In this appendix we explain in more detail the method for the evaluation of and the method we use to delete points from the grid.

In our program the coefficient vector is ordered with the
index running fastest. Thus, the application of
[whose
elements are defined just below Eq. (7)] is
straightforward.
is not a function of
and *j*,
therefore we store this matrix only once and apply it to the
coefficient vector for each
combination. Applying
to the coefficient vector is *not* straightforward, because
the coefficient vector is ordered with
running fastest and
not with
running fastest. Hence, prior to applying
the coefficient vector is reordered with the
index running
fastest. Afterwards, the resultant vector is reordered back to the
original ordering and added to the resultant vector of the
operation. Also here, we store the matrix
only
once. The same procedure is followed for the application of
and the rotational terms
from Eq. (7), but now we reorder with *j* running fastest.
Of course, actually resorting the vector each time is costly.
Therefore, we only sort the vectors once and store the permutations
needed to get from one ordering to another (and back again). This
procedure is taken from Ref. [57]. Any general
purpose routine can be used in this algorithm. We used the
HEAPSORT algorithm from Numerical Recipes[90]. After
calculating the kinetic energy terms, the potential terms and the
rotational terms, we then calculate the Coriolis terms for substate
by first getting the coefficient vector for substate
from the neighboring processor (using MPI[77]) and
subsequently calculating the Coriolis coupling between substates
and .
Then we do the same thing for the Coriolis
coupling with substate .
At this point we cycle back until we
have completed one time-step.

The calculation of the kinetic/potential terms of the Hamiltonian are the most expensive steps, computationally. The storage of the reduced potential matrix is most expensive memory-wise and we do not want to recalculate this matrix at each time-step, since its calculation is also time-consuming. With the parameters given in Table I, we need approximately 180 Mb to store the matrix. However, a significant number of grid points lie in regions, where the potential is such that the wave function will be zero. Therefore, we have devised a way to get rid of the points in those areas.

The algorithm is simple. We determine, for each combination, whether it lies the asymptotic region of 3 separate atoms or close to the repulsive wall. If either of the two situations is the case, we delete the point from our grid. The way we determine where a point lies on the PES is illustrated by Fig. 6. In this figure we plot the angular dependence of the PES at certain characteristic combinations. Curve 1 lies in the region of three separate atoms. There the potential is always lower than the asymptote and almost flat. For curve 2, which lies close to the repulsive wall, the potential is high for all angles. Thus, we discard points which have certain characteristics.

- 1.
- Separated atoms:

The potential for all angles at a certain combination is lower than the asymptote for 3 separate atoms*and*it is always higher than a cut-off energy*V*_{cut}^{(o)}, which is typically 4 eV above the H + O_{2}entrance channel (see curve 1 in Fig. 6). - 2.
- Repulsive wall:

For all angles the potential is higher than a cut-off energy*V*_{cut}^{(i)}which is typically 8 eV above the H + O_{2}entrance channel (see curve 2 in Fig. 6).

Using this algorithm only slightly complicates the algorithm for the
calculation of
.
In fact, it only has an impact on the
application of
and
to the coefficient
vector, because deleting a
point means deleting a
row and column from
and
.
However, this is
straightforward to implement. As a result of deleting points, one
could end up with a row of DVR-points which is cut in two (see
Fig. 7 around *R*=7 a.u.). This means in this case a kinetic
energy matrix
,
which consists of two diagonal blocks
corresponding to the two sections of the row of DVR-points and two
rectangular off-diagonal block connecting the two sections. At this
point we can choose to ignore those off-diagonal
blocks[57] or use them in the calculation. We examined
both options and conclude that there is not much difference between the
two. In the end, we choose to include the off-diagonal blocks in the
calculation with the rationale that although we assume that the wave
function is zero in the gap between the points, it might still be that
there is overlap between the wave function on one side and the second
derivative of the wave function on the other side of the gap.