INTRODUCTION

TOP ABSTRACT INTRODUCTION FORMULATION OF THE PROBLEM SIMULATION ALGORITHM ONE ANISOTROPIC MOLECULE TWO ANISOTROPIC MOLECULES CONCLUDING REMARKS APPENDIX A APPENDIX B REFERENCES

Anisotropic reactivity is a typical feature of the molecules of large size, especially biomolecules. Ligand binding to proteins, enzyme catalysis, self-aggregation of colloids, and polymerization are examples of diffusion-mediated processes involving anisotropically reactive species. Anisotropic reactivity is generally associated with a considerable slowing down of the reaction rate with respect to that for isotropic molecules of the same size (McCammon, 1984). To what extent the rate will be reduced is a nontrivial question that cannot be accounted for on the basis of simple geometrical arguments.

Different models of anisotropic reactivity have been suggested. The most popular one is the model of reactive patches (MRP) (Solc and Stockmayer, 1971, 1973; Schmitz and Schurr, 1972; Schurr and Schmitz, 1976; Hill, 1975; Berg and Purcell, 1977; Shoup et al., 1981; Northrup et al., 1984; Temkin and Yakobson, 1984; Pritchin and Salikhov, 1985; Berg, 1985; Lee and Karplus, 1987; Zwanzig and Szabo, 1991; Traytak, 1994, 1995). It assumes that a portion of the surface of a molecule is reactive while the rest of the surface is inert. The reaction occurs upon diffusive encounter of two molecules with favorable orientations. Usually, to make the problem tractable, the molecules are assumed to be spherical and the reactive patches are assumed to be circular. Although this formulation seems too idealistic from the practical point of view, it still poses quite a complicated mathematical task, even in the case where only one molecule is anisotropic while the other is uniformly reactive. Exact expression for the reaction rate in this case has been obtained in the limit of small angular size of the patch (Hill, 1975; Berg and Purcell, 1977; Shushin, 1986). In particular, it was found that the rate constant is proportional to the radius of the patch instead of its square, as one would expect from a naive geometrical consideration. To be able to analyze analytically the whole range of angular sizes and to consider even more complicated problems, such as the case of two anisotropic molecules, approximations have been introduced. The most actively used are the closure approximation (Wilemski and Fixman, 1973; Doi, 1975) and the constant flux approximation (Shoup et al., 1981), which are, in fact, proved to be equivalent. The expressions for the rate constant derived in these approximations are still fairly cumbersome and demand certain numerical efforts. Considerable simplification without essential loss of accuracy can be obtained by using further the so-called quasi-chemical approximation, proposed by Solc and Stockmayer (1973) and elaborated upon by Berg (1985).

The MRP offers a clear view on many important features of the manifestation of anisotropic reactivity in the kinetics of diffusion-controlled reactions. However, the assumption of the patch-like reactive regions is not always realistic. The reactivity typically depends on distance between certain reactive centers (atoms or molecular fragments) on surfaces of molecules. Therefore, the reactivity regions are actually expected to be convex rather than of patch-like shape. Importance of the effect of "thickness" of the reaction zone on the diffusion-controlled reaction rate is often overlooked. A particular realization of such a situation is found in the model of reactive hemispheres (MRH) (Shushin, 1986, 1988a). It assumes that a molecule is reactive within a certain distance from a point reactive center on its surface. In other words, the reactive region is hemispherical. The reaction occurs whenever the reactive hemispheres of two different molecules overlap.

A powerful method of local analysis has been developed allowing for rigorous analytical treatment of quite realistic models of anisotropic reactivity and interaction, including both the MRP and the MRH, in the practically most interesting limit of strong anisotropy, where the size of the reactive region is much smaller than the size of the molecule (Shushin, 1986, 1988a, 1988b, 1999). The method is based on a simple observation that, in this limit, the reactive flux is determined by the specific features of relative motion of the molecules in the vicinity of reactive orientations. It was found that the MRH predicts much higher reaction rates than the MRP for reactive regions of similar size. The difference was shown to be especially large in the case of two anisotropic molecules (as much as orders of magnitude for typical values of parameters) (Shushin, 1988a, 1988b, 1999, 2000).

