# Persistent effects of inertia on diffusion-influenced reactions: Theoretical
methods and applications

Sangyoub Lee<EMAIL_ADDRESS>Department of Chemistry, Seoul National
University, Seoul 08826, South Korea  Sergey D. Traytak<EMAIL_ADDRESS>Semenov Federal Research Center for Chemical Physics, Russian Academy of
Sciences, 4 Kosygina St., 119991 Moscow, Russian Federation  Kazuhiko Seki
<EMAIL_ADDRESS>National Institute of Advanced Industrial Science and
Technology (AIST), Onogawa 16-1 AIST West, Ibaraki, 305-8569, Japan

###### Abstract

The Cattaneo–Vernotte model has been widely studied to take momentum
relaxation into account in transport equations. Yet, the effect of reactions
on the Cattaneo–Vernotte model has not been fully elucidated. At present, it
is unclear how the current density associated with reactions can be expressed
in the Cattaneo–Vernotte model. Herein, we derive a modified Cattaneo–Vernotte
model by applying the projection operator method to the Fokker–Planck–Kramers
equation with a reaction sink. The same modified Cattaneo–Vernotte model can
be derived by a Grad procedure. We show that the inertial effect influences
the reaction rate coefficient differently depending on whether the intrinsic
reaction rate constant in the reaction sink term depends on the solute
relative velocity or not. The momentum relaxation effect can be expressed by a
modified Smoluchowski equation including a memory kernel using the
Cattaneo–Vernotte model. When the intrinsic reaction rate constant is
independent of the reactant velocity and is localized, the modified
Smoluchowski equation should be generalized to include a reaction term without
a memory kernel. When the intrinsic reaction rate constant depends on the
relative velocity of reactants, an additional reaction term with a memory
kernel is required because of competition between the current density
associated with the reaction and the diffusive flux during momentum
relaxation. The competition effect influences even the long-time reaction rate
coefficient.

††preprint: AIP/123-QED

##  I Introduction

In the conventional theory of diffusion-influenced reactions,Rice (1985)
diffusive motion under a potential is expressed by the Smoluchowski equation
and the reaction is taken into account by setting proper boundary conditions
or adding a reaction sink term to the Smoluchowski equation. Rice (1985) The
latter approach is inevitable for taking into account long-range reactions,
such as energy-transfer, and electron transfer in the Marcus inverted region.
Wilemski and Fixman (1973); Murata and Tachiya (1996); Seki, Barzykin, and
Tachiya (1999) For localized reactions, the effect of the reaction can be
taken into account either by setting a boundary condition or introducing a
reaction sink term with the reflecting boundary condition at the contact
distance. Collins and Kimball (1949); Shoup and Szabo (1982) Although the
conventional approach has been successful for studying most reactions in
condensed phases, the more detailed motion of solutes should be taken into
account to study reactions at shorter timescales compared with the solute
velocity correlation time or to study collision-induced reactions. Traytak
(2023); Lee and Traytak (2023)

The most straightforward method to assess the short time and/or the effect of
collision-induced reactions can be to perform Langevin dynamics simulations,
molecular dynamics simulations, or Monte-Carlo simulations. Tachiya (1986);
Dong, Baros, and Andre (1989); Zhou and Szabo (1991); Van Beijeren, Dong, and
Bocquet (2001); Yang, Han, and Lee (2001); Litniewski and Gorecki (2004); Lee
_et al._ (2004); Kim  _et al._ (2009); Piazza, Foffi, and Michele (2013) For
analytical approaches, reactions under ballistic transport at short times and
diffusive transport at long times can be most simply formulated using the
Fokker–Planck–Kramers equation (i.e., the Fokker–Planck equation in phase
space). Kramers (1940); Wang and Uhlenbeck (1945); Doi (1975); Northrup and
Hynes (1978); Naqvi, Mork, and Waldenstrøm (1982); Harris (1983a); Naqvi,
Waldenstrøm, and Mork (1983); Molski (1988); Ibuki and Ueno (1997, 2003,
2006); Kim  _et al._ (2009) The formulation using the Fokker–Planck–Kramers
equation is inevitable when a reaction is induced by collisions, where the
intrinsic reaction rate constant for the Fokker–Planck–Kramers equation
depends on the relative velocity of the reactants. Doi (1975); Northrup and
Hynes (1978); Naqvi, Mork, and Waldenstrøm (1982); Harris (1983a); Naqvi,
Waldenstrøm, and Mork (1983); Molski (1988); Ibuki and Ueno (1997, 2003,
2006); Kim  _et al._ (2009) However, setting the boundary conditions to
describe reactions in the Fokker–Planck–Kramers equation is not
straightforward. Harris (1983b, 1981, 1982, a); Naqvi, Mork, and Waldenstrøm
(1982); Naqvi, Waldenstrøm, and Mork (1983); Molski (1988); Ibuki and Ueno
(1997, 2003) Here, we note that a collision-induced reaction can be taken into
account by adding the delta-function reaction sink term into the
Fokker–Planck–Kramers equation.Northrup and Hynes (1978) Systematic
perturbation expansion in terms of the reaction sink term can be carried out.

We present two approaches to eliminating velocity variables from the
Fokker–Planck–Kramers equation with a reaction sink term. One approach is
based on decoupling between the configurational distribution and momentum
distribution, where an equilibrium Maxwell (Gaussian) distribution is assumed
for the momentum. The second approach is based on the projection operator
method with perturbation expansion. We show that the both approaches lead to
the same modified Smoluchowski equation that includes a memory kernel and
reaction sink terms. Using the decoupling approach, we also obtain
differential equations for the density and current density. The differential
equation for the current density can be regarded as the modified
Cattaneo–Vernotte equation, where the inertial effect is modeled as a
Markovian relaxation of the current density. The Cattaneo–Vernotte equation
was modified, with an additional term added to account for collision-induced
reactions.

In the absence of a potential, the standard Cattaneo–Vernotte equation
Cattaneo (1958); Vernotte (1958); Kubo, Toda, and Hashitsume (1991) reads

$\displaystyle\frac{\partial}{\partial t}\bm{j}(\bm{r},t)$
$\displaystyle=-\frac{1}{\tau_{\rm
D}}\left[\bm{j}(\bm{r},t)+D\frac{\partial}{\partial\bm{r}}\rho(\bm{r},t)\right].$
(1)

Hereafter, $\rho(\bm{r},t)$ denotes the non-equilibrium pair correlation
function (the concentration field of one reactant species at the position
$\bm{r}$ at time $t$ around a reactant partner at the coordinate origin
divided by the bulk concentration value) and $\bm{j}(\bm{r},t)$ is the
associated current density. Corresponding arguments are the mutual separation
vector $\bm{r}$ and time $t$ given in the augmented configuration space
${\mathbb{R}}_{\bm{r}}^{3}\times\\{t>0\\}$; $\tau_{\rm D}$ and $D$ denote the
momentum relaxation time and the mutual diffusion constant, respectively.

We show that the Cattaneo–Vernotte equation must be modified by the coupling
with reaction as,

$\displaystyle\frac{\partial}{\partial t}\bm{j}(\bm{r},t)$
$\displaystyle=-\frac{1}{\tau_{\rm
D}}\left[\bm{j}(\bm{r},t)+D\frac{\partial}{\partial\bm{r}}\rho(\bm{r},t)-\tau_{\rm
D}{\bf{R}}_{\rm v}(\bm{r},t)\right],$ (2)

where the explicit expression for the term $\tau_{\rm D}{\bf{R}}_{\rm
v}(\bm{r},t)$, representing the influence of a collision-induced reaction on
the relaxation of the current density, will be given below.

The modification of the Cattaneo–Vernotte equation may be ignored when the
intrinsic reaction rate constant in the reaction sink term of the
Fokker–Planck–Kramers equation is independent of the solute relative velocity.
Casting into the form of the Cattaneo–Vernotte equation clarifies physical
interpretations of the effect of a collision-induced reaction on the current
density. In addition to the Cattaneo–Vernotte equation, we show that the
continuity equation is modified by a reaction sink term,

$\displaystyle\frac{\partial}{\partial t}\rho(\bm{r},t)$
$\displaystyle=-\frac{\partial}{\partial\bm{r}}\cdot\bm{j}(\bm{r},t)-R(\bm{r},t),$
(3)

where the last term indicates an influence due to the reaction rate
$R(\bm{r},t)$.

Note that in our previous paper, Lee and Traytak (2023) we called the system
of Eqs. (2) and (3) "diffusive Cattaneo system". It is noteworthy that Eq. (3)
is an exact equation representing the conservation law, whereas the standard
Cattaneo–Vernotte equation [Eq. (1)] is essentially a constitutive equation,
being an approximation based on linear relaxation to local stationarity.
Performing the known Kac’s trick with the Cattaneo–Vernotte system of Eqs. (2)
and (3), Kac (1974) differentiating Eq. (3) with respect to $t$ and applying
the divergence operator to Eq. (2), we can eliminate $\bm{j}(\bm{r},t)$ and
derive the so-called telegraph equation for $\rho(\bm{r},t)$ alone:

$\displaystyle\tau_{\rm D}\frac{\partial^{2}}{\partial
t^{2}}\rho(\bm{r},t)+\frac{\partial}{\partial t}\rho(\bm{r},t)$
$\displaystyle=D\frac{\partial^{2}}{\partial\bm{r}^{2}}\rho(\bm{r},t)-R(\bm{r},t)-\frac{\partial}{\partial
t}R(\bm{r},t)-\tau_{\rm D}\frac{\partial}{\partial\bm{r}}\cdot{R}_{\rm
v}(\bm{r},t).$ (4)

If the last two terms are not included, only one reaction sink term remains;
the resultant equation has been called the hyperbolic reaction-diffusion
equation. Mendez, Fedotov, and Horsthemke (2010); Ghorai and Bairagi (2022);
Traytak (2023) If the last term is dropped, the equation is called the
reaction–telegraph equation and has been studied. Mendez, Fedotov, and
Horsthemke (2010); Wakou, Gorecki, and Kitahara (1998); Holmes (1993); Tilles
and Petrovskii (2019); Lee and Traytak (2023) We show that the
reaction–telegraph equation is appropriate when the intrinsic reaction rate
constant in the reaction sink term of the Fokker–Planck–Kramers equation is
independent of reactant velocity, such as in cases of electron transfer and
energy transfer. Our derivation can be regarded as a microscopic approach in
which a reaction at the contact distance is considered. The reaction–telegraph
equation has been studied in both mesoscopic and microscopic approaches.
Mendez, Fedotov, and Horsthemke (2010); Wakou, Gorecki, and Kitahara (1998);
Holmes (1993); Tilles and Petrovskii (2019); Lee and Traytak (2023) We
emphasize that the last term in Eq. (4) is inevitable for collision-induced
reactions and influences the long-time rate coefficient.

##  II Fokker–Planck–Kramers equation with a reaction sink term

Assume that the origin of the coordinates is located at the center of a
reactant and $\sigma$ indicates a contact reaction distance. Let us consider
the generalized non-equilibrium pair correlation function
$f(\bm{r},\bm{v},t)$, where $\bm{r}$ indicates the relative configuration
vector between reactants and $\bm{v}$ stands for the relative velocity. Thus
the function $f(\bm{r},\bm{v},t)$ is defined in the augmented phase space
$(\bm{r},\bm{v},t)\in\\{\bm{r}\in{\mathbb{R}}_{\bm{r}}^{3}:|\bm{r}|>\sigma\\}\times{\mathbb{R}}_{\bm{v}}^{3}\times\\{t>0\\}$.
The Fokker–Planck–Kramers equation with a reaction sink term introduced to
incorporate the effects of a reaction reads as

$\displaystyle\frac{D}{Dt}f(\bm{r},\bm{v},t):=\frac{\partial}{\partial
t}f(\bm{r},\bm{v},t)+$
$\displaystyle\left(\bm{v}\cdot\frac{\partial}{\partial\bm{r}}-\frac{1}{\mu}\frac{\partial
U}{\partial\bm{r}}\cdot\frac{\partial}{\partial\bm{v}}\right)f(\bm{r},\bm{v},t)$
$\displaystyle=\left[\frac{\partial}{\partial\bm{v}}\frac{1}{\tau_{\rm
D}}\cdot\left(\bm{v}+\frac{k_{\rm
B}T}{\mu}\frac{\partial}{\partial\bm{v}}\right)\right]f(\bm{r},\bm{v},t)-R(\bm{r},\bm{v},t)\,,$
(5)

where the Fokker–Planck–Kramers operator is separated into a streaming and a
collisional part in parentheses on the left hand side and in square brackets
on the right hand side, respectively. $D/Dt$ denotes the material derivative.
For the subsequent study it is convenient to rewrite Eq. (5) as follows:

$\displaystyle\frac{\partial}{\partial t}f(\bm{r},\bm{v},t)=({\cal
L}_{0}+{\cal L}_{1})f(\bm{r},\bm{v},t),$ (6)

where ${\cal L}_{0}$ is the Fokker-Planck collision operator,

$\displaystyle{\cal L}_{0}$
$\displaystyle=\frac{\partial}{\partial\bm{v}}\frac{1}{\tau_{\rm
D}}\cdot\left(\bm{v}+\frac{k_{\rm
B}T}{\mu}\frac{\partial}{\partial\bm{v}}\right)$ (7)

and ${\cal L}_{1}$ is defined by

$\displaystyle{\cal L}_{1}$ $\displaystyle={\cal L}_{\rm
L}-R(\bm{r},\bm{v},t),$ (8)

using ${\cal L}_{\rm L}$ representing the streaming term:

$\displaystyle{\cal L}_{\rm L}$
$\displaystyle=-\bm{v}\cdot\frac{\partial}{\partial\bm{r}}+\frac{1}{\mu}\frac{\partial
U}{\partial\bm{r}}\cdot\frac{\partial}{\partial\bm{v}}.$ (9)

In Eq. (5), $U(\bm{r})$ is an interacting potential, $\tau_{\rm
D}=\mu/\xi_{r}$ is the momentum relaxation time constant, where $\mu$ and
$\xi_{r}$ are the reduced mass and the reduced friction coefficient,
respectively Lee and Traytak (2023), $k_{\rm B}$ is the Boltzmann constant,
$T$ is the absolute temperature; $R(\bm{v},\bm{r},t)$ is the reaction sink
term, whose explicit expression is given below. We consider isotropic systems,
and introduce the unit normal vector to the reaction surface,
$\bm{n}=\bm{r}/r$.

Consider the case when the reaction does not occur; i.e. in Eq. (5)
$R(\bm{r},\bm{v},t)\equiv 0$. Posing the reflecting boundary condition for Eq.
(5) we integrate it over $\bm{v}$. In this way the continuity equation is
obtained by assuming that $f(\bm{r},\bm{v},t)$ tends to zero exponentially
fast as $|\bm{v}|\rightarrow\infty$. As shown in Appendix A, using the
divergence theorem for the continuity equation, Eq. (5) should be supplemented
with a reflecting boundary condition expressed by

$\displaystyle\int d\bm{v}\bm{n}\cdot\bm{v}f(\sigma-\epsilon,\bm{v},t)$
$\displaystyle=0,$ (10)

when the velocity components are integrated out, where $\sigma$ indicates a
contact distance and $\epsilon$ indicates a small positive value; later, we
take the limit of $\epsilon\rightarrow 0$. $\epsilon>0$ is introduced to avoid
interference between reflection and reaction sink term represented using a
delta-function set at $\sigma$ as shown below. The reflecting boundary is set
at the radius $\sigma-\epsilon$ to ensure that not a half but the whole
magnitude of the reactive sink term can be taken into account by the delta-
function set at $\sigma$; the reflecting boundary condition approaches to
$\sigma$ by taking the limit of $\epsilon\rightarrow 0$. Hence,
$f(\bm{r},\bm{v},t)$ is conserved by setting the reflecting boundary condition
when $R(\bm{r},\bm{v},t)\equiv 0$ and, conversely, the reaction term breaks
the conservation of $f(\bm{r},\bm{v},t)$.

For collision-induced reactions without reflection at a contact distance
$\sigma$ in isotropic systems, the reaction sink term can be written as

$\displaystyle R(\bm{r},\bm{v},t)=\frac{R_{\rm i}(r,\bm{v},t)}{4\pi
r^{2}}\delta(r-\sigma).$ (11)

$R_{\rm i}(r,\bm{v},t)$ is given by Molski (1988)

$\displaystyle R_{\rm
i}(r,\bm{v},t)=-4\pi\sigma^{2}\bm{n}\cdot\bm{v}\theta(-\bm{n}\cdot\bm{v})f(r,\bm{v},t),$
(12)

where $\theta(x)$ is the Heaviside step function; $\theta(x)=1$ for $x\geq 0$
otherwise zero. Equation (12) indicates that inward fluxes characterized by
$v_{z}=\bm{n}\cdot\bm{v}\leq 0$ are perfectly absorbed. Here, the perfectly
absorbing boundary condition is defined in phase space, which should be
distinguished from the perfectly absorbing boundary condition in configuration
space; the density at $\sigma$ is zero by imposing the perfectly absorbing
boundary condition in configuration space. Rice (1985) In phase space, we need
to take into account the ultimate escape from collision if the condition given
by $v_{z}>0$ is satisfied.

For collision-induced reactions with a partially absorbing boundary condition
at a contact distance $\sigma$, we introduce the fraction of reactive flux
$f_{\rm r}(|v_{z}|)$ with $0<f_{\rm r}(|v_{z}|)\leq 1$. Equation (12) should
be generalized to

$\displaystyle R_{\rm i}(r,\bm{v},t)=-4\pi\sigma^{2}\bm{n}\cdot\bm{v}f_{\rm
r}(|v_{z}|)\theta(-\bm{n}\cdot\bm{v})f(r,\bm{v},t).$ (13)

The boundary condition for collision-induced reaction without reflection can
be obtained when $f_{\rm r}(|v_{z}|)=1$. When $R_{\rm i}(r,\bm{v},t)$ is
independent of $\bm{v}$, Eq. (5) posed earlier is applied; this equation has
been studied explicitly for long-range reactions. Doi (1975)

##  III Cattaneo–Vernotte model

Here we shall introduce the Cattaneo–Vernotte model by defining
$\rho(\bm{r},t)$ and $\bm{j}(\bm{r},t)$ by the zeroth and first moments of the
infinite moment equations for $f(\bm{r},\bm{v},t)$ Davies (1954); Wilemski
(1976)

$\displaystyle\rho(\bm{r},t)$ $\displaystyle=\int d\bm{v}f(\bm{r},\bm{v},t),$
(14) $\displaystyle\bm{j}(\bm{r},t)$ $\displaystyle=\int
d\bm{v}\bm{v}f(\bm{r},\bm{v},t),$ (15)

respectively. Although only the first two moments are considered, our method
belongs to the Grad procedure, where the orthogonality relation of Hermite
functions is used as a moment closure method to an infinite hierarchy of
kinetic equations for the Hermite moments; Meyer and Schröter (1983); Grad
(1949) our results can be systematically extended to higher orders. Besides,
the other systematic method may also be applicable here. Traytak (2014) Using
Eq. (5), we find

$\displaystyle\frac{\partial}{\partial
t}\rho(\bm{r},t)+\frac{\partial}{\partial\bm{r}}\cdot\bm{j}(\bm{r},t)=-\int
d\bm{v}R(\bm{r},\bm{v},t).$ (16)

By multiplying $\bm{v}$ to Eq. (5) and integrating over $\bm{v}$, we obtain

$\displaystyle\frac{\partial}{\partial t}\bm{j}(\bm{r},t)$
$\displaystyle=-\frac{1}{\tau_{\rm D}}\left[\bm{j}(\bm{r},t)+\tau_{\rm
D}\left(\int
d\bm{v}\bm{v}\frac{\partial}{\partial\bm{r}}\cdot\bm{v}f(\bm{r},\bm{v},t)+\frac{1}{\mu}\frac{\partial
U}{\partial\bm{r}}\rho(\bm{r},t)+\int
d\bm{v}\bm{v}R(\bm{r},\bm{v},t)\right)\right].$ (17)

For isotropic systems, the reflecting boundary condition can be written in
terms of $\bm{j}(\bm{r},t)$ as

$\displaystyle\bm{n}\cdot\bm{j}(\sigma-\epsilon,t)$ $\displaystyle=0,$ (18)

from Eq. (10).

We define the diffusion coefficient [$D_{ij}(t)$] using, Davies (1954)

$\displaystyle\rho(\bm{r},t)D_{ij}(t)=\tau_{\rm D}\int
d\bm{v}\bm{v}_{i}\bm{v}_{j}f(\bm{r},\bm{v},t),$ (19)

where $\bm{v}_{i}=\bm{v}\cdot\bm{r}_{i}/|\bm{r}_{i}|$ and $\bm{v}_{j}$ is
defined in the similar manner. Equation (17) can be rewritten as

$\displaystyle\frac{\partial}{\partial t}\bm{j}_{i}(\bm{r},t)$
$\displaystyle=-\frac{1}{\tau_{\rm
D}}\left[\bm{j}_{i}(\bm{r},t)+\sum_{j}\frac{\partial}{\partial\bm{r}_{j}}D_{ij}(t)\rho(\bm{r},t)+\frac{\tau_{\rm
D}}{\mu}\frac{\partial U}{\partial\bm{r}_{i}}\rho(\bm{r},t)+\tau_{\rm D}\int
d\bm{v}\bm{v}_{i}R(\bm{r},\bm{v},t)\right].$ (20)

When the velocity distribution is in equilibrium at the initial time, we
introduce a decoupling approximation,

$\displaystyle f(\bm{r},\bm{v},t)=\rho(\bm{r},t)g_{\rm eq}(\bm{v}),$ (21)

where $g_{\rm eq}(\bm{v})$ indicates the Maxwell (Gaussian) distribution of
velocity. Ignoring the hydrodynamic interactions effects, we can express the
translational diffusion tensor [Eq. (19)] using the Einstein relation as

$\displaystyle D_{ij}=\tau_{\rm D}\int d\bm{v}v_{i}v_{j}g_{\rm
eq}(\bm{v})=D\delta_{ij},$ (22)

where

$D=\tau_{\rm D}\frac{k_{\rm B}T}{\mu}\,.$

Thus, Eq. (20) can be simplified as

$\displaystyle\frac{\partial}{\partial t}\bm{j}(\bm{r},t)$
$\displaystyle=-\frac{1}{\tau_{\rm
D}}\left[\bm{j}(\bm{r},t)+D\left(\frac{\partial}{\partial\bm{r}}\rho(\bm{r},t)+\frac{1}{k_{\rm
B}T}\frac{\partial U}{\partial\bm{r}}\rho(\bm{r},t)\right)+\tau_{\rm D}\int
d\bm{v}\bm{v}R(\bm{r},\bm{v},t)\right].$ (23)

Equations (16) and (23) constitute the Cattaneo–Vernotte differential model
extended to include a reaction.

For collision-induced reactions in isotropic systems, we obtain

$\displaystyle\frac{\partial}{\partial t}\rho(r,t)$
$\displaystyle+\frac{1}{r^{2}}\frac{\partial}{\partial
r}r^{2}j_{r}(r,t)=-\kappa\rho(r,t)\frac{\delta(r-\sigma)}{4\pi r^{2}},$ (24)
$\displaystyle\frac{\partial}{\partial t}j_{\rm r}(r,t)$
$\displaystyle=-\frac{1}{\tau_{\rm D}}\left[j_{\rm
r}(r,t)+D\left(\frac{\partial}{\partial r}\rho(r,t)+\frac{1}{k_{\rm
B}T}\frac{\partial U}{\partial r}\rho(r,t)\right)-\tau_{\rm D}\kappa_{\rm
r}\rho(r,t)\frac{\delta(r-\sigma)}{4\pi r^{2}}\right],$ (25)

where $\kappa$ is defined by

$\displaystyle\kappa=-4\pi\sigma^{2}\int d\bm{v}\bm{n}\cdot\bm{v}f_{\rm
r}(|v_{z}|)\theta(-\bm{n}\cdot\bm{v})g_{\rm eq}(\bm{v}),$ (26)

and $\kappa_{\rm r}$ satisfies $\vec{\bm{\kappa}}=\kappa_{\rm r}\bm{n}$, where
$\vec{\bm{\kappa}}$ is defined by

$\displaystyle\vec{\bm{\kappa}}=4\pi\sigma^{2}\int
d\bm{v}\bm{v}\bm{n}\cdot\bm{v}f_{\rm
r}(|v_{z}|)\theta(-\bm{n}\cdot\bm{v})g_{\rm eq}(\bm{v}).$ (27)

In Eq. (23), the term $\tau_{\rm D}\int d\bm{v}\bm{v}R(\bm{r},\bm{v},t)$
indicates the influence of a collision-induced reaction on the relaxation of
the current density. The current density relaxes toward the steady-state
quantity denoted by $\bm{j}_{\rm s}(\bm{r})$:

$\displaystyle\bm{j}_{\rm
s}(\bm{r})=-D\left(\frac{\partial}{\partial\bm{r}}\rho_{\rm
s}(\bm{r})+\frac{1}{k_{\rm B}T}\frac{\partial U}{\partial\bm{r}}\rho_{\rm
s}(\bm{r})\right)-\tau_{\rm D}\int d\bm{v}\bm{v}R_{\rm s}(\bm{r},\bm{v}),$
(28)

where the subscript denoted by “s” indicates the steady-state quantity. The
last term in Eq. (28) can be expressed as

$\displaystyle\tau_{\rm D}\int d\bm{v}\bm{v}R_{\rm
s}(\bm{r},\bm{v})=-\tau_{\rm D}\vec{\bm{\kappa}}\rho_{\rm
s}(r)\frac{\delta(r-\sigma)}{4\pi r^{2}},$ (29)

with $\vec{\bm{\kappa}}$ given by Eq. (27). $\vec{\bm{\kappa}}$ should be an
outward vector at $\sigma$ when $\bm{v}$ is the inward velocity induced by a
collision-induced reaction because $-\bm{n}\cdot\bm{v}$ should be negative to
have a non-zero value because of $\theta(-\bm{n}\cdot\bm{v})$. Later, in Eq.
(34), we show that $\kappa_{\rm r}>0$. The result $\kappa_{\rm r}>0$ indicates
a positive correlation of the velocity vector and the inward normal component
of the velocity vector at the contact distance. Equation (28) indicates that
the inward current density is reduced by the positive correlation at the
contact distance. In a steady state, the reduction of the inward current
density is given by the last term in Eq. (28), which shows the outward current
density. The sign of the last term in Eq. (25) is chosen accordingly.

When $f_{\rm r}(|v_{z}|)=f_{\rm r}$, where $f_{\rm r}$ is a constant
satisfying $0<f_{\rm r}\leq 1$, Eq. (26) is simplified to

$\displaystyle\kappa=-4\pi\sigma^{2}\int d\bm{v}\bm{n}\cdot\bm{v}f_{\rm
r}\theta(-\bm{n}\cdot\bm{v})g_{\rm eq}(\bm{v}).$ (30)

We choose the $z$-axis in the Cartesian coordinate for expressing $\bm{v}$
[$\bm{v}=(v_{x},v_{y},v_{z})$] in the direction of $\bm{n}$. In this
coordinate, we have $\bm{n}\cdot\bm{v}=v_{z}$ and
$\int_{-\infty}^{\infty}dv_{z}\theta(-v_{z})\cdots=\int_{-\infty}^{0}dv_{z}\cdots$.
Using the explicit expression

$\displaystyle g_{\rm eq}(\bm{v})=\left(\frac{\mu}{2\pi k_{\rm
B}T}\right)^{3/2}\exp\left(-\frac{\mu v^{2}}{2k_{\rm B}T}\right)$ (31)

with $v^{2}=v_{x}^{2}+v_{y}^{2}+v_{z}^{2}$, we find using Eq. (22) Kapral
(1978); Naqvi, Waldenstrøm, and Mork (1983); Molski (1988); Zhou and Szabo
(1991); Kim  _et al._ (2009) that

$\displaystyle\kappa$ $\displaystyle=2\sigma^{2}\sqrt{\frac{2\pi k_{\rm
B}T}{\mu}}f_{\rm r}=2\sigma^{2}\sqrt{\frac{2\pi D}{\tau_{\rm D}}}\,f_{\rm r}.$
(32)

We also calculate $\vec{\bm{\kappa}}$:

$\displaystyle\vec{\bm{\kappa}}=4\pi\sigma^{2}f_{\rm r}\int
d\bm{v}\bm{v}\left(\bm{n}\cdot\bm{v}\right)\theta(-\bm{n}\cdot\bm{v})g_{\rm
eq}(\bm{v}).$ (33)

When evaluating $\vec{\bm{\kappa}}$, we choose the $z$-axis in the Cartesian
coordinates for expressing $\bm{v}$ in the direction of $\bm{n}$ and find
$\vec{\bm{\kappa}}=\kappa_{\rm r}\bm{n}$ with

$\displaystyle\kappa_{\rm r}=\frac{2\pi\sigma^{2}k_{\rm B}T}{\mu}f_{\rm
r}=2\pi\sigma^{2}\frac{D}{\tau_{\rm D}}f_{\rm r},$ (34)

where Eq. (22) is used. Equation (33) can be interpreted as the correlation
between the velocity vector and the inward normal component of the velocity
vector at the contact distance. Equation (34) with Eq. (25) indicates that the
correlation at the contact distance is positive and reduces the inward current
density.

The perfectly absorbing boundary condition for a collision-induced reaction
implies that all inward fluxes lead to a reaction and that the reflected
outward fluxes are zero at the reactive boundary. Harris (1981); Burschka and
Titulaer (1981) In this case, a half-Gaussian distribution, where the velocity
component in the outward direction at the reactive boundary is set to zero,
has been introduced. In general, the velocity vector projected to $\bm{n}$ in
the negative direction is absorbed and the Gaussian distribution can be
distorted. Determining the proper distribution requires another approach. For
the moment, we introduce the Gaussian distribution in Eq. (21). In the
subsequent section, we apply the projection operator method to derive the
corresponding terms perturbatively.

We denote the Laplace transform of $\rho(r,t)$ and $j_{\rm r}(r,t)$ with
respect to $t$ by $\hat{\rho}(r,s)$ and $\hat{j}_{\rm r}(r,s)$, respectively.
By the Laplace transform of Eq. (25), we obtain,

$\displaystyle\hat{j}_{\rm r}(r,s)-j_{\rm r}(r,0)$
$\displaystyle=-\left[D(s)\left(\frac{\partial}{\partial
r}\hat{\rho}(r,s)+\frac{1}{k_{\rm B}T}\frac{\partial U}{\partial
r}\hat{\rho}(r,s)\right)-\tau_{\rm D}(s)\kappa_{\rm
r}\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi r^{2}}\right].$ (35)

By substituting Eq. (35) into the Laplace transform of Eq. (24), we find

$\displaystyle
s\hat{\rho}(r,s)-\rho(r,0)+\frac{1}{r^{2}}\frac{\partial}{\partial
r}r^{2}j_{\rm r}(r,0)$ $\displaystyle=\frac{1}{r^{2}}\frac{\partial}{\partial
r}D(s)r^{2}\frac{\partial}{\partial r}\hat{\rho}(r,s)$
$\displaystyle-\kappa\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi
r^{2}}-\tau_{\rm D}(s)\frac{1}{r^{2}}\frac{\partial}{\partial
r}r^{2}\kappa_{\rm r}\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi r^{2}},$ (36)

where $D(s)=D/(1+s\tau_{\rm D})$ and $\tau_{\rm D}(s)=\tau_{\rm
D}/(1+s\tau_{\rm D})$ are introduced and the potential function is set equal
to zero for simplicity.

Equations (24) and (25) can be expressed as the modified telegraph equation
given by a function of $\rho(r,t)$ alone:

$\displaystyle\tau_{\rm D}$ $\displaystyle\frac{\partial^{2}}{\partial
t^{2}}\rho(r,t)+\frac{\partial}{\partial
t}\rho(r,t)=\frac{1}{r^{2}}\frac{\partial}{\partial
r}Dr^{2}\left(\frac{\partial}{\partial r}\rho(r,t)+\frac{1}{k_{\rm
B}T}\frac{\partial U}{\partial r}\rho(r,t)\right)$
$\displaystyle-\left(\kappa+\tau_{\rm D}\frac{\partial}{\partial
t}\kappa\right)\rho(r,t)\frac{\delta(r-\sigma)}{4\pi r^{2}}-\tau_{\rm
D}\frac{1}{r^{2}}\frac{\partial}{\partial r}r^{2}\kappa_{\rm
r}\rho(r,t)\frac{\delta(r-\sigma)}{4\pi r^{2}},$ (37)

where the boundary condition must be given by Naqvi, Mork, and Waldenstrøm
(1982); Naqvi, Waldenstrøm, and Mork (1983)

$\displaystyle\left.\left(\frac{\partial}{\partial r}\rho(r,t)+\frac{1}{k_{\rm
B}T}\frac{\partial U}{\partial r}\rho(r,t)\right)\right|_{r=\sigma-\epsilon}$
$\displaystyle=0.$ (38)

In the Laplace domain, Eq. (37) can be expressed as

$\displaystyle s\hat{\rho}(r,s)-\rho(r,0)-$ $\displaystyle\tau_{\rm
D}(s)\left[\left.\frac{\partial}{\partial
t}\rho(r,t)\right|_{t=0}+\kappa\rho(r,0)\frac{\delta(r-\sigma)}{4\pi
r^{2}}\right]=\frac{1}{r^{2}}\frac{\partial}{\partial
r}D(s)r^{2}\frac{\partial}{\partial r}\hat{\rho}(r,s)$
$\displaystyle-\kappa\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi
r^{2}}-\tau_{\rm D}(s)\frac{1}{r^{2}}\frac{\partial}{\partial
r}r^{2}\kappa_{\rm r}\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi r^{2}},$ (39)

where $U=0$ is assumed and we take into account that the right-hand side of
Eq. (37) contains the term

$\displaystyle-\tau_{\rm D}\frac{\partial}{\partial
t}\kappa\rho(r,t)\frac{\delta(r-\sigma)}{4\pi r^{2}}.$

We consider the case of the homogeneous initial conditions

$\displaystyle\rho(r,0)=H(r-\sigma)\,,\quad\left.\frac{\partial}{\partial
t}\rho(r,t)\right|_{t=0}=0\,,\quad j_{\rm r}(r,0)=0\,,$ (40)

where $H(x)=1$ for $x>0$ otherwise zero; [$H(x)=1-\theta(-x)$]. When the
initial conditions are given by Eq. (40), both Eqs. (36) and (39) reduce to,

$\displaystyle s\hat{\rho}(r,s)-\rho(r,0)$
$\displaystyle=\frac{1}{r^{2}}\frac{\partial}{\partial
r}D(s)r^{2}\frac{\partial}{\partial r}\hat{\rho}(r,s)$
$\displaystyle-\kappa\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi
r^{2}}-\tau_{\rm D}(s)\frac{1}{r^{2}}\frac{\partial}{\partial
r}r^{2}\kappa_{\rm r}\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi r^{2}},$ (41)

which will be also derived using the projection operator method in the
subsequent section. Therefore, Eq. (41) can be derived either from Eqs. (24)
and (25) (the Cattaneo–Vernotte differential model) with the initial condition
given by $j_{\rm r}(r,0)=0$ or from Eq. (37) (the modified telegraph equation)
for the initial condition given by Eq. (40). It should be reminded that the
first initial condition of Eq. (40) is not required to derive Eq. (41) from
Eqs. (24) and (25) (the Cattaneo–Vernotte differential model); the initial
condition for $\rho(r,0)$ can be chosen arbitrary. Moreover, the first initial
condition of Eq. (40) should be given by $H(r-\sigma)$ rather than
$\theta(r-\sigma)$ because of the 4th term in Eq. (39). In these senses, it
might be preferable to use the Cattaneo–Vernotte differential model rather
than the modified telegraph equation as already concluded previously. Lee and
Traytak (2023)

By multiplying $\int_{\sigma-\epsilon}^{\infty}4\pi r^{2}dr$ on both sides of
Eq. (41) and introducing $\hat{p}(s)=\int_{\sigma-\epsilon}^{\infty}4\pi
r^{2}dr\hat{\rho}(r,s)$ with $p(0)=1$, we obtain

$\displaystyle s\hat{p}(s)-1$ $\displaystyle=-\kappa\hat{\rho}(\sigma,s),$
(42)

where the partial integration was performed to evaluate the last term
associated with $\kappa_{\rm r}$ in Eq. (41). By applying the inverse Laplace
transform, we obtain from Eq. (42),

$\displaystyle\frac{d}{dt}p(t)=-\kappa\rho(\sigma,t),$ (43)

and the rate coefficient can be given by (see Appendix B),

$\displaystyle k(t)=\kappa\rho(\sigma,t),$ (44)

which can be expressed in the Laplace domain as
$\hat{k}(s)=\kappa\hat{\rho}(\sigma,s)$. In principle, the reaction rate
coefficient should be obtained by taking the limit of $\epsilon\rightarrow 0$
after the inverse Laplace transform of $\kappa\hat{\rho}(\sigma,s)$, where the
reflecting boundary condition is set at $\sigma-\epsilon$. Here, we assume
that the limits are exchangeable.

It is instructive to express Eq. (41) in the time domain as