Simultaneously with intensive analytical analysis of different models of anisotropic reactivity, much effort has been devoted to the development of efficient numerical methods aimed to test the existing theories and to explore realistic systems, where analytical solution is not feasible. Among a variety of simulation techniques, the Brownian dynamics (BD) approach is found to be particularly suitable for studying the kinetics of diffusion-influenced reactions. Northrup et al. suggested a BD method to extract the stationary rate constant using the capture probability of one reactant molecule started on a spherical surface enclosing the other reactant molecule (Northrup et al., 1984; McCammon et al., 1986; Luty et al., 1992). Lee and Karplus proposed an alternative BD method based on calculating the time-dependent probability that a pair of reactant molecules, having started in the reaction zone, will be found again in this zone in the absence of reaction (Lee and Karplus, 1987; Yang et al., 1999). Zhou developed an efficient algorithm to calculate the time-dependent rate coefficient via the survival probability of the pair of reactants started in the reaction zone (Zhou, 1990, 1993). He also suggested a method based on the mean residence time in the reaction zone (Zhou, 1998). A weighted-ensemble BD algorithm has been recently proposed by Huber and Kim (1996), which enables one to more efficiently sample the reactive trajectories. Another algorithm is based on averaging the reaction functional over nonreactive BD trajectories in the path-integral sense (Lamm and Schulten, 1983; Zhou and Szabo, 1996a, 1996b). All these methods have recently been compared in practical applications, and their merits are now well understood (Zhou and Szabo, 1996b; Zhou, 1998; Yang et al., 1999).

This paper is concerned with a comparative analysis of the MRP and the MRH. In particular, the accuracy of the known and newly derived simple analytical expressions for the stationary diffusion-controlled reaction rate constant K of neutral spherical molecules is tested by comparison with numerical results obtained using the original BD method. The simulation scheme is based on the exact asymptotic formula for the long-time behavior of the pair survival probability S(t) (Shushin, 1988b, 1999). We found this method to be most convenient for the purpose of calculating K.

Comparison with the simulation results for the MRP demonstrates that the quasi-chemical approximation provides sufficient accuracy, and there is no practical need to resort to more elaborate approaches. It is also shown that simple formulas derived for the MRH, using the method of local analysis, are quite accurate in the limit of strong anisotropy corresponding to small size of the reactive regions. To reproduce the numerically calculated MRH rates in a wide range of model parameters, simple interpolation formulas in the spirit of the quasi-chemical approximation are proposed. Numerical calculations confirm earlier predictions that the MRP rates are much smaller than the MRH rates for the same radii of the reactive regions, especially in the case where both reacting molecules are anisotropic.

	FORMULATION OF THE PROBLEM

TOP ABSTRACT INTRODUCTION FORMULATION OF THE PROBLEM SIMULATION ALGORITHM ONE ANISOTROPIC MOLECULE TWO ANISOTROPIC MOLECULES CONCLUDING REMARKS APPENDIX A APPENDIX B REFERENCES

General equations

We consider a bimolecular reaction between two spherical molecules A and B of radii R_A and R_B undergoing translational diffusion with a relative diffusion coefficient D = D_A + D_B and rotational diffusion with the coefficients D_rA and D_rB, respectively. Configuration of the reacting pair can be fully described by the relative position r and individual orientations _A and _B, combined into a general configuration vector Q = (r, _A, _B). The steady-state kinetics and, in particular, the rate of the considered bimolecular reaction is determined by the pair distribution function (Q), which satisfies the Smoluchowski equation

(1)

Here the operator  represents multidimensional gradient in the curvilinear coordinates Q, is the diffusion matrix expressed in terms of the diffusion coefficients D, D_rA and D_rB (see below), and u(Q) = U(Q)/k_BT is the intermolecular mean force potential. In this work, we consider reactions without interaction, i.e., with u(Q) = 0.

The function (Q) is defined outside the surface of closest approach S(Q), whose small part S_r(Q) corresponds to reactive orientations. For spherical molecules, S(Q) is determined by r = |r| = R = R_A + R_B. Eq. 1 should be solved with the following boundary conditions:

(2)

(3)

where S_n = S  S_r denotes the nonreactive part of S and n is the outward vector normal to S. Hereafter, the reactivity at contact (Q) is assumed to be strong, corresponding to the diffusion-controlled limit.

The reaction rate constant K is expressed in terms of the steady-state reaction flux J,

(4)