$\displaystyle\frac{\partial}{\partial t}\rho(r,t)$
$\displaystyle=\int_{0}^{t}dt_{1}\exp\left(-\frac{t-t_{1}}{\tau_{\rm
D}}\right)\left[\frac{1}{r^{2}}\frac{\partial}{\partial r}\frac{D}{\tau_{\rm
D}}r^{2}\frac{\partial}{\partial
r}\rho(r,t_{1})-\frac{1}{r^{2}}\frac{\partial}{\partial r}r^{2}\kappa_{\rm
r}\rho(r,t_{1})\frac{\delta(r-\sigma)}{4\pi r^{2}}\right]$
$\displaystyle-\kappa\rho(r,t)\frac{\delta(r-\sigma)}{4\pi r^{2}}.$ (45)

The term multiplied by $\kappa_{\rm r}$ is coupled to the memory kernel,
whereas the term multiplied by $\kappa$ is not coupled to the memory kernel.
The term multiplied by $\kappa_{\rm r}$ indicates a reduction of the current
density by the positive correlation between the velocity vector and the inward
normal component of the velocity vector at the contact distance, as explained
below Eq. (29). The memory kernel multiplied to $\kappa_{\rm r}$ represents
the competition of the current density associated with collision-induced
reaction and the diffusive flux during momentum relaxation (i.e., an inertial
effect).

##  IV Projection operator method

Here, we apply the projection operator method to Eq. (6) to obtain a closed
equation for $\rho(\bm{r},t)$. In Eq. (6), an explicit expression of
$R(\bm{r},\bm{v},t)$ is substituted in ${\cal L}_{1}$; using Eq. (13), we
first consider ${\cal L}_{1}={\cal L}_{\rm L}+{\cal L}_{\rm R}$, where ${\cal
L}_{\rm R}$ is given by,

$\displaystyle{\cal L}_{\rm R}$ $\displaystyle=\bm{n}\cdot\bm{v}f_{\rm
r}(|v_{z}|)\theta(-\bm{n}\cdot\bm{v})\delta(r-\sigma).$ (46)

We define the projection operator for any function ${\cal O}(\bm{r},\bm{v})$
as Kubo, Toda, and Hashitsume (1991); Northrup and Hynes (1978); Seki,
Barzykin, and Tachiya (1999); Bandyopadhyay, Seki, and Tachiya (2000)

$\displaystyle{\cal P}{\cal O}(\bm{r},\bm{v})=g_{\rm
eq}(\bm{v})\int\mbox{d}\bm{v}{\cal O}(\bm{r},\bm{v}),$ (47)

and the complementary projection operator by

$\displaystyle{\cal Q}=1-{\cal P}.$ (48)

We note that ${\cal L}_{0}g_{\rm eq}(\bm{v})=0$ and

$\displaystyle{\cal P}{\cal L}_{0}\cdots$ $\displaystyle=$ $\displaystyle 0,$
(49)

where we used the fact that $g_{\rm eq}(\bm{v})$ exponentially goes to zero as
$|\bm{v}_{i}|\rightarrow\infty$. Following the conventional projection
operator formalism using ${\cal L}={\cal L}_{0}+{\cal L}_{1}$, we first
introduce Kubo, Toda, and Hashitsume (1991)

$\displaystyle\frac{\partial}{\partial t}{\cal{P}}f(\bm{r},\bm{v},t)$
$\displaystyle=$ $\displaystyle{\cal P}({\cal L}{\cal P}f+{\cal L}{\cal Q}f),$
(50) $\displaystyle\frac{\partial}{\partial t}{\cal Q}f(\bm{r},\bm{v},t)$
$\displaystyle=$ $\displaystyle{\cal Q}({\cal L}{\cal Q}f+{\cal L}{\cal P}f).$
(51)

The solution of Eq. (51) can be formally expressed as

$\displaystyle{\cal
Q}f(\bm{r},\bm{v},t)=\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal Q}{\cal
L}\right]{\cal Q}{\cal L}{\cal P}f(\bm{r},\bm{v},t_{1})+\exp\left(t{\cal
Q}{\cal L}\right){\cal Q}f_{0},$ (52)

where $f_{0}=f(\bm{r},\bm{v},0)$. By substituting Eq. (52) into Eq. (50), we
can express the formally closed equation ${\cal{P}}f(\bm{r},\bm{v},t)$ as

$\displaystyle\frac{\partial}{\partial t}{\cal{P}}f(\bm{r},\bm{v},t)=$
$\displaystyle{\cal P}{\cal L}{\cal P}f(\bm{r},\bm{v},t)+{\cal P}{\cal
L}\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal Q}{\cal L}\right]{\cal Q}{\cal
L}{\cal P}f(\bm{r},\bm{v},t_{1})$ $\displaystyle+{\cal P}{\cal
L}\exp\left(t{\cal Q}{\cal L}\right){\cal Q}f_{0}.$ (53)

We consider the case of an initial equilibrium velocity distribution
$f_{0}=\rho(\bm{r},0)g_{\rm eq}(\bm{v})$ and find ${\cal Q}f_{0}=0$ because
${\cal P}g_{\rm eq}(\bm{v})=g_{\rm eq}(\bm{v})$. By noting that ${\cal
L}_{0}{\cal P}=0$ using ${\cal L}_{0}g_{\rm eq}(\bm{v})=0$ and Eq. (49), we
can simplify Eq. (53) to

$\displaystyle\frac{\partial}{\partial t}{\cal{P}}f(\bm{r},\bm{v},t)$
$\displaystyle={\cal P}{\cal L}_{1}{\cal P}f(\bm{r},\bm{v},t)+{\cal P}{\cal
L}_{1}\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal Q}{\cal L}\right]{\cal
Q}{\cal L}_{1}{\cal P}f(\bm{r},\bm{v},t_{1}).$ (54)

Because $g_{\rm eq}(\bm{v})$ exponentially goes to zero as
$|\bm{v}_{i}|\rightarrow\infty$, we have

$\displaystyle{\cal P}\frac{\partial}{\partial\bm{v}}\cdots{\cal P}f=0;$ (55)

we also have ${\cal{P}}\bm{v}\cdots{\cal P}f=0$. Equation (54) is further
simplified to

$\displaystyle\frac{\partial}{\partial t}{\cal{P}}f(\bm{r},\bm{v},t)$
$\displaystyle={\cal P}{\cal L}_{\rm R}{\cal P}f+{\cal
P}\left(-\bm{v}\cdot\frac{\partial}{\partial\bm{r}}+{\cal L}_{\rm
R}\right)\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal Q}{\cal L}\right]{\cal
Q}{\cal L}_{1}{\cal P}f(\bm{r},\bm{v},t_{1}).$ (56)

By applying perturbation expansion for the reaction sink term and the
streaming term, we obtain the lowest-order approximation by changing
$(t-t_{1}){\cal Q}{\cal L}$ to $(t-t_{1}){\cal L}_{0}$, where Eq. (49) is
used; we also obtain

$\displaystyle\frac{\partial}{\partial t}\rho(\bm{r},t)$
$\displaystyle\approx\int\mbox{d}\bm{v}{\cal L}_{\rm R}g_{\rm
eq}\rho+\int\mbox{d}\bm{v}\left(-\bm{v}\cdot\frac{\partial}{\partial\bm{r}}+{\cal
L}_{\rm R}\right)\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal
Q}{\cal L}_{1}g_{\rm eq}\rho(\bm{r},t_{1}),$ (57)

where ${\cal P}f=g_{\rm eq}\rho$. Using Eq. (26), we can express
$\int\mbox{d}\bm{v}{\cal L}_{\rm R}g_{\rm
eq}=-\kappa/(4\pi\sigma^{2})\delta(r-\sigma)$ and

$\displaystyle\frac{\partial}{\partial t}\rho(\bm{r},t)=$
$\displaystyle-\kappa\rho(\bm{r},t)\frac{\delta(r-\sigma)}{4\pi\sigma^{2}}+$
$\displaystyle\int\mbox{d}\bm{v}\left(-\bm{v}\cdot\frac{\partial}{\partial\bm{r}}+{\cal
L}_{\rm R}\right)\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal
Q}{\cal L}_{1}g_{\rm eq}\rho(\bm{r},t_{1}).$ (58)

We here further study the last term in Eq. (58). By introducing

$\displaystyle\frac{\partial}{\partial{\bm{v}}}g_{\rm
eq}(\bm{v})=-\frac{\mu{\bm{v}}}{k_{\rm B}T}g_{\rm eq}(\bm{v}),$ (59)

we obtain

$\displaystyle{\cal L}_{\rm L}g_{\rm
eq}\rho(\bm{r},t_{1})=-{\bm{v}}\cdot\left(\frac{\partial}{\partial\bm{r}}+\frac{1}{k_{\rm
B}T}\frac{\partial U}{\partial\bm{r}}\right)g_{\rm eq}\rho(\bm{r},t_{1}).$
(60)

Because ${\bm{v}}g_{\rm eq}$ is the eigen vector of the operator ${\cal
L}_{0}$,

$\displaystyle{\cal L}_{0}{\bm{v}}g_{\rm eq}=-\frac{1}{\tau_{\rm
D}}{\bm{v}}g_{\rm eq},$ (61)

and ${\cal P}{\bm{v}}g_{\rm eq}=0$; we therefore obtain

$\displaystyle\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal Q}{\cal L}_{\rm
L}g_{\rm eq}\rho(\bm{r},t_{1})=\exp\left[-(t-t_{1})/\tau_{\rm D}\right]{\cal
L}_{\rm L}g_{\rm eq}\rho(\bm{r},t_{1}).$ (62)

Part of Eq. (57) can be expressed as

$\displaystyle-\int\mbox{d}\bm{v}\bm{v}\cdot\frac{\partial}{\partial\bm{r}}\exp\left[(t-t_{1}){\cal
L}_{0}\right]{\cal Q}{\cal L}_{\rm L}g_{\rm
eq}\rho(\bm{r},t_{1})=-\int\mbox{d}\bm{v}\bm{v}\cdot\frac{\partial}{\partial\bm{r}}\exp\left[-(t-t_{1})/\tau_{\rm
D}\right]{\cal L}_{\rm L}g_{\rm eq}\rho(\bm{r},t_{1}).$ (63)

Using Eq. (22), we obtain

$\displaystyle-\int\mbox{d}\bm{v}\bm{v}\cdot\frac{\partial}{\partial\bm{r}}$
$\displaystyle\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal
Q}{\cal L}_{\rm L}g_{\rm eq}\rho(\bm{r},t_{1})=$
$\displaystyle\int_{0}^{t}dt_{1}\,\exp\left[-(t-t_{1})/\tau_{\rm
D}\right]\frac{\partial}{\partial\bm{r}}\cdot\left[\frac{D}{\tau_{\rm
D}}\left(\frac{\partial}{\partial\bm{r}}+\frac{1}{k_{\rm B}T}\frac{\partial
U}{\partial\bm{r}}\right)\right]\rho(\bm{r},t_{1}).$ (64)

Similarly, we consider

$\displaystyle-\int\mbox{d}\bm{v}\bm{v}\cdot\frac{\partial}{\partial\bm{r}}$
$\displaystyle\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal
Q}{\cal L}_{\rm R}g_{\rm eq}\rho(\bm{r},t_{1}).$ (65)

By introducing Eq. (59), we have

$\displaystyle{\cal L}_{\rm R}g_{\rm eq}\rho(\bm{r},t_{1})$
$\displaystyle=-\bm{n}\cdot\bm{v}f_{\rm
r}(|v_{z}|)\theta(-\bm{n}\cdot\bm{v})\delta(r-\sigma)g_{\rm
eq}\rho(\bm{r},t_{1}).$ (66)

Using Eq. (61) and ${\cal P}{\bm{v}}g_{\rm eq}=0$, we obtain

$\displaystyle\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal
Q}{\cal L}_{\rm R}g_{\rm
eq}\rho(\bm{r},t_{1})=\int_{0}^{t}dt_{1}\exp\left[-(t-t_{1})/\tau_{\rm
D}\right]{\cal L}_{\rm R}g_{\rm eq}\rho(\bm{r},t_{1}).$ (67)

By introducing Eq. (27), we find

$\displaystyle\int d\bm{v}\bm{v}{\cal L}_{\rm R}g_{\rm
eq}\rho(\bm{r},t_{1})=\vec{\bm{\kappa}}\rho(\bm{r},t_{1})\frac{\delta(r-\sigma)}{4\pi\sigma^{2}}$
(68)

and obtain