where K_s = 4RD is the Smoluchowski rate constant for uniformly reacting molecules, and f is the so-called steric factor. It is important to note that only the rate constant itself is a real physically sensible kinetic parameter, whereas the steric factor is generally not. The latter can be uniquely defined only in some particular (MRP-like) models of reactivity and for spherical shape of the molecules (see below). The ambiguity of f results from that of K_s, which, for nonspherical molecules, depends not only on translational but also on rotational diffusion coefficients. Nevertheless, the steric factor defined in Eq. 4 appears to be very helpful in the analysis of the problem at hand, especially in constructing interpolation formulas.

Models of reactivity

The shape of the reactive surface S_r, which determines the value of the reaction rate, depends on the model of reactivity. Here we consider two models: the model of reactive patches (MRP) and the model of reactive hemispheres (MRH), as illustrated in Fig. 1.

View larger version (18K): [in this window] [in a new window]   FIGURE 1   Schematic picture of two models of anisotropic reactivity, the MRP and the MRH. In the MRP, the surface of a molecule is reactive over a circular patch with the rest of the surface being inert. In the MRH, a molecule is reactive within a certain distance from a point reaction center on its surface.

The model of reactive patches

The model of reactive hemispheres

	SIMULATION ALGORITHM

TOP ABSTRACT INTRODUCTION FORMULATION OF THE PROBLEM SIMULATION ALGORITHM ONE ANISOTROPIC MOLECULE TWO ANISOTROPIC MOLECULES CONCLUDING REMARKS APPENDIX A APPENDIX B REFERENCES

Here we propose a simulation scheme that we found to be most convenient for calculating the stationary diffusion-controlled reaction rate constant K. The method is based on the exact asymptotic formula for the long-time behavior of the pair survival probability S(t) (Shushin, 1988b, 1999),

(5)

that is valid for t R²/D for any model of anisotropic reactivity, any shape of molecules, and any model of relative diffusion (except for the case of frozen rotation of both molecules). The parameter S() defines the total escape probability and depends on the initial condition. It should be emphasized that S(t) is independent of R. This fact can be considered as an additional evidence that Eq. 5 is valid for molecules of any shape. The above universal asymptotic behavior of the survival probability entails the same universal behavior of the time-dependent rate coefficient (Shushin, 1988b, 1999; Zhou and Szabo, 1996b).

Numerically, it is more efficient to initiate the trajectories near the reaction region and to take larger time steps as the separation between the reactants is increased. The molecules were moved using standard BD algorithms. The trajectories were terminated either when the reaction criteria were satisfied or when a sufficiently large cutoff time was exceeded. For each set of parameters, at least 2 · 10⁴ trajectories were run and the asymptotic decay kinetics of the resulting average survival probability was fit to Eq. 5, as illustrated in Fig. 2. Technical details can be found in Appendix A.

View larger version (18K): [in this window] [in a new window]   FIGURE 2   Normalized time-dependent survival probability S(t)/S()  1 as a function of R/(Dt)1/2 obtained by averaging 2 · 104 BD trajectories for the MRH with RA = RB, DrAR/D = DrBR/D = 0.75, and lA/RA = lB/RB = 0.2 (solid curve). Dotted line is the result of an asymptotic least-squares fit to Eq. 5.

	ONE ANISOTROPIC MOLECULE

TOP ABSTRACT INTRODUCTION FORMULATION OF THE PROBLEM SIMULATION ALGORITHM ONE ANISOTROPIC MOLECULE TWO ANISOTROPIC MOLECULES CONCLUDING REMARKS APPENDIX A APPENDIX B REFERENCES

For definiteness, we assume that molecule A is anisotropic, whereas molecule B is uniformly reactive.

The model of a reactive patch

The reactive patch on the surface of A is assumed to be circular with the angular radius _A = ₀. For simplicity of notation we also denote D_rA  D_r. The resulting mixed boundary-value problem can be reduced to a linear matrix equation of infinite order (see Appendix B), which can be solved numerically to any degree of accuracy by iteration or truncation (Traytak, 1994, 1995). We used this method as an exact reference. A relatively simple and quite accurate expression for the rate was derived within the closure (constant flux) approximation (Shoup et al., 1981). An even simpler but still sufficiently accurate formula suitable for practical applications is provided by the quasi-chemical approximation (Solc and Stockmayer, 1973; Berg, 1985). We are not going to analyze this approximation in detail but restrict ourselves to outlining its idea and presenting the final result.

The idea of the quasi-chemical approximation is based on the assumption that the reaction kinetics is exponential and that reactive and non-reactive orientations can be treated as two different states of the system. Diffusive reorientation is described as transitions between these two states with the rates proportional to the statistical weights of the corresponding types of orientations. The expression for the combination of the transition rates in the final formula for the overall reaction rate is obtained by comparison with that derived in the closure (constant flux) approximation. An especially simple and accurate formula is obtained with the use of the interpolation expression proposed for this combination by Berg (1985). In so doing, one arrives at the following result for the steric factor:

(6)

where  = . In what follows, we will use this formula for comparison with the MRH predictions and simulations.

In the limit of small reactive patch, ₀ 1, Eq. 6 reduces to

(7)

which agrees quite well with the exact formula (Berg and Purcell, 1977; Shushin, 1986),

(8)

In the opposite limit of isotropic reactivity, ₀ = , we obtain from Eq. 6 f = 1, as expected.

The model of a reactive hemisphere

The case of one anisotropic molecule in the MRH has not been analyzed so far. Therefore, we present some details of this analysis here. The general problem is mathematically very complicated. Fortunately, in the most interesting limit of strongly anisotropic reactivity (l_A R_A) it can be markedly simplified by using the method of local analysis (Shushin, 1986, 1988a). This method allows one to reduce the general Eq. 1 in curvilinear coordinates to that in the local Cartesian coordinates in the vicinity of reactive orientations. In the case of one anisotropic molecule, the local Cartesian coordinates are the distance between the surfaces of the molecules, x = (r  R)/R, and two angular deviations ₁ and ₂ of the center of the reactive hemisphere (on the surface of molecule A) from the intermolecular axis r (see Fig. 3) (Shushin, 1999). Here we define the reaction radius as l_r = l_A to have a unified notation with the case of two anisotropic molecules.

View larger version (43K): [in this window] [in a new window]   FIGURE 3   The coordinates describing relative orientation of molecules in the case of one anisotropic molecule (A).

First, we need to obtain the shape of the reaction surface S_r. The formula for S_r is represented as a dependence of x on  =  (the surface S_r is axially symmetric, i.e., x depends only on the absolute value  of the angular deviation vector). An equation for x() can be easily written when taken into account that the surfaces of B and the reactive hemisphere on A are in contact if the distance l between the center of B and the reactive center on A is given by l = R_B + l_r (see Fig. 3), i.e.,

(9)

where

(10)

For relatively large l_r  R_A the shape of S_r is quite peculiar, but it can be characterized by two length parameters, the height _n (corresponding to  = 0) and the radius _t of the circular base (corresponding to x = 0),

(11)

(12)

In the limit of strongly anisotropic reactivity, l_r R_A (but for arbitrary R_B), the shape reduces to ellipsoidal,

(13)

S_r is a segment of this ellipsoid (x > 0) with the height _n and the radius _t of the circular base given by

(14)

Eq. 14 shows that, for molecules B of small radius R_B < l_r R_A, the surface S_r is a hemisphere with _t  n = l_r/R. However, in the opposite limit of large molecule B, R_B  R_A, the surface S_r is of a patch-like shape with the base radius much larger than the height, _t  [(2R_B/R_A)(l_r/R)]^1/2 _n = l_r/R. As the radius R_B is increased, the height _n decreases, tending to 0 while _t monotonically grows approaching the asymptotic value of  = . These estimations show that the characteristic size of the reaction region S_r increases from l_r/R to l_r/R as R_B/R_A is increased from R_B/R_A 1 to R_B/R_A 1.

It is already clear from the above analysis that, for molecules of comparable size (R_A ~ R_B), the MRH rate will be much larger than the MRP rate corresponding to the reactive patch of angular size of the reactive hemisphere (~l_r/R_A). This interesting and important effect is similar to that predicted for the case of two anisotropic molecules (Shushin, 1988a, 1988b) (see below). Although the assumption of uniform reactivity for a large molecule is, perhaps, not quite realistic, discussion of this case is very instructive and useful for further analysis because the origin of the rate enhancement effect is especially clear.

Relative diffusion of molecules A and B is described as diffusion of a point particle in the local coordinates (₁, ₂, x). The corresponding diffusion matrix is diagonal in these coordinates with

(15)

To obtain a fairly accurate estimation for the rate, valid for any relation between _n and _t, we can approximate S_r by an oblate ellipsoid with semiaxes _n and _t. In so doing, we get, for the steric factor,

(16)

where

(17)

and  = , as before. In the limit of _{n t}, corresponding to the patch, one can put    and obtain f = _t/ (Shushin, 1986, 1999).