$\displaystyle\int\mbox{d}\bm{v}\bm{v}\cdot\frac{\partial}{\partial\bm{r}}$
$\displaystyle\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal
Q}{\cal L}_{\rm R}g_{\rm eq}\rho(\bm{r},t_{1})=$
$\displaystyle\int_{0}^{t}dt_{1}\,\exp\left[-(t-t_{1})/\tau_{\rm
D}\right]\frac{\partial}{\partial\bm{r}}\cdot\vec{\bm{\kappa}}\rho(\bm{r},t_{1})\frac{\delta(r-\sigma)}{4\pi\sigma^{2}}.$
(69)

Finally, we consider

$\displaystyle\int\mbox{d}\bm{v}{\cal L}_{\rm R}\int_{0}^{t}dt_{1}$
$\displaystyle\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal Q}{\cal L}_{1}g_{\rm
eq}\rho(\bm{r},t_{1})\approx$ $\displaystyle\int\mbox{d}\bm{v}{\cal L}_{\rm
R}\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal Q}{\cal
L}_{\rm L}g_{\rm eq}\rho(\bm{r},t_{1}),$ (70)

in Eq. (57), where we substituted ${\cal L}_{\rm L}$ for ${\cal L}_{1}$ as the
lowest order in the perturbation expansion. By substituting Eq. (60) and using
Eq. (62), we obtain

$\displaystyle\int\mbox{d}\bm{v}$ $\displaystyle{\cal L}_{\rm
R}\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal Q}{\cal
L}_{\rm L}g_{\rm eq}\rho(\bm{r},t_{1})=$
$\displaystyle-\int\mbox{d}\bm{v}{\cal L}_{\rm
R}\int_{0}^{t}dt_{1}\exp\left[-(t-t_{1})/\tau_{\rm
D}\right]{\bm{v}}\cdot\left(\frac{\partial}{\partial\bm{r}}+\frac{1}{k_{\rm
B}T}\frac{\partial U}{\partial\bm{r}}\right)g_{\rm eq}\rho(\bm{r},t_{1}).$
(71)

Using Eq. (68), we can express Eq. (71) as

$\displaystyle\int\mbox{d}\bm{v}$ $\displaystyle{\cal L}_{\rm
R}\int_{0}^{t}dt_{1}\exp\left[(t-t_{1}){\cal L}_{0}\right]{\cal Q}{\cal
L}_{\rm L}g_{\rm eq}\rho(\bm{r},t_{1})=$
$\displaystyle-\int_{0}^{t}dt_{1}\exp\left[-(t-t_{1})/\tau_{\rm
D}\right]\frac{\delta(r-\sigma)}{4\pi\sigma^{2}}\vec{\bm{\kappa}}\cdot\left(\frac{\partial}{\partial\bm{r}}+\frac{1}{k_{\rm
B}T}\frac{\partial U}{\partial\bm{r}}\right)\rho(\bm{r},t_{1})$
$\displaystyle=0,\mbox{ for }\epsilon\rightarrow 0,$ (72)

where the reflecting boundary condition given by Eq. (38) is introduced in the
limit of $\epsilon\rightarrow 0$. In principle, the limit of
$\epsilon\rightarrow 0$ should be taken after obtaining the reaction rate
coefficient; however, this term can be shown to not contribute to the final
expression. Note that the derivative with respect to ${\bm{r}}$ is not applied
to the delta-function. If the derivative with respect to ${\bm{r}}$ is applied
to the delta-function, as in Eq. (69), we should not take this limit.

By collecting Eqs. (57), (64), (69), and (72), we obtain

$\displaystyle\frac{\partial}{\partial t}\rho(\bm{r},t)=$
$\displaystyle\int_{0}^{t}dt_{1}\exp\left(-\frac{t-t_{1}}{\tau_{\rm
D}}\right)\frac{\partial}{\partial\bm{r}}\cdot\left[\frac{D}{\tau_{\rm
D}}\left(\frac{\partial}{\partial\bm{r}}+\frac{1}{k_{\rm B}T}\frac{\partial
U}{\partial\bm{r}}\right)\right]\rho(\bm{r},t_{1})-\kappa\rho(\bm{r},t)\frac{\delta(r-\sigma)}{4\pi\sigma^{2}}$
$\displaystyle-\int_{0}^{t}dt_{1}\exp\left(-\frac{t-t_{1}}{\tau_{\rm
D}}\right)\frac{\partial}{\partial\bm{r}}\cdot\vec{\bm{\kappa}}\rho(\bm{r},t_{1})\frac{\delta(r-\sigma)}{4\pi\sigma^{2}},$
(73)

where $\kappa$ is given by Eq. (32) and $\vec{\bm{\kappa}}$ can be expressed
as $\vec{\bm{\kappa}}=\kappa_{\rm r}\bm{n}$ using $\kappa_{\rm r}$ given by
Eq. (34). In this way, Eq. (45) is reproduced using the projection operator
method. When the time convolution is decoupled, the time integration gives the
factor $\tau_{\rm D}$ multiplied to $D/\tau_{\rm D}$ and $\vec{\bm{\kappa}}$
in Eq. (73). In isotropic systems, we can express $\tau_{\rm D}\kappa_{\rm
r}=2\pi\sigma^{2}D$ using Eqs. (22) and (34).

Thus far, we have considered the case of collision-induced reactions, where
the intrinsic reaction rate constant depends on the velocity, as shown by the
multiplication factor of $f(\sigma,\bm{v},t)$ in Eq. (12). For electron
transfer and energy transfer, the intrinsic reaction rate constant can be
independent of the reactant velocity. When the intrinsic reaction rate
constant is independent of the reactant velocity and is localized, the
intrinsic reaction rate constant can be expressed as

$\displaystyle R(\bm{v},\bm{r},t)=\frac{\kappa_{\rm
i}}{4\pi\sigma^{2}}f(\sigma,\bm{v},t)\delta(r-\sigma).$ (74)

We can express ${\cal L}_{\rm R}$ as

$\displaystyle{\cal L}_{\rm R}$ $\displaystyle=-\kappa_{\rm
i}f(\sigma,\bm{v},t)\delta(r-\sigma)$ (75)

instead of Eq. (46) and obtain $\int\mbox{d}\bm{v}{\cal L}_{\rm R}g_{\rm
eq}=-[\kappa_{\rm i}/(4\pi\sigma^{2})]\delta(r-\sigma)$. We then have

$\displaystyle\int d\bm{v}\bm{v}{\cal L}_{\rm R}g_{\rm
eq}\rho(\bm{r},t_{1})=0,$ (76)

instead of Eq. (68) because $\vec{\bm{\kappa}}=0$ when the intrinsic reaction
rate constant is given by Eq. (74). In this case, we obtain

$\displaystyle\frac{\partial}{\partial t}\rho(\bm{r},t)$
$\displaystyle=\int_{0}^{t}dt_{1}\exp\left(-\frac{t-t_{1}}{\tau_{\rm
D}}\right)\frac{\partial}{\partial\bm{r}}\cdot\left[\frac{D}{\tau_{\rm
D}}\left(\frac{\partial}{\partial\bm{r}}+\frac{1}{k_{\rm B}T}\frac{\partial
U}{\partial\bm{r}}\right)\right]\rho(\bm{r},t_{1})-\kappa_{\rm
i}\rho(\bm{r},t)\frac{\delta(r-\sigma)}{4\pi\sigma^{2}}$ (77)

instead of Eq. (73).

##  V Analytical solution

We solve Eq. (73) in the Laplace domain for $U=0$ expressed as

$\displaystyle s\hat{\rho}(r,s)-\rho(r,0)$
$\displaystyle=\frac{1}{r^{2}}\frac{\partial}{\partial
r}D(s)r^{2}\frac{\partial}{\partial r}\hat{\rho}(r,s)$
$\displaystyle-\kappa\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi
r^{2}}-\tau_{\rm D}(s)\frac{1}{r^{2}}\frac{\partial}{\partial
r}r^{2}\kappa_{\rm r}\hat{\rho}(r,s)\frac{\delta(r-\sigma)}{4\pi r^{2}},$ (78)

where $D(s)=D/(1+s\tau_{\rm D})$ and $\tau_{\rm D}(s)=\tau_{\rm
D}/(1+s\tau_{\rm D})$. The long-time limit of the rate coefficient as well as
the initial rate coefficient will be derived from the analytical expression of
the rate coefficient in the Laplace domain.

The Green’s function for

$\displaystyle s\hat{\rho}(r,s)-\frac{\delta(r-r_{\rm i})}{4\pi r^{2}}$
$\displaystyle=\frac{1}{r^{2}}\frac{\partial}{\partial
r}r^{2}D(s)\frac{\partial}{\partial r}\hat{\rho}(r,s),$ (79)

is given by

$\displaystyle\hat{G}_{\epsilon}\left(r,r_{\rm i},s\right)=\frac{1}{8\pi
D(s)\zeta(s)rr_{\rm i}}\left[\exp\left[-|r-r_{\rm
i}|\zeta(s)\right]-\frac{1-\sigma_{\rm i}\zeta(s)}{1+\sigma_{\rm
i}\zeta(s)}\exp\left[-(r+r_{\rm i}-2\sigma_{\rm i})\zeta(s)\right]\right],$
(80)

where $\zeta(s)=\sqrt{s/D(s)}$, and $\sigma_{\rm i}=\sigma-\epsilon$ is
introduced because of the reflecting boundary condition at
$r=\sigma-\epsilon$. The formal solution can be expressed using the Green
function as

$\displaystyle\hat{\rho}(r,s)=$
$\displaystyle\frac{1}{s}-\kappa\hat{G}_{\epsilon}(r,\sigma,s)\hat{\rho}(\sigma,s)-$
$\displaystyle 4\pi\tau_{\rm
D}(s)\int_{\sigma-\epsilon}^{\infty}dr_{0}\hat{G}_{\epsilon}(r,r_{0},s)\frac{\partial}{\partial
r_{0}}r_{0}^{2}\kappa_{\rm
r}\hat{\rho}(r_{0},s)\frac{\delta(r_{0}-\sigma)}{4\pi r_{0}^{2}},$ (81)

where we used

$\displaystyle\lim_{\epsilon\rightarrow 0}\int_{\sigma-\epsilon}^{\infty}4\pi
r_{\rm i}^{2}dr_{\rm i}\hat{G}_{\epsilon}(\sigma,r_{\rm i},s)\rho(r_{\rm
i},0)=1,$ (82)

for the initial uniform distribution with a reflecting boundary condition at
$r=\sigma-\epsilon$. By performing the partial integration, the last term in
Eq. (81) can be rewritten as

$\displaystyle 4\pi\tau_{\rm
D}(s)\int_{\sigma-\epsilon}^{\infty}dr_{0}\hat{G}_{\epsilon}(r,r_{0},s)$
$\displaystyle\frac{\partial}{\partial r_{0}}r_{0}^{2}\kappa_{\rm
r}\hat{\rho}(r_{0},s)\frac{\delta(r_{0}-\sigma)}{4\pi r_{0}^{2}}$
$\displaystyle=-\tau_{\rm D}(s)\left.\kappa_{\rm
r}\hat{\rho}(\sigma,s)\frac{\partial}{\partial
r_{0}}\hat{G}_{\epsilon}(r,r_{0},s)\right|_{r_{0}=\sigma}.$ (83)

By substituting Eq. (83) into Eq. (81), a closed equation for
$\hat{\rho}(\sigma,s)$ can be obtained,

$\displaystyle\hat{\rho}(\sigma,s)=$
$\displaystyle\frac{1}{s}-\kappa\hat{G}_{\epsilon}(\sigma,\sigma,s)\hat{\rho}(\sigma,s)+\tau_{\rm
D}(s)\left.\kappa_{\rm r}\hat{\rho}(\sigma,s)\frac{\partial}{\partial
r_{0}}\hat{G}_{\epsilon}(\sigma,r_{0},s)\right|_{r_{0}=\sigma}.$ (84)

There remains a subtlety in the boundary condition. Depending on the limit, we
have

$\displaystyle\lim_{\epsilon\rightarrow
0}\lim_{r\rightarrow\sigma}\frac{\partial}{\partial
r}\hat{G}_{\epsilon}\left(\sigma,r,s\right)$ $\displaystyle=0\mbox{ for
}\sigma>r,$ (85) $\displaystyle\lim_{\epsilon\rightarrow
0}\lim_{r\rightarrow\sigma}\frac{\partial}{\partial
r}\hat{G}_{\epsilon}\left(\sigma,r,s\right)$
$\displaystyle=-\frac{1}{4\pi\sigma^{2}D(s)}\mbox{ for }r>\sigma.$ (86)

$\hat{G}_{\epsilon}\left(\sigma,r,s\right)$ is not differentiable at the point
$r=\sigma$ where the first derivative of $|x|$ with respect to $x$ is used
when performing the derivative with respect to $r_{i}$ for Eq. (80) by
applying the chain rule for calculating the derivative of composition of $|x|$
and $x=r-r_{\rm i}$. We, therefore, consider Eq. (85) or (86). In the limit of
$\epsilon\rightarrow 0$, $r$ approaches $\sigma$ from the side satisfying
$r>\sigma$ rather than the side satisfying $r<\sigma$; by considering that the
physical meaningful limit towards the reaction surface is
$r\rightarrow\sigma+$, we use Eq. (86) rather than Eq. (85). Equation (86) is
also consistent with the previous theoretical model, where a division into
inner and outer spatial regions separated by a boundary is introduced.
Northrup and Hynes (1978) By adopting Eq. (86) and using Eq. (84), we obtain

$\displaystyle\hat{\rho}(\sigma,s)$ $\displaystyle=\lim_{\epsilon\rightarrow
0}\frac{1}{s}\left(1+\kappa\hat{G}_{\epsilon}(\sigma,\sigma,s)+\frac{\tau_{\rm
D}}{4\pi\sigma^{2}D}\kappa_{\rm r}\right)^{-1},$ (87)

Using Eq. (42), we obtain the Laplace transform of the reaction rate
coefficient, $k(t)$, as

$\displaystyle\hat{k}(s)=\lim_{\epsilon\rightarrow
0}\kappa\hat{\rho}(\sigma,s)=\frac{\kappa}{s}\left(1+\kappa\hat{G}_{0}(\sigma,\sigma,s)+\frac{\tau_{\rm
D}}{4\pi\sigma^{2}D}\kappa_{\rm r}\right)^{-1},$ (88)

where the limit of $\epsilon\rightarrow 0$ is taken and
$\hat{G}_{0}(\sigma,\sigma,s)$ can be expressed as

$\displaystyle\hat{G}_{0}(\sigma,\sigma,s)$ $\displaystyle=\frac{1+s\tau_{\rm
D}}{4\pi\sigma D}\frac{1}{1+\sigma\sqrt{s(1+s\tau_{\rm D})/D}}.$ (89)

Assuming that in Eq. (87) the limits as $r_{0},r\to\sigma$ and inverse Laplace
transformation are commutative, we can conclude: Equation (88) is the exact
reaction rate coefficient in the Laplace domain for the Cattaneo–Vernotte
differential model.

When $\kappa_{\rm r}=0$, Eq. (88) can be expressed as

$\displaystyle\hat{k}(s)=\lim_{\epsilon\rightarrow
0}\kappa\hat{\rho}(\sigma,s)=\frac{1}{s}\left(\frac{1}{\kappa}+\hat{G}_{0}(\sigma,\sigma,s)\right)^{-1}.$
(90)

The long-time reaction rate coefficient can be obtained as

$\displaystyle k_{\infty}$ $\displaystyle=\lim_{\epsilon\rightarrow
0}\lim_{s\rightarrow
0}\left(\frac{1}{\kappa}+\hat{G}_{\epsilon}(\sigma,\sigma,s)\right)^{-1}$ (91)
$\displaystyle=\left[\frac{1}{\kappa}\left(1+\frac{\tau_{\rm
D}}{\sigma}\frac{\kappa_{\rm r}}{4\pi\sigma D}\right)+\frac{1}{4\pi\sigma
D}\right]^{-1}$ (92) $\displaystyle=\left[\frac{1}{\kappa}\left(1+\frac{f_{\rm
r}}{2}\right)+\frac{1}{4\pi\sigma D}\right]^{-1},$ (93)

where Eq. (34) is substituted and $\kappa$ is given by Eq. (32). (The result
is the same as that obtained from $k_{\infty}=\lim_{s\rightarrow
0}s\hat{k}(s)$.) For collision-induced reactions without reflection at the
contact distance, we have $f_{\rm r}=1$. For $f_{\rm r}=1$, $\kappa$ is
effectively reduced to $2\kappa/3$ in Eq. (93). As explained below Eq. (29),
the inward current density is reduced by the positive correlation between the
velocity vector and the inward normal component of the velocity vector at the
contact distance. This effect reduces $\kappa$ to $2\kappa/3$ for collision-
induced reactions without reflection at the contact distance in the lowest
order of the perturbation expansion for the reaction sink term and the
streaming term. The persistence of inertial effects on the reduction of the
reaction rate coefficient has been also observed by Langevin dynamic
simulations. Yang, Han, and Lee (2001)

The aforementioned results can be compared to the long-time rate coefficient
obtained by assuming $\kappa_{\rm r}=0$ given by

$\displaystyle k_{\infty}^{(0)}$
$\displaystyle=\left[\frac{1}{\kappa}+\frac{1}{4\pi\sigma D}\right]^{-1}.$
(94)

In the opposite limit of $t\rightarrow 0$, $k_{\rm f}=k(0)$ can be estimated
from $k_{\rm f}=\lim_{s\rightarrow\infty}sk(s)$. We obtain

$\displaystyle k_{\rm f}=\left[\frac{1}{\kappa}\left(1+\frac{\tau_{\rm
D}}{\sigma}\frac{\kappa_{\rm r}}{4\pi\sigma
D}\right)+\frac{1}{4\pi\sigma^{2}}\sqrt{\frac{\tau_{\rm D}}{D}}\right]^{-1}$
(95)

using Eq. (88) and obtain

$\displaystyle k_{\rm f}^{(0)}$
$\displaystyle=\left[\frac{1}{\kappa}+\frac{1}{4\pi\sigma^{2}}\sqrt{\frac{\tau_{\rm
D}}{D}}\right]^{-1}$ (96)

by assuming $\kappa_{\rm r}=0$. The physical origin of this sudden drop of the
initial rate coefficient from the equilibrium value $\kappa$ was discussed in
Ref. Lee and Traytak, 2023.

##  VI Numerical results

The Fokker–Planck–Kramers equation is the continuous description in phase
space under Markovian momentum relaxation. Kramers (1940); Doi (1975);
Northrup and Hynes (1978); Naqvi, Mork, and Waldenstrøm (1982); Harris
(1983a); Molski (1988); Ibuki and Ueno (1997, 2003, 2006); Kim  _et al._
(2009) The assumption of Markovian momentum relaxation might be justified if
the mean free path is sufficiently smaller than the contact distance.
$\sqrt{k_{\rm B}T/(2\mu)}$ indicates the mean velocity of reactants in thermal
equilibrium, and $\ell_{v}=\sqrt{k_{\rm B}T/(2\mu)}\,\tau_{\rm D}$ can be
interpreted as the mean free path that reactants are able to travel during the
momentum relaxation time ($\tau_{\rm D}$) when a reaction is absent. The
condition that the mean free path is smaller than the contact distance can be
expressed by

$\displaystyle\tau_{\rm D}\sqrt{\frac{k_{\rm B}T}{2\mu}}<\sigma.$ (97)

If Eq. (97) is satisfied, momentum can be relaxed by successive collision-
induced events and the memory of successive collision-induced events can be
lost before the reactant moves at a distance characterized by $\sigma$. The
mean free path, $\ell_{v}$, divided by the contact distance corresponds to the
Knudsen number in hydrodynamics, which is required for a continuum description
of reactants. Under the condition $\ell_{v}/\sigma<1$, we take into account
the effect of a reaction in the Fokker–Planck–Kramers equation. In this
section, we present the results under $\ell_{v}/\sigma<1$.

For simulations, we rewrite Eqs. (24)–(25) with $U=0$. By multiplying both
sides of Eq. (24) by $4\pi\int_{\sigma-\epsilon}^{\sigma+\epsilon}drr^{2}$ and
taking the limit of $\epsilon\rightarrow 0$, we obtain

$\displaystyle\frac{\partial}{\partial t}\rho(r,t)$
$\displaystyle+\frac{1}{r^{2}}\frac{\partial}{\partial r}r^{2}j_{r}(r,t)=0,$
(98)

with the boundary condition given by

$\displaystyle j_{\rm r}(\sigma,t)$
$\displaystyle=-\frac{\kappa\rho(\sigma,t)}{4\pi\sigma^{2}},$ (99)

which is equivalent to Eq. (24) with the boundary condition given by $j_{\rm
r}(\sigma-\epsilon,t)=0$ if the limit of $\epsilon\rightarrow 0$ limit is
taken when calculating the reaction rate coefficient. The kinetic equation for
$j_{\rm r}(r,t)$ is given by Eq. (25) with $U=0$ as

$\displaystyle\frac{\partial}{\partial t}j_{\rm r}(r,t)$
$\displaystyle=-\frac{1}{\tau_{\rm D}}\left[j_{\rm
r}(r,t)+D\frac{\partial}{\partial r}\rho(r,t)-\tau_{\rm D}\kappa_{\rm
r}\rho(r,t)\frac{\delta(r-\sigma)}{4\pi r^{2}}\right].$ (100)

We here study the effect of the momentum relaxation time ($\tau_{\rm D}$) on
the long-time rate coefficient. For this purpose, we use $\tau_{\rm
c}=\sigma\sqrt{2\mu/(k_{\rm B}T)}$ as the time unit. We introduce
dimensionless radial coordinate $\zeta=r/\sigma$, time $\tau=t/\tau_{\rm c}$,
$J_{\rm r}=j_{\rm r}\tau_{\rm c}/\sigma$, $D_{\rm n}=D\tau_{\rm
c}/\sigma^{2}=2\tau_{\rm D}/\tau_{\rm c}$, $\kappa_{\rm n}=\kappa\tau_{\rm
c}/\sigma^{3}=4\sqrt{\pi}f_{\rm r}$, and $\kappa_{\rm rn}=\kappa_{\rm
r}\tau_{\rm D}\tau_{\rm c}/\sigma^{4}=2\pi(\tau_{\rm D}/\sigma)\sqrt{2k_{\rm
B}T/\mu}\,f_{\rm r}$, which is proportional to $\tau_{\rm D}$. By considering
that $D=\tau_{\rm D}k_{\rm B}T/\mu$ is also proportional to $\tau_{\rm D}$, we
express $\kappa_{\rm rn}=2\pi D_{\rm n}f_{\rm r}$ for numerical evaluation.

Figure 1: (Color online) Dimensionless rate coefficient
[$(k_{\infty}/\sigma^{2})\sqrt{2\mu/(k_{\rm B}T)}$] in the long-time limit
obtained using Eqs. (101)–(103) is shown against $D_{\rm
n}=(D/\sigma)\sqrt{2\mu/(k_{\rm B}T)}$ for
$(\kappa/\sigma^{2})\sqrt{2\mu/(k_{\rm B}T)}=4\sqrt{\pi}f_{\rm r}$, where
$\zeta_{\rm max}=100$. (a) and (b) show the results for $f_{\rm r}=1$ and
$f_{\rm r}=0.2$, respectively. Circles and crosses indicate the numerical
results obtained using the Lorentzian representation of the delta-function
$2\eta/[\pi(\zeta^{2}+\eta^{2})]$ with $\eta=10^{-6}$ and the Lorentzian
representation of the delta-function $\eta/[\pi(\zeta^{2}+\eta^{2})]$ located
at distance $5\eta$ shifted outward from the reflecting boundary,
respectively. The black short-dashed line is the diffusion-controlled limit of
$4\pi D_{\rm n}$. The black long-dashed line indicates the results for
$\kappa_{\rm r}=0$ given by Eq. (94). The red thick line indicates the results
obtained using Eq. (93).

Using Eqs. (24) and (100), we obtain the rate coefficient in the long-time
limit from $k_{\infty\rm{n}}=4\pi J_{\rm rs}(1)$ by numerically solving

$\displaystyle\frac{\partial}{\partial\zeta}\zeta^{2}J_{\rm rs}(\zeta)$
$\displaystyle=0,$ (101) $\displaystyle J_{\rm rs}(\zeta)+D_{\rm
n}\frac{\partial}{\partial\zeta}\rho_{\rm s}(\zeta)-\kappa_{\rm rn}\rho_{\rm
s}(\zeta)\frac{\delta(\zeta-1)}{4\pi}$ $\displaystyle=0,$ (102)

with the boundary condition given by

$\displaystyle J_{\rm rs}(1)=-\frac{\kappa_{\rm n}}{4\pi}\rho_{\rm s}(1),$
(103)