For not very strong anisotropy of reactivity (l_r  R_A, R_B) the size _t of the base of the reaction region appears to be fairly large: _t ~ 1. In this case, one can estimate the rate by using the quasi-chemical approximation for an effective patch of radius r_eff = _tR. It is expected to slightly underestimate the rate because it does not take into account the thickness of the reaction region. However, the effect of thickness turns out to be not so strong for realistic values of D_rR²/D ~ 1 (as long as R_B is large).

Another simple way to extend the theory beyond the strict limit of strong anisotropy for arbitrary R_B is by drawing a correspondence between the MRP and the MRH. Thus, the effective patch should be chosen to give the same steric factor as predicted by the method of local analysis (MLA) for the MRH, i.e.,

(18)

where f_MLA is given by Eq. 16. _eff is then used in Eq. 6 instead of ₀.

By definition of the steric factor, the rate K is referred to K_s = 4RD, which is the reaction rate for uniformly reacting molecules. This formula for K_s is evidently valid in the MRP but it is generally not in the MRH. Clearly, as soon as l_r > R, the rate constant K_s  4Dl_r. The simplest expression that combines the two limits is K_s = 4D(R + l_r). It should be noted, however, that, in principle, the rate K_s also depends on the ratio D_rR²/D. We are not going to discuss this problem in more detail because, first, this dependence does not seem to be very strong, and, second, in the most realistic case of D_rR²/D  1 and l_r  0.5R considered here, the deviation of K_s from 4RD is fairly small and can be neglected. The simplest approximation which takes this effect into account yields,

(19)

where f_QC is given by the quasi-chemical Eq. 6.

Comparison with numerical results

The first thing we have to do is to make sure that our BD algorithm actually works. Exact solution for the MRP can serve as a reference. Figure 4 shows that numerical simulations indeed reproduce quite well the dependence of the steric factor on the angular size of the patch for different values of D_rR²/D. The parameter D_rR²/D describes the influence of the rotational motions of A. Assuming that friction is of Stokes type (with stick boundary conditions) we have

(20)

i.e., this is just a geometric factor. Thus, for molecules of comparable size, D_rR²/D ~ 1, whereas, for R_B R_A, such as in the case of ligand binding to proteins, D_rR²/D  0. Our choice of parameters will be based on these arguments.

View larger version (16K): [in this window] [in a new window]   FIGURE 4   BD simulation results (circles) versus the exact solution (lines) for the effective steric factor in the case of one anisotropic molecule bearing an MRP site as a function of angular size of the patch for DrR2/D = 3 (upper curve) and 0 (lower curve).

The main problem in the comparison of the MRP and the MRH is to relate their size parameters. Although the reaction radius l_r in the MRH can be readily estimated for certain models of distance-dependent reactivity (see the preceding section), the angular size ₀ of the patch in the MRP is not directly related to the reaction mechanism. In what follows, we assume a reasonable relation,

(21)

i.e., the reaction radius l_r in the MRH is assumed to be equal to the radius r_A of the patch in the MRP.

Figure 5 compares the steric factors f = K/(4RD) calculated for the MRP and the MRH against the numerical results for R_A = R_B and R_A = 10³R_B (D_rR²/D = 1.5 and 0, respectively). It is seen that the MRH rate becomes significantly larger than the MRP rate as the ratio of R_B/R_A is increased. This effect can be attributed solely to the difference in the area of the reactive sites only in the limit of R_B/R_A 1. 

View larger version (19K): [in this window] [in a new window]   FIGURE 5   A plot of the BD simulation results for the effective steric factor as a function of the reaction radius in the case of one anisotropic molecule, comparing the MRP (triangles) versus the MRH (circles) for RA = RB, DrR2/D = 1.5 (top) and RA = 103 RB, DrR2/D = 0 (bottom). Dashed lines correspond to Eqs. 7 and 16 obtained in the limit of strong anisotropy for the MRP and the MRH, respectively. The correspondence between the MRH- and MRP-size parameters of Sr is established by Eq. 21. Solid lines represent the quasi-chemical approximation, Eqs. 6 and 19. In the case of the MRH, the effective patch size is used as explained in the text.