and $\rho_{\rm s}(\zeta_{\rm max})=1$; steady states are labeled with the
subscript “s”. We use the Lorentzian representation of the delta-function
given by $f_{\rm d1}(\zeta)=2\eta/[\pi(\zeta^{2}+\eta^{2})]$ with
$\eta=10^{-6}$, where the normalization is given by $\int_{0}^{\infty}dxf_{\rm
d1}(\zeta)=1$. We also use the Lorentzian representation of the delta-function
given by $f_{\rm d2}(\zeta)=\eta/[\pi(\zeta^{2}+\eta^{2})]$ by changing the
delta-function in Eq. (102) by $\delta(\zeta-1-5\eta)$ by noticing
$\int_{1}^{\infty}d\zeta f_{\rm d2}(\zeta)\approx 0.94$. The factor $5\eta$
with $\eta=10^{-6}$ corresponds to $\epsilon$ representing the difference
between the contact distance ($\sigma$), where the reaction sink term is
located, and $\sigma-\epsilon$, where the reflecting boundary is located;
here, the reflecting boundary is set at $\zeta=1$. As shown in Fig. 1, the
numerical results for $k_{\infty\rm{n}}$ obtained using Eqs. (101)–(103) and
$\zeta_{\rm max}=100$ are not influenced by the choice of Lorentzian
representation nor by changing $\eta$ to $\eta=10^{-4}$ (results not shown)
for $f_{\rm r}=1$ and $f_{\rm r}=0.2$. The results of Eq. (94), where
$\kappa_{\rm r}$ is ignored, can be regarded as the upper bound of the long-
time rate coefficient. The numerical results are close to the exact results
obtained from Eq. (93); they are consistent with each other.

By decreasing $D$, which is proportional to $\tau_{\rm D}$, we can approximate
the long-time rate coefficient by the diffusion-controlled rate. In Fig. 1,
the diffusion-controlled limit is shown by the black short-dashed line. The
diffusion-controlled rate coefficient is given by $4\pi\sigma D$, which is
proportional to the contact distance. In the reaction-controlled limit, the
long-time rate coefficient is limited by the rate constant expressed by
$\kappa=2\sigma^{2}\sqrt{2\pi k_{\rm B}T/\mu}\,f_{\rm r}$ [Eq. (32)], which is
proportional to the square of the contact distance because of the scattering
cross-sectional area for the ballistic motion of reactants. With an increase
in $D$, the long-time rate coefficient is more limited by the reaction.

Figure 2: (Color online) The dimensionless rate coefficient
[$(k_{\infty}/\sigma^{2})\sqrt{2\mu/(k_{\rm B}T)}$] in the long-time limit
obtained using Eqs. (101)-(103) for $\kappa_{\rm r}=0$ is shown against
$D_{\rm n}=(D/\sigma)\sqrt{2\mu/(k_{\rm B}T)}$ for
$(\kappa/\sigma^{2})\sqrt{2\mu/(k_{\rm B}T)}=4\sqrt{\pi}f_{\rm r}$, where
$\zeta_{\rm max}=100$. The upper red solid lines and circles and the lower
black solid lines and circles indicate $f_{\rm r}=1$ and $f_{\rm r}=0.2$,
respectively. Circles indicate the numerical solutions of Eqs. (101)–(103) for
$\kappa_{\rm r}=0$. The black short-dashed line is the diffusion-controlled
limit of $4\pi D_{\rm n}$. The thick lines indicate the results for
$\kappa_{\rm r}=0$ given by Eq. (94).

When the intrinsic reaction rate constant is independent of the reactant
velocity and is localized, we obtain $\kappa_{\rm r}=0$, as shown by Eq. (76).
In this case, the long-time rate constants are not influenced by the inertial
effect and can be obtained from Eq. (94), which is not influenced by
$\tau_{\rm D}$ except through the diffusion constant. In Fig. 2, the solid
lines are the same as the black long-dashed lines in Fig. 1. The numerical
solutions of Eqs. (101)–(103) using $f_{\rm d1}(\zeta)$ for $\kappa_{\rm r}=0$
are close to the values obtained from Eq. (94). Conversely, the influence of
the inertial effect remains in the values of the long-time rate constants, as
evidenced by the results in Fig. 1, when the intrinsic reaction rate constant
depends on the reactant velocity.

##  VII Conclusion

Although the Cattaneo–Vernotte model has been widely studied to take into
account momentum relaxation in transport equations, the effect of reactions on
the Cattaneo–Vernotte model has not yet been fully elucidated. How current
density associated with reactions can be expressed in the Cattaneo–Vernotte
model is unclear. We derived the effect of a reaction on the Cattaneo–Vernotte
model using the Fokker–Planck–Kramers equation.

We took into account momentum relaxation by applying the projection operator
method to the Fokker–Planck–Kramers equation with a reaction sink term
describing reactions. In the absence of a reaction sink term, a modified
Smoluchowski equation including a memory kernel could be derived. In the
lowest order of the perturbation expansion for the reaction sink term and the
streaming term, we obtain the modified Smoluchowski equation generalized to
include two reaction terms [Eq. (73) ] for collision-induced reactions, where
the intrinsic reaction rate constant depends on the relative velocity of
reactants. The term multiplied by $\kappa_{\rm r}$ is coupled to the memory
kernel, whereas the term multiplied by $\kappa$ is not coupled to the memory
kernel; the term multiplied by $\kappa_{\rm r}$ represents the competition of
the current density associated with a collision-induced reaction and the
diffusive flux during momentum relaxation. We showed that the current density
is reduced by the positive correlation between the velocity vector and the
inward normal component of the velocity vector at the contact distance.
Without $\kappa_{\rm r}$, the reaction rate coefficient is overestimated.
Equation (73) can be rewritten in the form of a generalized reaction–telegraph
equation [Eq. (37)], where a reaction–telegraph equation is generalized to
include the reduction effect of the current density by the positive
correlation of the current density at the contact distance. We also derived
the same equation as Eq. (73) by introducing decoupling between the
configurational distribution and the momentum distribution, where an
equilibrium Maxwell (Gaussian) distribution is assumed for the momentum.

The long-time rate coefficient turned out to be influenced by $\kappa_{\rm
r}$; the intrinsic reaction rate constant in the expression of the long-time
rate coefficient reduces to $2/3$ of the original value for a collision-
induced reaction without reflection at the contact distance in the lowest
order of the perturbation expansion for the reaction sink term and the
streaming term. The persistent inertial effects that reduce the reaction rate
coefficient are consistent with the results of Langevin dynamic simulations.
Yang, Han, and Lee (2001) Moreover, we examined the long-time rate coefficient
by changing the momentum relaxation time denoted by $\tau_{\rm D}$. The
diffusion constant is proportional to $\tau_{\rm D}$, whereas the reaction
rate constants such as $\kappa$ and $\kappa_{\rm r}$ are independent of
$\tau_{\rm D}$. When $\tau_{\rm D}$ is small, the long-time rate coefficient
can be approximated by the diffusion controlled rate, which is given by
$4\pi\sigma D$ and is proportional to the contact distance. When $\tau_{\rm
D}$ is increased, the long-time rate coefficient is influenced by the rate
constant given by $\kappa=2\sigma^{2}\sqrt{2\pi k_{\rm B}T/\mu}$ [Eq. (32)]
and is proportional to the square of the contact distance; in the reaction-
controlled limit, reactions proceed by ballistic collisions.

For electron transfer and energy transfer, the intrinsic reaction rate
constant in the reaction sink term of the Fokker–Planck–Kramers equation can
be independent of the reactant velocity. When the intrinsic reaction rate
constant is independent of the reactant velocity and is localized, we have
$\kappa_{\rm r}=0$; the aforementioned reduction of the long-time rate
coefficient is absent in this case.

For brevity, we ignored the liquid structure factor and the hydrodynamic
effect, which might influence reaction kinetics. Kapral (1978); Shin and
Kapral (1978); Zhou and Szabo (1991); Lee  _et al._ (2004); Northrup and Hynes
(1979) The hydrodynamic effect can be taken into account by the relative
distance dependence in the diffusion constant. The structure factor can be
taken into account by introducing the potential of mean force. Although we
considered a simplified model, we found that the inertial effect reduces the
long-time reaction rate coefficient for a collision-induced reaction.

###### Acknowledgements.

SL notes that this work was supported by the National Research Foundation of
Korea (NRF) grant funded by the Korea government (MSIT) (No.
2020R1F1A1071933). S.D.T. notes that this work was performed in accordance
with STATE TASK 45.12 (Grant No. 0082-2019-0009), registration number in GZ
122040500058-1. Scientific basis for the designing of new materials with
desired properties and functions, including highpurity and nanomaterials.
Topic 1.1 “Physics and chemistry of new nanostructured systems and composite
materials with prescribed properties.”

## Appendix A. Derivation of Eq. (10)

The Fokker–Planck–Kramers equation without the reaction sink term must be
conservative, and we infer

$\displaystyle\int_{r\geq\sigma-\epsilon}d\bm{r}\int
d\bm{v}\,\left[-\frac{\partial}{\partial\bm{r}}\cdot\bm{v}+\frac{\partial}{\partial\bm{v}}\frac{1}{\mu}\cdot\frac{\partial
U}{\partial\bm{r}}+\frac{\partial}{\partial\bm{v}}\frac{1}{\tau_{\rm
D}}\cdot\left(\bm{v}+\frac{k_{\rm
B}T}{\mu}\frac{\partial}{\partial\bm{v}}\right)\right]f(\bm{r},\bm{v},t)=0,$
(A.1)

where a reflecting boundary condition is set at $r=\sigma-\epsilon$ by
assuming isotropy. Provided function $f(\bm{r},\bm{v},t)$ vanishes
exponentially as $|\bm{v}|\rightarrow\infty$, we arrive at

$\displaystyle\int d\bm{v}\,\frac{\partial}{\partial\bm{v}}\frac{1}{\tau_{\rm
D}}\cdot\left(\bm{v}+\frac{k_{\rm
B}T}{\mu}\frac{\partial}{\partial\bm{v}}\right)f(\bm{r},\bm{v},t)=0,$ (A.2)
$\displaystyle\int
d\bm{v}\,\frac{\partial}{\partial\bm{v}}\frac{1}{\mu}\cdot\frac{\partial
U}{\partial\bm{r}}f(\bm{r},\bm{v},t)=0.$ (A.3)

We then must have

$\displaystyle\int_{r\geq\sigma-\epsilon}d\bm{r}\int
d\bm{v}\,\frac{\partial}{\partial\bm{r}}\cdot\bm{v}f(\bm{r},\bm{v},t)=0.$
(A.4)

We assume the initial equilibrium reactant distribution. Using the divergence
theorem and assuming that the influence of a reaction on the current density
at the infinity distance is negligible,

$\displaystyle-\lim_{r\rightarrow\infty}\int
d\bm{v}\,\bm{n}\cdot\bm{v}f(\bm{r},\bm{v},t)=0,$ (A.5)

we obtain

$\displaystyle\left.\int
d\bm{v}\,\bm{n}\cdot\bm{v}f(\bm{r},\bm{v},t)\right|_{r=\sigma-\epsilon}=0,$
(A.6)

for the isotropic distribution of reactants.

## Appendix B. Derivation of the reaction rate coefficient [Eq. (44)]

We consider bimolecular reaction between the reactants denoted by A and the
reactants denoted by B, where the volume number density of B denoted by
$c_{\rm B}$ is larger than the volume number density of A. The survival
probability of A is denoted by $S(t)$, which is the probability of finding A
survived from the bimolecular reaction with B at time $t$. We also introduce
the number of B reactants in the spherical volume ($V$) by $N_{\rm B}=c_{\rm
B}V$, where the volume is expressed as $V=4\pi[r_{\rm
o}^{3}-(\sigma-\epsilon]^{3})/3$ using $r_{\rm o}$ representing the distance
to an outer-spherical boundary from the center of A. Note that the pair
correlation function approaches $1$ as $r\rightarrow\infty$; $r_{\rm o}$
should be large enough to ensure $\int_{\sigma-\epsilon}^{r_{\rm o}}4\pi
r^{2}dr\rho(r,0)\approx V$. When A is surrounded by B reactants, the survival
probability of A decays by reaction to one of $N_{\rm B}$ B-reactants, which
can be expressed as

$\displaystyle S(t)=\left(\int_{\sigma-\epsilon}^{r_{\rm o}}4\pi
r^{2}dr\rho(r,t)/\int_{\sigma-\epsilon}^{r_{\rm o}}4\pi
r^{2}dr\rho(r,0)\right)^{N_{\rm B}},$ (A.7)

where the pair correlation function associated with each B reactant is assumed
to be independent and identical, therefore is given by using the same
expression, $\rho(r,t)$. From the definition, $S(0)=1$ is satisfied. By taking
the thermodynamic limit, $S(t)$ is related to $\rho(r,t)$ as Tachiya (1983)

$\displaystyle S(t)$ $\displaystyle=\lim_{r_{\rm
o}\rightarrow\infty}\left(1-\frac{1}{V}\int_{\sigma-\epsilon}^{r_{\rm o}}4\pi
r^{2}dr\left[\rho(r,0)-\rho(r,t)\right]\right)^{c_{\rm B}V}$ (A.8)
$\displaystyle=\exp\left[-c_{\rm B}\int_{\sigma-\epsilon}^{\infty}4\pi
r^{2}dr\left[\rho(r,0)-\rho(r,t)\right]\right],$ (A.9)

where $S(0)=1$ can be confirmed. Note that the original derivation is
formulated in terms of the survival probability of B; Tachiya (1983)
obviously, the derivation can be directly applied to the pair correlation
function. Rosspeintner  _et al._ (2007) The first-order reaction rate
coefficient of A can be given by,

$\displaystyle k_{\rm s}(t)=-\left[\frac{d}{dt}S(t)\right]/S(t)=-c_{\rm
B}\int_{\sigma}^{\infty}4\pi r^{2}dr\frac{d}{dt}\rho(r,t)=c_{\rm
B}\frac{d}{dt}p(t).$ (A.10)

Here, we consider the second-order reaction rate coefficient given by
$k(t)=k_{\rm s}(t)/c_{\rm B}$.

## References

  * Rice (1985) S. Rice,  _Diffusion-limited Reactions_ , Comprehensive chemical kinetics (Elsevier, 1985).
  * Wilemski and Fixman (1973) G. Wilemski and M. Fixman, “General theory of diffusion-controlled reactions,” J. Chem. Phys. 58, 4009–4019 (1973).
  * Murata and Tachiya (1996) S. Murata and M. Tachiya, “Transient effect in fluorescence quenching by electron transfer. 3. distribution of electron transfer distance in liquid and solid solutions,” J. Phys. Chem. 100, 4064–4070 (1996).
  * Seki, Barzykin, and Tachiya (1999) K. Seki, A. V. Barzykin, and M. Tachiya, “Diffusion-assisted long-range reactions in confined systems: Projection operator approach,” J. Chem. Phys. 110, 7639–7649 (1999).
  * Collins and Kimball (1949) F. C. Collins and G. E. Kimball, “Diffusion-controlled reaction rates,” J. Colloid Sci. 4, 425–437 (1949).
  * Shoup and Szabo (1982) D. Shoup and A. Szabo, “Role of diffusion in ligand binding to macromolecules and cell-bound receptors,” Biophys. J. 40, 33–39 (1982).
  * Traytak (2023) S. D. Traytak, “Accurate analytical calculation of the rate coefficient for the diffusion-controlled reactions due to hyperbolic diffusion,” J. Chem. Phys. 158, 044104 (2023).
  * Lee and Traytak (2023) S. Lee and S. D. Traytak, “Inertial dynamic effects on diffusion-influenced reactions: Approach based on the diffusive Cattaneo system,” J. Chem. Phys 158, 204111 (2023).
  * Tachiya (1986) M. Tachiya, “Influence of the mean free path of reactant particles on the kinetics of diffusion-controlled reactions. II. Rate of bulk recombination,” J. Chem. Phys. 84, 6178–6181 (1986).
  * Dong, Baros, and Andre (1989) W. Dong, F. Baros, and J. C. Andre, “Diffusion-controlled reactions. I. Molecular dynamics simulation of a noncontinuum model,” J. Chem. Phys. 91, 4643–4650 (1989).
  * Zhou and Szabo (1991) H. Zhou and A. Szabo, “Comparison between molecular dynamics simulations and the Smoluchowski theory of reactions in a hard-sphere liquid,” J. Chem. Phys. 95, 5948–5952 (1991).
  * Van Beijeren, Dong, and Bocquet (2001) H. Van Beijeren, W. Dong, and L. Bocquet, “Diffusion-controlled reactions: A revisit of Noyes’ theory,” J. Chem. Phys. 114, 6265–6275 (2001).
  * Yang, Han, and Lee (2001) S. Yang, H. Han, and S. Lee, “An efficient brownian dynamics method for evaluating inertial dynamic effects on diffusion-influenced reactions,” J. Phys. Chem. B 105, 6017–6024 (2001).
  * Litniewski and Gorecki (2004) M. Litniewski and J. Gorecki, “Molecular dynamics tests of the smoluchowski–collins–kimball model for fluorescence quenching of spherical molecules,” Phys. Chem. Chem. Phys. 6, 72–83 (2004).
  * Lee  _et al._ (2004) J. Lee, S. Yang, J. Kim, and S. Lee, “An efficient molecular dynamics simulation method for calculating the diffusion-influenced reaction rates,” J. Chem. Phys. 120, 7564–7575 (2004).
  * Kim  _et al._ (2009) J.-H. Kim, S. Lee, J. Lee, and S. Lee, “Kinetics of collision-induced reactions between hard-sphere reactants,” J. Chem. Phys. 131, 164503 (2009).
  * Piazza, Foffi, and Michele (2013) F. Piazza, G. Foffi, and C. D. Michele, “Irreversible bimolecular reactions with inertia: from the trapping to the target setting at finite densities,” J. Phys. Condens. Matter 25, 245101 (2013).
  * Kramers (1940) H. Kramers, “Brownian motion in a field of force and the diffusion model of chemical reactions,” Physica 7, 284–304 (1940).
  * Wang and Uhlenbeck (1945) M. C. Wang and G. E. Uhlenbeck, “On the theory of the brownian motion ii,” Rev. Mod. Phys. 17, 323–342 (1945).
  * Doi (1975) M. Doi, “Theory of diffusion-controlled reaction between non-simple molecules. ii,” Chem. Phys. 11, 115–121 (1975).
  * Northrup and Hynes (1978) S. H. Northrup and J. T. Hynes, “Coupling of translational and reactive dynamics for a Fokker–Planck model,” J. Chem. Phys. 68, 3203–3216 (1978).
  * Naqvi, Mork, and Waldenstrøm (1982) K. R. Naqvi, K. J. Mork, and S. Waldenstrøm, “Reduction of the fokker–planck equation with an absorbing or reflecting boundary to the diffusion equation and the radiation boundary condition,” Phys. Rev. Lett. 49, 304–307 (1982).
  * Harris (1983a) S. Harris, “The rate of absorption of Brownian particles by a sphere,” J. Phys. Chem. 78, 4698–4700 (1983a).
  * Naqvi, Waldenstrøm, and Mork (1983) K. R. Naqvi, S. Waldenstrøm, and K. J. Mork, “On the rate of absorption of Brownian particles by a black sphere: The connection between the Fokker–Planck equation and the diffusion equation,” J. Chem. Phys. 78, 2710–2712 (1983).
  * Molski (1988) A. Molski, “A source term formalism for the reactive fokker–planck dynamics,” Chem. Phys. Lett. 148, 562–566 (1988).
  * Ibuki and Ueno (1997) K. Ibuki and M. Ueno, “A generalized fokker–planck equation treatment of inertia and non-markovian effects on the short-time dynamics of a collision-induced reaction,” Bull. Chem. Soc. Japan 70, 543–553 (1997).
  * Ibuki and Ueno (2003) K. Ibuki and M. Ueno, “Fokker–planck–kramers equation treatment of dynamics of diffusion-controlled reactions using continuous velocity distribution in three dimensions,” J. Chem. Phys. 119, 7054–7064 (2003).
  * Ibuki and Ueno (2006) K. Ibuki and M. Ueno, “Analysis of short-time transient dynamics of a diffusion-controlled reaction in a hard-sphere fluid based on fokker–planck–kramers equation,” Bull. Chem. Soc. Jpn. 79, 1509–1518 (2006).
  * Harris (1983b) S. Harris, “The rate of absorption of brownian particles by a sphere,” J. Chem. Phys. 78, 4698–4700 (1983b).
  * Harris (1981) S. Harris, “Steady, one-dimensional brownian motion with an absorbing boundary,” J. Chem. Phys. 75, 3103–3106 (1981).
  * Harris (1982) S. Harris, “Steady absorption of brownian particles by a sphere,” J. Chem. Phys. 77, 934–938 (1982).
  * Cattaneo (1958) C. Cattaneo, C.R.Acad.Sci Paris Ser I , 431 (1958).
  * Vernotte (1958) P. Vernotte, “Les paradoxes de la theorie continue de l’equation de la chaleur,” Comptes Rendus 246, 3154 (1958).
  * Kubo, Toda, and Hashitsume (1991) R. Kubo, M. Toda, and N. Hashitsume,  _Statistical Physics II: Nonequilibrium Statistical Mechanics_ (Springer, Berlin, 1991).
  * Kac (1974) M. Kac, “A stochastic model related to the telegrapher’s equation,” Rocky Mt. J. Math. 4, 497 – 510 (1974).
  * Mendez, Fedotov, and Horsthemke (2010) V. Mendez, S. Fedotov, and W. Horsthemke,  _Reaction-Transport Systems: Mesoscopic Foundations, Fronts, and Spatial Instabilities_ , Springer Series in Synergetics (Springer Berlin Heidelberg, 2010).
  * Ghorai and Bairagi (2022) S. Ghorai and N. Bairagi, “Instabilities in hyperbolic reaction–diffusion system with cross diffusion and species–dependent inertia,” Chaos Solitons Fractals 165, 112800 (2022).
  * Wakou, Gorecki, and Kitahara (1998) J. Wakou, J. Gorecki, and K. Kitahara, “On the growth of nonequilibrium spatial correlations in a model reaction diffusion system: The effect of the diffusive flow relaxation,” Acta Phys. Pol. B 29, 1691–1705 (1998).
  * Holmes (1993) E. E. Holmes, “Are diffusion models too simple? a comparison with telegraph models of invasion,” The American Naturalist 142, 779–795 (1993).
  * Tilles and Petrovskii (2019) P. F. C. Tilles and S. V. Petrovskii, “On the consistency of the reaction-telegraph process within finite domains,” J. Stat. Phys. 177, 569–587 (2019).
  * Davies (1954) R. W. Davies, “The connection between the smoluchowski equation and the kramers-chandrasekhar equation,” Phys. Rev. 93, 1169–1170 (1954).
  * Wilemski (1976) G. Wilemski, “On the derivation of Smoluchowski equations with corrections in the classical theory of Brownian motion,” J. Stat. Phys. 14, 153–169 (1976).
  * Meyer and Schröter (1983) J. Meyer and J. Schröter, “Comments on the grad procedure for the Fokker-Planck equation,” J. Stat. Phys. 32, 53–69 (1983).
  * Grad (1949) H. Grad, “On the kinetic theory of rarefied gases,” Commun. Pure and Appl. Math. 2, 331–407 (1949).
  * Traytak (2014) S. D. Traytak, “Asymptotic solution of the diffusion equation in slender impermeable tubes of revolution. I. The leading-term approximation,” J. Chem. Phys. 140, 224102 (2014).
  * Kapral (1978) R. Kapral, “Kinetic theory of chemical reactions in dense fluids,” J. Chem. Phys. 68, 1903–1912 (1978).
  * Burschka and Titulaer (1981) M. A. Burschka and U. M. Titulaer, “The kinetic boundary layer for the fokker–planck equation: Selectively absorbing boundaries,” J. Stat. Phys. 26, 59–71 (1981).
  * Bandyopadhyay, Seki, and Tachiya (2000) T. Bandyopadhyay, K. Seki, and M. Tachiya, “Theoretical analysis of the influence of stochastic gating on the transient effect in fluorescence quenching by electron transfer,” J. Chem. Phys. 112, 2849–2862 (2000).
  * Shin and Kapral (1978) K. J. Shin and R. Kapral, “Kinetic theory of reactive pair dynamics in liquids,” J. Chem. Phys. 69, 3685–3696 (1978).
  * Northrup and Hynes (1979) S. H. Northrup and J. T. Hynes, “Short range caging effects for reactions in solution. I. Reaction rate constants and short range caging picture,” J. Chem. Phys. 71, 871–883 (1979).
  * Tachiya (1983) M. Tachiya, “Theory of diffusion-controlled reactions: Formulation of the bulk reaction rate in terms of the pair probability,” Radiat. Phys. Chem. 21, 167–175 (1983).
  * Rosspeintner  _et al._ (2007) A. Rosspeintner, D. Kattnig, G. Angulo, S. Landgraf, G. Grampp, and A. Cuetos, “On the coherent description of diffusion-influenced fluorescence quenching experiments,” Chem. Eur. J. 13, 6474–6483 (2007).