The accuracy of our simple limiting expression, Eq. 16, for the MRH proves to be quite reasonable in the limit of strongly anisotropic reactivity, l_r R. Further improvement of the theory is achieved by using the quasi-chemical interpolation formula, Eq. 19. Good accuracy of the quasi-chemical approximation in the MRP is, in principle, well established (Berg, 1985; Zhou, 1993). Here we just emphasize once again its practical power as applied to the MRH. The effective angular radius _eff (see Eq. 18), which determines the MRP representation of the MRH steric factor f_MRH in Eq. 19, can significantly exceed ₀ = 2 arcsin(l_r/2R) for molecules of comparable size. Note that the steric factor for the MRH can be larger than 1. This is because the effective contact radius becomes larger than R as soon as l_r and R become comparable, as discussed in the preceding subsection.

	TWO ANISOTROPIC MOLECULES

TOP ABSTRACT INTRODUCTION FORMULATION OF THE PROBLEM SIMULATION ALGORITHM ONE ANISOTROPIC MOLECULE TWO ANISOTROPIC MOLECULES CONCLUDING REMARKS APPENDIX A APPENDIX B REFERENCES

The model of reactive patches

No exact analytical expressions are available for the rate constant in the case of two anisotropic molecules. Nevertheless, quite accurate formulas are derived in the closure (constant flux) approximation (Temkin and Yakobson, 1984; Lee and Karplus, 1987). Unfortunately, these formulas are fairly cumbersome and thus not very much convenient for applications. Considerable simplification without essential loss of accuracy is obtained by using the quasi-chemical approximation (Solc and Stockmayer, 1973; Berg, 1985). The basic idea underlying this approximation was briefly discussed in the preceding section. It allows one to express the rate through the parameters obtained for the case of one reactive patch (Berg, 1985):

(22)

where F₀ = F_AF_B, F_v = sin²(_v/2) (v = A, B) are the geometric steric factors,

(23)

where _v = , and

(24)

In the limit of strongly anisotropic reactivity, _A,B 1, Eq. 22 reduces to a simple formula:

(25)

In the opposite limit of isotropic reactivity, _A = _B = , Eq. 22 gives f = 1, as it should.

In the limit when one of the molecules is much larger than the other (R_B R_A, i.e., the surface of B approaches a plane), f is not an appropriate kinetic parameter because K_s = 4DR becomes infinite, though the rate constant K naturally retains its physical meaning. In this case, the reduced rate constant is conveniently defined in terms of a new dimensionless steric factor f_A,

(26)

By taking the limit R_B  , one obtains, in the quasi-chemical approximation,

(27)

where _B = r_B/R_A is the dimensionless parameter, determined by the radius r_B = R_BB of the reactive patch on the (plane-like) surface of the molecule B, and _A = (D_rAR/D_A)^1/2. For strongly anisotropic molecule A, i.e., for _A 1, Eq. 27 simplifies into

(28)

The model of reactive hemispheres

Approximate formulas for the reaction rate in the MRH have been derived so far only in the limit of strongly anisotropic reactivity (l_r R_A, R_B) using the method of local analysis (Shushin, 1986, 1988a). The details of the derivation have been thoroughly considered in a number of articles (Shushin, 1988a, 1999, 2000). Therefore, here we restrict ourselves to a brief discussion of the principal ideas of the method and present only the final analytical expression, which is to be compared with the simulation results.

In the MLA, one considers relative diffusive motion of the molecules in close vicinity of the reactive orientations, where the system configuration can be represented as a point in five-dimensional space ₅ of Cartesian coordinates. These coordinates are the normalized distance between the surfaces of the molecules, x_N = (r  R)/R, and the angular deviations of the reactive centers from the intermolecular axis r, i.e., the normalized projection vectors x_v = z_v/R_v (v = A, B) of the reactive centers onto the plane perpendicular to r, as shown in Fig. 6.

View larger version (21K): [in this window] [in a new window]   FIGURE 6   The coordinates describing relative orientation in the case of two anisotropic molecules.

The shape of the reaction region S_r is determined by the condition that the distance l between the reactive centers is equal to the reaction radius, l = l_r  l_A + l_B, leading to

(29)

where _r and µ_v are defined in Eq. 10.

In the limit of strong anisotropy of reactivity (_r 1), a convenient approximate expression for S_r can be obtained in coordinates

(30)

(31)

(32)

where  = . We have

(33)

Analysis of this expression shows that the shape of S_r is close to ellipsoidal with semiaxes

(34)

and
The projection of S_r described by Eq. 33 is shown schematically in Fig. 7.

View larger version (18K): [in this window] [in a new window]   FIGURE 7   Schematic picture of the reaction regions in the MRP (small rectangle) and the MRH (ellipsoid). Big square represents the projection of the boundary in the subspace {xA, xB}.

In X-coordinates the diffusion matrix is nondiagonal with the following diagonal elements:

(35)

(36)

(37)

It can be rigorously shown that nondiagonal elements only weakly affect the absolute value of the reaction rate and thus can be neglected (Shushin, 2000). In fact, nondiagonal diffusion matrix elements are also neglected in the quasi-chemical approximation. The limit of strong anisotropy corresponds to _r 1, and thus _{2 1}, _N. Therefore, one can use the so-called adiabatic approximation and neglect the effect of diffusion along the coordinates X₂, i.e., put D₂₂ = 0 (Shushin, 1988a, 1999). By further assuming, for simplicity, that the reaction region S_r is spherical in the subspace {X₁, X_N} with the radius  = _r  1 (this approximation leads to a minor underestimation of the rate) one finally obtains (Shushin, 1988a, 1999, 2000)

(38)

where

(39)

Eq. 38 predicts the principally different dependence of the reaction rate on the size of the reactive sites than the MRP. In the MRP we had f ~ l for equal reactive sites, while f ~ l in the MRH.

As in the case of one anisotropic molecule, we can extend the theory beyond the limit of strong anisotropy by using the quasi-chemical type of approximation. First of all, we notice the similarity of the reaction regions and the diffusion matrices in the MRH and the MRP. The reaction region S_r in the MRP is a patch (the height is _N = 0) of cylinder shape in ₅ space whose projections in the subspaces x_A and x_B are circles of radii _v = _v (v = A, B), whereas the projections onto the planes {x, x} (the superscripts i, k = 1, 2 correspond to the coordinates in the two-dimensional subspaces x_A and x_B) are rectangles, as pictured schematically in Fig. 7 (Shushin, 2000). The diffusion matrix used in the quasi-chemical approximation is diagonal in x-coordinates and the diagonal elements are given by

(40)

It is easily seen that, by changing S_r of Eq. 33 to the cylinder patch with the corresponding size parameters of Eq. 34 (but with _N = 0), the MRH formulation of the problem becomes identical to that in the MRP quasi-chemical approximation. The only difference is in the orthogonal transformation (rotation) of variables in Eqs. 30 and 31. This kind of modification of S_r is known to affect the absolute value of the rate only slightly (Shushin, 1999, 2000). Therefore, the MRH rate should be well approximated by the quasi-chemical expression of Eq. 22 after the corresponding transformation of the size parameters of S_r and the elements of the diffusion matrix,

(41)

in which D_jj are defined in Eqs. 35-37. The rate is expected to be slightly underestimated in this way because the thickness of the reaction region was neglected.

Again, as in the case of one anisotropic molecule, the reaction rate constant for two uniformly reacting MRP molecules is given by K_s = 4RD only in the limit of strong anisotropy, l_v R_v (v = A, B), when the area of the reactive surface S_r in the four-dimensional angular space {x_A, x_B} is much smaller than (4)². The correction is expected to be ~l_r/R. In the opposite limit of l_v > 2R_v, the molecules are uniformly reactive, though the centers of the reactive spheres are displaced from the centers of (spherical) molecules. This displacement affects the rate only for not very large l_v  2R_v. In the limit of l_v 2R_v this effect is negligible and K_s  4l_rD. These two estimations for K_s can be combined in a simple interpolation expression K_s = 4D(R + l_r), which is expected to be fairly accurate for relatively slow rotational diffusion, i.e., for D_rA,B  D/R², typical of real systems.

Comparison with numerical results

Figure 8 compares the steric factors f = K/(4DR) for R_A = R_B and f_A = K/(4DR_A) for R_B R_A calculated in the MRP and the MRH against the numerical results for l_A = l_B = l_r/2 (recall that, in the MRH, the rate depends only on l_r = l_A + l_B and the relation l_A = l_B is taken only for definiteness). We assume, as before, that the reaction radius in the MRH is equal to the radius of the patch in the MRP,

(42)

It is seen that the MRH rate is orders of magnitude larger than the MRP rate when l_r R. Analytically predicted scaling behaviors, i.e., f ~ l for the MRH and f ~ l for the MRP are well confirmed numerically. In fact, the accuracy of the simple limiting expressions proves to be quite reasonable ev