PHYSICAL REVIEW E 101, 012610 (2020)

Shear viscosity of two-state enzyme solutions
Yuto Hosaka and Shigeyuki Komura*

Department of Chemistry, Graduate School of Science, Tokyo Metropolitan University, Tokyo 192-0397, JapanDavid Andelman†
Raymond and Beverly Sackler School of Physics and Astronomy, Tel Aviv University, Ramat Aviv, Tel Aviv 69978, Israel(Received 16 September 2019; published 28 January 2020)
We discuss the shear viscosity of a Newtonian solution of catalytic enzymes and substrate molecules. The
enzyme is modeled as a two-state dimer consisting of two spherical domains connected with an elastic spring.
The enzymatic conformational dynamics is induced by the substrate binding and such a process is represented byan additional elastic spring. Employing the Boltzmann distribution weighted by the waiting times of enzymatic
species in each catalytic cycle, we obtain the shear viscosity of dilute enzyme solutions as a function of substrateconcentration and its physical properties. The substrate affinity distinguishes between fast and slow enzymes,
and the corresponding viscosity expressions are obtained. Furthermore, we connect the obtained viscosity withthe diffusion coefficient of a tracer particle in enzyme solutions.
DOI: 10.1103/PhysRevE.101.012610

I. INTRODUCTION

Molecular enzymes are nanometer-size proteins that catalyze chemical reactions in the presence of substratemolecules. Here substrates are chemical species that reactwith enzymes and generate product molecules. Catalytic processes that are carried out by molecular enzymes in the cytoplasm and the membrane are essential for cellular metabolism
and homeostasis [1]. In the presence of a substrate, enzymes
undergo conformational changes in each turnover cycle ofthe chemical reaction [2]. In order to mimic actual enzymes,
these conformational dynamics have been simulated usingelastic network models [3–5], and the relationship between
conformational dynamics and the chemical reaction stages has
been studied recently [6].
One of the long-standing and interesting questions in thefield is whether a single enzyme exhibits a motile behavior[7]. Thanks to recent developments of experimental techniques, diffusion phenomena in enzyme solutions have beenstudied by several groups. Using fluorescence correlation
spectroscopy, Muddana et al. [8] reported that diffusion ofa single enzyme is enhanced in the presence of a substrate.
Later on, Riedel et al. [9] showed that the heat released
during turnovers also enhances the enzyme diffusion. Illien
et al. [10], however, revealed experimentally that not only
exothermic enzymes but also endothermic ones contribute to
the diffusion enhancement. In the presence of a gradient in
substrate concentrations, enzymes exhibit collective motionsin the direction of higher or lower concentrations [11,12].
Moreover, the enhanced diffusion of passive objects in enzymatic solutions have been observed independently [13,14].

*†

komura@tmu.ac.jp
andelman@tauex.tau.ac.il

2470-0045/2020/101(1)/012610(11)To understand these experimental findings, several models
have been proposed using equilibrium as well as nonequilibrium approaches. Illien et al. [15] modeled an enzyme
consisting of hydrodynamically coupled subunits, and introduced two discrete equilibrium states corresponding to afree enzyme and a substrate-enzyme complex. They showedthat diffusion of an enzyme is enhanced due to equilibrium
fluctuations [15,16]. Within a nonequilibrium framework,
Golestanian [17] proposed four possible mechanisms leading
to diffusion enhancement by enzymes. They included selfthermophoresis, boost in kinetic energy, stochastic swimming,and collective heating. Mikhailov and Kapral [18,19] modeled
an enzyme as an active force dipole that exerts forces on
the surrounding fluid. When such dipoles are immersed inaqueous fluids, hydrodynamic collective effects due to force
dipoles can lead to diffusion enhancement [18–20].
In spite of these extensive studies on enzyme diffusion,a recent experimental work pointed out the difficulty of accounting quantitatively for the observed enhanced diffusion
within such models as above [21]. Moreover, recent experiments did not observe any change in the diffusion behavior for a specific enzyme that was previously reported to
exhibit enhanced diffusion [22,23]. It was also noticed thatthe viscosity of enzyme solutions is locally reduced while a
specific enzymatic reaction is taking place [7,24]. However,
the effect of enzyme conformational changes on the solution
shear viscosity has not been considered theoretically despite
its importance.
In this paper, we present an analytical study on the shear
viscosity of a dilute enzyme solution under steady shear flow.
As a coarse-grained model of catalytic enzymes, we use thetwo-state dimer model in which conformational changes are
induced by substrate binding and product release [18]. Our
two-state dimer model consists of two hard spheres representing enzymatic domains, which are connected by a harmonic

012610-1©2020 American Physical Society

HOSAKA, KOMURA, AND ANDELMANPHYSICAL REVIEW E 101, 012610 (2020)

In the presence of potential forces, the equation of motionof an overdamped dimer can be written as [27,28]
2
2kB T ∂ ln ψ∂rα
= fα −+ dαβ rβ ,
∂t
ζ
ζ
∂rα

(2)

where ζ is the friction coefficient of the sphere, kB is theBoltzmann constant, T is the temperature, ψ (r, t ) is the timedependent configurational distribution of a dimer, and thevelocity gradient tensor is given bydαβ =FIG. 1. A dilute solution of two-state dimers under steady shear
flow with shear rate γ̇ . Dimers consist of two green spheres of radiusa connected with an elastic spring, and immersed in a Newtonian
fluid having viscosity ηs . The enzymatic reaction, in which a dimer,
a substrate (red circle) and a product (blue square) participated, is
explained in Fig. 2.

spring [18,25,26]. Assuming that the conformational distribution is given by the Boltzmann distribution function, weighted
by the waiting time of an enzyme, we obtain analytically theshear viscosity of a two-state dimer solution as a function of
the substrate concentration. As a result of the competition
between the energy difference of the enzyme two internal
states and the substrate concentration, we find that the enzyme
solution viscosity exhibits a nonmonotonic behavior that depends on the physical properties of the binding substrates. We
shall also connect the obtained viscosity with the diffusion
coefficient of a tracer particle in enzyme solutions.
The outline of our manuscript is the following. In Sec. II,
we review the derivation of the shear viscosity of dimer
solutions originally used to describe polymer solutions. In
Sec. III, we discuss the shear viscosity of a two-state dimer
solution that represents enzyme solutions. We first introduce
the two-state dimer model and discuss the conformational
distribution function of dimers. Analytical results for the
shear viscosity due to dimers and its limiting expressions arepresented. Finally, some discussions and a summary are given
in Sec. IV.
II. VISCOSITY OF DIMER SOLUTIONS
A. Shear viscosity

We consider a dilute solution of dimers under steady shear
flow as schematically depicted in Fig. 1. Here the solvent viscosity is ηs and each dimer is composed of two rigid spheresof radius a, which are connected by an elastic spring. The
positions of two spheres are denoted by the three-dimensional
vectors r1 and r2 . Then, the force acting between the two
spheres within the dimer is given by
fα = −∂U (r),∂rα(1)

where U (r) is the elastic potential energy, r = |r| = |r2 − r1 |is the distance between the two spheres, and rα is the αcomponent of the vector r = (rx , ry , rz ).

∂vα
.
∂rβ

(3)

Notice that vα is the α component of the velocity v =(vx , vy , vz ). Throughout this work, we assume summation
over repeated indices. The second and third terms on theright-hand side of Eq. (2) represent the velocity due to thermal
motion of the solvent and that imposed by the flow field,
respectively.
Such models of dimers have been used extensively tomodel polymer solutions. For polymer solutions, the stresstensor due to the presence of dimers is given [27,28]σαβ = nrα fβ ,(4)where n is the number density (per unit volume) of dimers,and · · ·  denotes the thermal average over all dimer configurations. To calculate the statistical average in Eq. (4),
we introduce the following Fokker-Planck equation for theconformational distribution ψ (r, t ),


∂ 2
∂ψ
2kB T ∂ψ=−fα ψ −+ dαβ rβ ψ .
(5)
∂t
∂rα ζ
ζ ∂rα
In the above, the continuity equation,


∂ψ
∂r
= −∇ ·
ψ ,
∂t
∂t

(6)

where ∇ = (∂rx , ∂ry , ∂rz ) and Eq. (2) have been used. From
the time evolution of rα rβ  in a steady state, the stress tensor
in Eq. (4) can be written as [27,28]nζ[dαγ rβ rγ  + dβγ rα rγ ].
(7)4For simple shear flow whose velocity components are
given by vx = γ̇ ry , vy = vz = 0, where γ̇ is the shear rate (seeFig. 1), the viscosity due to dimers has a simple form,
σxynζ  2 
r .
(8)
η=
=
γ̇
4 y
σαβ = nkB T δαβ +In order to calculate the average ry2 , we need to specify the
conformational distribution function ψ (r).
B. Fraenkel dimer model

Let us first discuss a dimer consisting of two spheres thatare connected by a harmonic spring having an elastic constantK0 , and a natural length 0 . Its potential energy is then given
byK0(r − 0 )2 .
(9)
2This is the “Fraenkel dimer model” [29], and is different thanother polymer dynamic models, such as the Hookean dimer012610-2

U0 (r) =SHEAR VISCOSITY OF TWO-STATE ENZYME SOLUTIONSPHYSICAL REVIEW E 101, 012610 (2020)

model. For Fraenkel dimers, the conformational distributionfunction ψ0 is given by


K0
(10)
(r − 0 )2 ,
ψ0 (r) = C exp −2kB T
where C is the normalization constant. Here, we assume thatthe characteristic relaxation time of a dimer is much smaller
than that of a shear flow, i.e., ζ 20 γ̇ /(kB T )  1. The physicalmeaning of this condition will be separately explained inSec. IV.
Although the shear viscosity of the Fraenkel dimer model
was discussed in Ref. [30], its explicit expression was notderived. By calculating ry2  in Eq. (8) using Eq. (10), we
obtain the shear viscosity for a Fraenkel dimer solution η0 as
η0 ( )
Gτ
√√
2 2 (5+2 )e− + π (3 + 12 +4 2 )[1 + erf ( )]
,
=
√
√
3
4 2 e− + 2 π ( + 2 2 )[1 + erf ( )](11)
where = K0 20 /(2kB T ) is the dimensionless elastic energy,
G = nkB T is the relaxation modulus, τ = ζ /(4K0 ) is the√2xrelaxation time, and erf (x) = (2/ π ) 0 dt e−t is the error
function [31]. Notice that Gτ corresponds to the viscosity of adimer solution when the natural length of the spring vanishes,i.e., = 0 [28,30].
The limiting behaviors of η0 for the Hookean,  1, andstiff Fraenkel dimers,  1, are given by [27,30]
⎧
4
⎪
⎪
 1,
η 0 ( ) ⎨1 + 3 π
=
(12)
⎪
2
Gτ
⎪
⎩
 1.
3
For  1, the viscosity is almost constant, indicating that
thermal energy dominates over elastic energy. For  1, on
the other hand, the viscosity increases linearly with .
III. TWO-STATE DIMER SOLUTIONS
A. Two-state dimer model

Catalytic enzymes undergo conformational changes in thepresence of substrate molecules. To model such situations, weuse a previously proposed two-state dimer model with a stateparameter that can get two values, s = 0 or 1 [18,25,26]. In
Fig 2(a), we schematically illustrate an enzymatic cycle thatis driven by binding a substrate to an enzyme. In the s = 0
state, i.e., the state of the dimer with the elastic constant K0
and the natural length of the spring 0 , this model corresponds
to the Fraenkel dimer model.
When a substrate is supplied to a dimer enzyme whose size
is r = 0 , a transition from s = 0 to s = 1 occurs with the reaction rate k1 . At the same time, the reverse reaction, namely,
the substrate dissociation process, can occur also when r = 0
with the reaction rate k−1 . For the state s = 1, the substrateadds another intradimer interaction, which is modeled as anadditional spring, whose elastic constant and natural length
are K1 and 1 , respectively. Then, the dimer relaxes to a new
equilibrium conformation having the size r = ∗ , as will be

FIG. 2. (a) The enzymatic cycle of the two-state dimer model. A
substrate (red circle) binds to a free enzyme (s = 0) with the reaction
rate k1 (A → B), while its dissociation also occurs with the reaction
rate k−1 (B → A). Once the substrate-enzyme complex (s = 1) isformed, it starts to contract until the equilibrium conformation is
attained (B → C). Then, the product (blue square) is irreversibly
released with the reaction rate kcat , and the bare enzyme comesback to its initial conformation. (b) The schematic illustration of theenergy for a two-state dimer as described by Eq. (13). There are twoenergy branches U (r, 0) and U (r, 1). The transition between themtakes place at r = 0 and r = ∗ , which are the equilibrium valuesof U (r, 0) and U (r, 1), respectively, as are indicated by black circles.
This transition is followed by the downhill relaxational motion alongeach branch. The forward and reverse transition rates, (s = 0) 
(s = 1), are given by k1 , k−1 , respectively, and (s = 1) → (s = 0)
is given by kcat .

explicitly given after Eq. (13). Once the substrate molecule isirreversibly converted to a product molecule with the reaction
rate kcat , a transition from s = 1 to s = 0 takes place at r = ∗ .
Finally, the product is released from the enzyme.
Notice that the reaction rates, k1 , k−1 , and kcat , are the barerate constants that do not depend on the energy difference
between any two states. This also holds for the reaction rates
in the cascade reactions discussed in Appendix A. Moreover,
the transition of a dimer occurs only when r = 0 or r = ∗ ;hence, the reaction rates k1 , k−1 , and kcat are simply taken to
be constant in our model.
The state-dependent total potential energy of this two-state
dimer can be written as
K0
sK1(r − 0 )2 +(r − 1 )2 ,U (r, s) =(13)22
which gives the equilibrium length for s = 1 as ∗ =(K0 0 + K1 1 )/(K0 + K1 ). In Fig. 2(b), we schematically illustrate the energy of a two-state dimer given by Eq. (13)when 0 > 1 . Under this condition, the substrate-enzyme
complex shrinks as compared to the bare enzyme [18,25,26].

012610-3

HOSAKA, KOMURA, AND ANDELMANPHYSICAL REVIEW E 101, 012610 (2020)

In this work, however, we do not require such a condition. In
physiological conditions, the sizes of actual substrate-enzyme
complexes either decrease ( 0 > 1 ) or increase ( 1 > 0 )upon substrate binding [7]. Hereafter, the subscripts “0” and“1” denote physical values for the enzyme and the substrateenzyme complex, respectively.
As represented by the second term in the right-hand sideof Eq. (2), a dimer in our model undergoes conformational
fluctuations due to thermal energy. In other words, a free
enzyme (or a substrate-enzyme complex) fluctuates aroundr = 0 (or r = ∗ ) during turnover cycles. This corresponds tothe situation in which enzymes are subject to thermal motionof solvent molecules. Notice, however, that conformational
fluctuations between multistate enzymes [32,33] are not considered. This is because the original dimer model [18,25,26]
that we employ follows the simple Michaelis-Menten kinetics[see Eq. (14)] with the advantage that the problem becomes
tractable.

where ps (t ) is the time-dependent probability distributionfunction of an enzyme in state s, which will be explicitly
given in Eq. (18). The case of a cascade reaction containing
N substrate-enzyme complexes is discussed in Appendix A asa generalization, and Eq. (15), hence, corresponds to the caseN = 1.
C. Waiting times

Since we consider a dilute solution of two-state dimers, weemploy a single enzyme kinetics to obtain the waiting timethat an enzyme spends at each catalytic step (see also Appendix B). The validity of using a single enzyme kinetics foran enzyme solution will be discussed later in this subsection.
For two-state dimers, the corresponding kinetic equations are
written in terms of the probability functions as [32,33,37,38]
d p0= k−1 p1 − k1 p0 ,dt
d p1
= k1 p0 − (k−1 + kcat )p1 ,dt
d p∗= kcat p1 .
dt

B. Conformational distribution function

The above two-state dimer model describes a chemical
equation following the standard Michaelis-Menten reaction
[34]:
k1

kcat

E + S  ES −→ E∗ + P.

(14)

k−1

This chemical reaction equation describes the enzymatic cycle
composed of three states of an enzyme: a free enzyme (E),a substrate-enzyme complex (ES), and a free enzyme after
the reaction (E∗ ), as depicted in Fig. 2. Furthermore, S and P
stand for the substrate and product, respectively. When dimers
are connected by elastic springs, the time spent during thetransition between these chemical states can be characterized
by a relaxation time τ = ζ /(4K0 ) as introduced after Eq. (11).
For a two-state dimer, we assume that the characteristic
relaxation time is much smaller than that of a shear flow, i.e.,ζ 20 γ̇ /(kB T )  1 as adopted for the Fraenkel dimer modelin Sec. II. We further assume that the transition time spent
between enzymatic states is much smaller than the waitingtime in each of the states, s = 0, 1, i.e., τ /Ws  1, where
the waiting time Ws will be defined later in Eq. (16). This
assumption is justified for enzymes such as adenylate kinase
having a relatively large waiting time, τ /W1 ≈ 0.1 [6]. For
completeness, however, the general case of arbitrary waitingtimes is discussed in Sec. IV. Under these conditions, we
can introduce the Boltzmann distribution function that isweighted only by the waiting time in the respective enzymaticstates. The validity of this assumption has been confirmed
by numerical solutions of the Langevin equation for a single
two-state dimer [26].
The distribution function for the two-state dimer model foran enzyme is then given byψe (r) = 

W0 e−βU (r,0) + W1 e−βU (r,1),
dr [W0 e−βU (r,0) + W1 e−βU (r,1) ]Here, p0 (t ), p1 (t ), and p∗ (t ) are the probability distributionfunctions for the two-state dimer in one of the two states, s =0, 1, and the free enzyme after the catalysis (E∗ ), respectively.
In the above, we have introduced the pseudo-first-order rate
constant k1 = k1 cS , where cS is the time-independent substrateconcentration. Such an assumption is justified when cE  cS
is satisfied, where cE is the enzyme concentration.
By solving the above coupled kinetic equations using theinitial conditions, p0 (0) = 1 and p1 (0) = p∗ (0) = 0, under
the normalization condition p0 (t ) + p1 (t ) + p∗ (t ) = 1, thetime-dependent probability distributions are obtained [37]:1[(a + b − k1 )e(a−b)t + (a − b + k1 )e−(a+b)t ],2akp1 (t ) = 1 [e(a−b)t − e−(a+b)t ],2a
1 (a−b)t
k kcat1 −(a+b)t+ 1, (18)eep∗ (t ) = 1+
2a a − ba+bp0 (t ) =where
a = [(k1 + k−1 + kcat )2 /4 − k1 kcat ]1/2 ,b = (k1 + k−1 + kcat )/2.

(19)

Because a − b < 0 and a + b > 0, both p0 (t ) and p1 (t ) decay
exponentially for t → ∞, and consequently p∗ → 1.
Substituting p0 (t ) and p1 (t ) of Eq. (18) into Eq. (16), we
obtain the waiting times for s = 0 and 1 as
W0 =

(15)

where β = 1/(kB T ). Here the waiting time in the state s isdefined by [35,36]
 ∞dt ps (t ),
(16)
Ws =

(17)

k−1 + kcat,
k1 kcat

W1 =

1
.
kcat

(20)

As a result, the distribution function in Eq. (15) can be written
as0012610-4

ψe (r) = 

e−βU (r,0) + νe−βU (r,1),
dr [e−βU (r,0) + ν e−βU (r,1) ](21)

SHEAR VISCOSITY OF TWO-STATE ENZYME SOLUTIONSPHYSICAL REVIEW E 101, 012610 (2020)

where we have introduced the dimensionless parameter ν,
ν=

cS
k1
cS =
,
k−1 + kcat
KM

(22)

and KM is the Michaelis constant [1],KM =k−1 + kcat
.
k1

(23)

Physically, ν represents the fraction of the s = 1 state during
one turnover cycle of the enzymatic reaction. It depends only
on the substrate concentration and the bare rate constants.
In the following analyses, we vary this state parameter ν toinvestigate the shear viscosity of enzyme solutions. Somenumerical estimates of ν are given in the end of this section.
We discuss here the validity of using a single-enzyme kinetics. In our model, we have assumed that the concentrationof enzymes is small enough so that hydrodynamic interactions
between enzymes are negligible [28]. Such a dilute condition
corresponds to having only a single enzyme in the system,
leading to a renewal process [37]. In the renewal process,the probability distribution function is identically and independently distributed [39]. This means that in every turnovercycle, waiting times follow the same probability distribution,and hence these times can be uniquely determined as shownin Eq. (20).
For systems containing mesoscopic numbers of enzymes,
however, stochasticity in enzymatic reactions plays more important roles as discussed in Refs. [39,40]. Enzyme stochasticity leads to nonrenewal processes and causes breakdownof the Michaelis-Menten equation in steady state [39,40].
Since the waiting time distributions depends on the number ofenzymes for nonrenewal processes, one needs to derive masterequations for waiting time distributions when a solution ofmultiple enzymes is considered [39]. This is beyond the scopeof the present work.

FIG. 3. Contour plot of ηe /(Gτ ) as a function of the parametersν = cS /KM [see Eq. (22)] and κ = K1 /K0 for = K0 20 /(2kB T ) = 1
and λ = 1 / 0 = 1.

and its explicit expression is obtained in Appendix C [seeEq. (C6)]. Specifically, the functions g2 (p, q) (m = 2) andg4 (p, q) (m = 4) are given by√√π (1 + 2pq2 )[1 + erf ( pq)]
q −pq2eg2 (p, q) =
+
,
2p
4p3/2
(28)
and
g4 (p, q) =q(5 + 2pq2 ) −pq2
e
4p2
√
√
π (3 + 12pq2 + 4p2 q4 )[1 + erf ( pq)]+
, (29)
8p5/2

D. Viscosity of two-state dimer solutions

To calculate the shear viscosity of a two-state enzyme
solution, we introduce the following notations: κ = K1 /K0 ,λ = 1 / 0 , and λ∗ = ∗ / 0 = (1 + κλ)/(1 + κ ), where ∗ =(K0 0 + K1 1 )/(K0 + K1 ) is the effective natural length for a
dimer in the s = 1 state. In Appendix C, we show that theviscosity of a two-state enzyme solution is given by
zν
,
(24)
ηe (ν, , κ, λ) = η0 + (η1 − η0 )1 + zν
where the quantity η1 (η0 ) corresponds to the viscosity when
all the enzymes are in the s = 1 (s = 0) stateη1 ( , κ, λ)2 g4 ( (1 + κ ), λ∗ ),=Gτ3 g2 ( (1 + κ ), λ∗ )andz( , κ, λ) = exp −(25)


g2 ( (1 + κ ), λ∗ )κ(λ − 1)2. (26)1+κ
g2 ( , 1)

See also Eq. (11) for the Fraenkel dimer viscosity η0 ( ). In
the above, gm (p, q) is given by an integral, ∞2
dr r m e−p(r−q) ,(27)
gm (p, q) =0

respectively. Equations (24)–(29) for the viscosity are themain result of this work.
In Eq. (26), the factor κ (λ − 1)2 /(1 + κ ) in the exponential function corresponds to the dimensionless energy difference, U ( ∗ , 1) − U ( 0 , 0), between the two equilibrium statesof a two-state dimer with 0 and ∗ , as shown in Fig. 2(b).
Although only the bare reaction rates are taken into account,
the above energy difference naturally emerges by defining theweighted distribution function as in Eq. (21).
When ν = 0, ηe of Eq. (24) simply reduces to η0 , theviscosity of the Fraenkel dimer solution [s = 0, see Eq. (11)].
For ν = 0, the enzyme solution viscosity ηe is determined by
the ratio between the two viscosities η0 and η1 . Due to thefactor z, however, ηe also depends on the energy difference
between the two states of the enzyme. This effect causes anonmonotonic behavior of the viscosity as we will show later.
Before analyzing the behavior of ηe , we estimate typicalvalues of = K0 20 /(2kB T ). The enzyme size can be taken as
0 ≈ 10 nm [1]. Moreover, considering typical forces, 1 pN,
generated by a two-state dimer with size 0 , we estimate thespring constant as K0 ≈ 10−4 N/m [18]. Using these values
and kB T ≈ 4 × 10−21 J in physiological conditions, we obtain≈ 1. Hence, we fix the value hereafter to = 1.

012610-5

HOSAKA, KOMURA, AND ANDELMANPHYSICAL REVIEW E 101, 012610 (2020)

FIG. 4. Contour plot of ηe /(Gτ ) as a function of the parametersν = cS /KM [see Eq. (22)] and λ = 1 / 0 for = K0 20 /(2kB T ) = 1
and κ = K1 /K0 = 1. The white region corresponds to larger absolutevalues of ηe .

In Fig. 3, we present the contour plot of the rescaled viscosity due to two-state dimers, ηe /(Gτ ), as a function of ν and κ
for = λ = 1. One can see that ηe becomes smaller for large
ν and κ, implying that the viscosity decreases when enzymatic
reactions occur more frequently and substrates are stiffer(large K1 ). Notice that stiff dimers lead to a decrease of ηebecause its stiffness suppresses the enzyme size fluctuation. In
Fig. 4, we plot the rescaled viscosity, ηe /(Gτ ), as a function ofν and λ for = κ = 1. Here we see a nonmonotonic behaviorof the viscosity in λ characterized by a peak around λ ≈ 3.2.
Note that for larger λ values, ηe becomes independent of ν.
To see more detailed behavior, we plot in Fig. 5 the rescaled
viscosity, ηe /(Gτ ), as a function of ν for κ = 0.1, 1, and10, while keeping = λ = 1. The dashed line corresponds
to the constant viscosity for a Fraenkel dimer solution, i.e.,

FIG. 6. Plot of ηe /(Gτ ) as a function of the parameter ν for
λ = 1 / 0 = 0.1, 1, 4, and 5.3. The other parameter values are =K0 20 /(2kB T ) = 1 and κ = K1 /K0 = 1. The black dashed line represents η0 in Eq. (11). The blue dotted lines represent the two limiting
expressions in Eq. (30) for λ = 4.

η0 /(Gτ ) ≈ 2.13. We see that ηe decreases with increasing ν
for all the κ values. The decrease of ηe is more enhanced forlarger κ values.
In Fig. 6, we plot ηe as a function of ν for λ = 0.1, 1, 4,
and 5.3, while keeping = κ = 1. We see that ηe shows
both increasing and decreasing dependency as a function ofν depending on the value of λ. When λ = 0.1, 1, and 4, theviscosity ηe increases with λ, reflecting the fact that larger
enzymes lead to higher viscosity. For larger λ such as λ = 5.3,
however, ηe becomes smaller, and as λ is further increased, theviscosity approaches the value of η0 as indicated by the dashedline. In this limit, both Fraenkel dimer solutions and two-state
enzyme solutions exhibit the same viscosity even when ν isvery large.
We discuss now the nonmonotonic behavior of ηe thatis seen in Fig. 6. Such a behavior occurs because z in
Eq. (26) increases for smaller λ, but strongly decreases forlarger λ due to the Gaussian function of Eq. (26). The factorκ (λ − 1)2 /(1 + κ ) in the Gaussian function corresponds to
the rescaled energy difference between the s = 0 and s = 1
states. Hence, it can be regarded as an Arrhenius’ equation that
determines the transition rate from the s = 0 to s = 1 state.
E. Limiting expressions

FIG. 5. Plot of ηe /(Gτ ) as a function of the parameter ν forκ = K1 /K0 = 0.1, 1, and 10. The other parameter values are =K0 20 /(2kB T ) = 1 and λ = 1 / 0 = 1. The black dashed line represents η0 in Eq. (11). The red dotted lines represent the two limiting
expressions in Eq. (30) for κ = 1.

Next, we present the limiting expressions of ηe for small
and large values of the ν parameter, ν  1 and ν  1. The
viscosity of two-state dimer solution in Eq. (24) becomes⎧ν  1,⎨η0 + C1 ν(30)
ηe (ν, , κ, λ) ≈C2
⎩η1 +ν  1,
ν
where C1 ( , κ, λ) = (η1 − η0 )z and C2 ( , κ, λ) = (η0 − η1 )/z. In Figs. 5 and 6, we have plotted the above limits by the red(for κ = 1) and blue (for λ = 4) dotted line, respectively.

012610-6

SHEAR VISCOSITY OF TWO-STATE ENZYME SOLUTIONSPHYSICAL REVIEW E 101, 012610 (2020)λ ≈ 2. Since ηe is inversely proportional to ν in Eq. (30),
positive C2 leads to a decreasing behavior of ηe , whereas
negative C2 results in an increasing behavior.
F. Numerical estimates

FIG. 7. Contour plot of C1 /(Gτ ) [see Eq. (30)] as a functionof κ = K1 /K0 and λ = 1 / 0 for = K0 20 /(2kB T ) = 1 under thecondition ν  1. The quantity C1 changes its sign from negative topositive around λ ≈ 2.

In Fig. 7, we study the ν  1 behavior and plot thecoefficient C1 = (η1 − η0 )z of ν in Eq. (30) as a function ofκ and λ for = 1. The behavior of C1 is nonmonotonic, having a minimum and a maximum around (κ, λ) ≈ (1, 1) and(κ, λ) ≈ (1, 2.5), respectively. The quantity C1 vanishes forlarge λ values, because the Gaussian function in z, Eq. (26),dominates over the viscosity difference, η0 − η1 . Notice thatC1 changes its sign from negative to positive around λ ≈ 2,
where the switching from decreasing to increasing behaviorof ηe as a function of ν occurs.
In Fig. 8, we study the ν  1 behavior and plot the coefficient C2 = (η0 − η1 )/z of ν −1 in Eq. (30) as a function of κand λ when = 1. Here C2 exhibits a monotonic behavior inκ and λ, and changes its sign from positive to negative aroundTo end this section, we give some numerical estimates ofthe parameter ν = cS /KM in Eq. (22). The experimentally
accessible substrate concentration is 10−6 M < cS < 10−3 M[12,14]. On the other hand, the value of the Michaelis constant
KM differs between fast and slow enzymes. For fast enzymes,
such as urease and catalase, it is given by KM ≈ 10−3 M[9,12]. For slow enzymes, such as aldolase and adenylatekinase, it is KM ≈ 10−6 M [6,10]. Hence, the ν range isestimated as 10−3 < ν < 1 and 1 < ν < 103 , respectively,
for fast and slow enzymes. These estimates imply that thelimiting expressions derived for ν  1 and ν  1 in Eq. (30)correspond to these two types of enzymes for cS < 10−4 M
and cS > 10−5 M, respectively.
Next we discuss the values of κ and λ in order to estimate
the viscosity ηe for typical physiological conditions. Since an
enzyme consists of a large complex of macromolecules, thesize of substrate molecules is typically smaller than that of
enzymes [1]. Due to this size difference, the condition λ < 1holds generally. Noncovalent bonds, such as hydrogen bonds,
van der Waals attractions, and hydrophobic forces, are responsible for the formation of macromolecular assemblies. On theother hand, covalent bonds are responsible for the formationof substrate molecules. Then, the molecular flexibilities forthe substrates compared with the enzymes are different, whichleads to the condition κ > 1.
From the above argument, we choose λ = 0.1 and κ = 10.
Using these values and setting = 1, we obtain ηe /(Gτ ) ≈2.11 and ηe /(Gτ ) ≈ 0.39 for fast and slow enzymes, respectively, assuming that the maximum substrate concentrationcS = 10−3 M is attained. Since η0 /(Gτ ) ≈ 2.13 for = 1,
the difference between the enzyme solution with substrates ηeand that without substates η0 is negligible for fast enzymes,
whereas the viscosity ηe is approximately five times smaller
than η0 for slow enzymes.
IV. DISCUSSION AND CONCLUSION

FIG. 8. Contour plot of C2 /(Gτ ) [see Eq. (30)] as a functionof κ = K1 /K0 and λ = 1 / 0 for = K0 20 /(2kB T ) = 1 under thecondition ν  1. The quantity C2 changes its sign from positive
to negative around λ ≈ 2. The white region corresponds to larger
absolute values of C2 .

In this paper, we have investigated the viscosity of dilute two-state enzyme solutions under steady shear flow. We
have obtained the shear viscosity by taking into account the
enzyme conformational changes in a solution with a supplyof substrates. The waiting times, which correspond to therespective conformations of the enzyme, are connected to
the reaction rates in the enzymatic cycle by using the single
enzyme kinetics [37]. In our approach, the two-state dimer
model [18,25,26] and the polymer dimer model [27–29] arecombined.
When the enzyme has the same structural properties as
the substrate, the shear viscosity decreases as the substrateconcentration becomes higher (see Fig. 5). For a substratelarger than the enzyme, the viscosity increases with substrate concentrations (see Fig. 6). When the substrate is large
enough, however, the viscosity reduces to that of a Fraenkel
dimer solution. Furthermore, we have obtained the limiting

012610-7

HOSAKA, KOMURA, AND ANDELMANPHYSICAL REVIEW E 101, 012610 (2020)

expressions of the viscosity for fast and slow enzymes [seeEq. (30)]. For slow enzymes, the coefficient shows only amonotonic behavior. For fast enzymes, on the other hand, thecoefficient of the substrate concentration exhibits a nonmonotonic behavior as functions of the stiffness and size of thesubstrate.
Next, we comment on the connection between the viscosityof a two-state dimer solution and the diffusion coefficient of a
tracer particle in such a solution. By following the discussion
in Refs. [41,42], the diffusion coefficient of a passive sphericalparticle of radius R can be given by Einstein’s relation,
De =

kB T
,
6π (ηs + ηe )R(31)

where we have assumed R  0 . In terms of the enzyme
volume fraction φ = 4π ( 0 /2)3 n/3, De can be expanded up
to first order in φ as


9aηe
kB T
1−
φ .
(32)
De ≈
6π ηs R
2 0 Gτ
Hence, the relative change of the diffusion coefficient withrespect to that of a Fraenkel dimer solution (denoted by D0 ) is
3kB T a(η0 − ηe )φ.
δD = De − D0 =4π ηs R
0 Gτ

(33)

Since η0 > ηe holds for both fast and slow enzymes asestimated before, catalytic enzymes give rise to the diffusion
enhancement under physiological conditions. Moreover, we
see that δD increases as cS is increased in the limits offast and slow enzymes (see Figs. 7 and 8). This behavior
qualitatively agrees with experiments for both tracers andenzymes [8,14,21]. More specifically, using values such as
cS = 10−3 M, a/ 0 = 0.2, φ = 0.1, we obtain that the diffusion increases for slow enzymes as δD/D0 ≈ 0.15. In existing
experiments, however, φ is typically of the order of 10−5 , andhence experimental measurements using higher cE concentration are needed for a more accurately checking of the validityof our model.
Here we discuss how the obtained viscosity is modified
by hydrodynamic effects that have been neglected so far. In
the presence of hydrodynamic interactions, the equation of
motion, Eq. (2), can be rewritten as [43]


∂rα
2
2kB T ∂ ln ψ
+ dαβ rβ , (34)= (δαβ − ζ Gαβ )
fβ −
∂t
ζ
ζ
∂rβ
where Gαβ (r) = (δαβ + rα rβ /r 2 )/(8π ηs r) is the hydrodynamic Oseen tensor [44]. If we assume all orientations tobe equally probable, an equilibrium-averaged hydrodynamic
interaction can be defined by taking the average of Gαβ (r)over all orientations [45],
dr ψ (r)Gαβ (r)
1

,
(35)
h = Tr3
dr ψ (r)where Tr denotes the trace operation. This is called the preaveraging approximation [44]. Then, the equation of motion
can be approximated as∂rα2(1 − ζ h)∂ ln ψ
+ dαβ rβ .
≈fα − kB T
(36)
∂t
ζ
∂rα

Comparing Eqs. (2) and (36), one finds that the change overfrom negligible hydrodynamic interactions to equilibriumaveraged ones can be accomplished by replacing ζ withζ /(1 − ζ h). Hence, for a single-state dimer as in Eq. (8), thehydrodynamic interaction modifies the viscosity by a factorof 1/(1 − ζ h). In Appendix D, we derive h for the Fraenkel
dimer model. When a/ 0 = 0.2 and = 1, for example, wefind that the viscosity is about 20% larger as compared to
the negligible hydrodynamic case. For the two-state dimers,
hydrodynamic effects do not affect the ν dependence of ηealthough some geometrical factors such as κ and λ can enterin h.
In this study, we have assumed that the distribution functions do not depend on shear flow [see Eqs. (10) and (15)].
Here we discuss how these distribution functions are modified
by an external flow and the regime where the flow does notaffect the distributions as assumed in this paper. For a steadystate homogeneous potential flow, Eq. (5) has an analyticalsolution [27],




ζ
U (r)
exp
rα dαβ rβ ,(37)
ψ (r) = C exp −
kB T
kB T
where C is the normalization constant.
For a simple shear flow characterized by a shear rate γ̇ , thedistribution function becomes

 2


ζ r γ̇
U (r)ψ (r, θ , φ, γ̇ ) = C exp −exp
sin2 θ sin 2φ ,
kB T
2kB T
(38)
where rx = r sin θ cos φ and ry = r sin θ sin φ. When thelength of a dimer is r = 0 , the characteristic relaxation time
is given by ζ 20 /(kB T ) [29]. Hence, the shear flow does notaffect the distribution functions when ζ 20 γ̇ /(kB T )  1.
We have assumed that the transition time spent from one
enzymatic species to another is much smaller than the waiting
time, i.e., τ /Ws  1. Here, we consider the general case ofarbitrary waiting time. Because the total times in state s = 0
and s = 1 are given by W0 + τ and W1 + τ1 , respectively, themodified parameter ν becomesν=

k1 (1 + kcat τ1 )cS,
k−1 + kcat (1 + k1 τ cS )(39)

where τ = ζ /(4K0 ) as before and τ1 = ζ /(4K1 ). Since thereverse reaction rate k−1 is negligible in general but may havea finite value, we set it to be a constant. There are only four−1
relevant time scales, namely, kcat
, (k1 cS )−1 , τ , and τ1 , andEq. (39) has four limiting expressions. When the transition
rates are vanishingly small, the modified parameter coincides
with ν in Eq. (22) as it should. For the two intermediate
regimes, Eq. (39) shows linear and inverse dependencies on
the transition time. When the transition time is infinitely large,
we have ν ∼ κ −1 , indicating that the transition dynamics isgoverned only by the relative stiffness between the enzyme
and substrate.
The transition rates can depend on κ and/or λ for general
enzymatic solutions although these effects were not considered in this work. Using Kramers’ reaction-rate theory [46],Aviram et al. [6] obtained free-energy profiles of enzymes byexperimentally measuring the transition rates. In the presenceof such an effect, the enzyme solution viscosity may exhibit

012610-8

SHEAR VISCOSITY OF TWO-STATE ENZYME SOLUTIONSPHYSICAL REVIEW E 101, 012610 (2020)

more complicated dependences on κ and/or λ. Finally, we
have assumed that the viscosity due to enzymes does notdepend on the shear rate. Since the dimer model with finite
natural lengths predicts a viscosity that depends on the shear
rate [27,30], one can extend the present model to a nonNewtonian enzymatic fluid.
ACKNOWLEDGMENTS

We thank R. M. Adar, Y. Avni, V. Démery, K. K. Dey, M.
Doi, T. Kato, A. S. Mikhailov, and K. Yasuda for fruitful discussions and helpful suggestions. Y.H. acknowledges support
by a Grant-in-Aid for JSPS Fellows (Grant No. 19J20271)
from the Japan Society for the Promotion of Science (JSPS).
Y.H. also thanks the hospitality of Tel Aviv University, wherepart of this research was conducted under the TMU-TAU
co-tutorial program. S.K. acknowledges the support by Grantin-Aid for Scientific Research (C) (Grants No. 18K03567 and
No. 19K03765) from the JSPS. D.A. acknowledges support
from the Israel Science Foundation (ISF) under Grant No.
213/19.
APPENDIX A: PROBABILITY DISTRIBUTION FUNCTION
FOR MULTIPLE-STATE ENZYMESIn this Appendix, we generalize the dimer enzyme into a
N-mer one. We derive the probability distribution function
for a single enzyme that has multiple intermediate statesin catalytic chemical reactions. We consider the following
cascade reaction containing N intermediate substrate-enzyme
complexes:
k1

k2

kN

k−1

k−2

k−N

kcat

E + S  (ES)1  · · · (ES)s · · ·  (ES)N −→ E∗ + P.
(A1)
Here (ES)s denotes the sth intermediate complex in the reaction, and ks and k−s are the forward and backward reaction
rates to the states s and s − 1, respectively. At the final step,
the complex is irreversibly converted to an enzyme and a
product with the reaction rate kcat . The enzyme after thecatalysis is denoted by E∗ .
Since we assume that a substrate having the energy Es
binds to (ES)s−1 with the reaction rate ks , the energy of an
enzyme in the state s can be written as
U (r, s) = E0 +

s


Es ,

(A2)

s =1

where E0 is the energy of the free enzyme. Then, the waiting
time-weighted distribution function is given byN
Ws e−βU (r,s)ψN (r) = N s=0 .
(A3)−βU (r,s)s=0 Ws dr e
Here Ws is the waiting time in the state s, which is defined in
Eq. (16).
In order to obtain the viscosity of dimer solutions using
Eq. (8), we need to calculate the second moment ry2 . In
general, the average of any function f (r) over the distributionfunction, Eq. (A3), can be written as f (r)N =  f (r)0 +N
[ f (r)s −  f (r)0 ]s=1

zs0 ws0,×N
1 + s =1 zs 0 ws 0

(A4)

where  f (r)s denotes the average of f (r) over all configurations in the state s,
dr f (r)e−βU (r,s)
 f (r)s = ,(A5)
dr e−βU (r,s)
while zss and wss are defined by∞
dt ps (t )dr e−βU (r,s)


zss =
.
, wss = 0∞
dr e−βU (r,s )dt
ps (t )
0

(A6)

Notice that the quantity z in Eq. (26) corresponds to z10 in theabove notation.
APPENDIX B: MICHAELIS-MENTEN KINETICS ANDSINGLE ENZYME KINETICSIn this Appendix, we briefly review the Michaelis-Menten
kinetics [34] and the single-enzyme kinetics. In the two-state
dimer model, the cascade reaction in Eq. (A1) reduces to theMichaelis-Menten reaction [see Eq. (14)]. In the ensemble ofenzymatic experiments, the corresponding kinetic equations
becomedcE= k−1 cES − k1 cE cS ,dt
dcES
= k1 cE cS − (k−1 + kcat )cES ,
dt
dcP
(B1)
= kcat cES ,dt
where cE and cS were defined before, whereas cES and cP are
the concentrations of substrate-enzyme complex and product,
respectively. By replacing the concentrations of the chemical species with the probability distributions, we obtain the
kinetic equations for a single enzyme as in Eq. (17). In thesteady sate, dcES /dt = 0, the enzymatic velocity is given by
V =dcP
Vmax cS,=dt
KM + cS(B2)

where Vmax = kcat (cE + cES ) is the maximum enzymatic velocity and KM = (k−1 + kcat )/k1 is the Michaelis constantdefined in Eq. (23).
For a single enzyme, the corresponding reaction velocity
can be obtained from the inverse of the total waiting time
during one catalytic cycle. With the use of Eq. (20), thisvelocity becomes
kcat cS
1
1
=
,
(B3)
=
W
W0 + W1
KM + cS
which is termed the single-molecule Michaelis-Menten equation [32]. Comparison of Eqs. (B2) and (B3) yields therelation,
V
1
= .
(B4)
cE + cESW

012610-9

HOSAKA, KOMURA, AND ANDELMANPHYSICAL REVIEW E 101, 012610 (2020)

This relation originates from the equivalence between the average over a single molecule’s long-time trace and that over a large
ensemble of identical molecules, i.e., the ergodicity [32,33].
APPENDIX C: DERIVATION OF ηe

In this Appendix, we present the derivation of ηe in Eq. (24). Using Eq. (21), we calculate ry2  in Eq. (8) asnζ dr ry2 e−βU (r,0) + νe−βU (r,1)

.
ηe =
4
dr [e−βU (r,0) + νe−βU (r,1) ](C1)With the use of Eq. (A4) for N = 1, we obtain 2 
 2   zνnζ  2 ry 0 + r y 1 − ry 0.
ηe =4
1 + zν

(C2)

Since nζ ry2 0 /4 = η0 and nζ ry2 1 /4 = η1 , we obtain Eq. (24). The viscosity of a Fraenkel dimer solution η0 is given by Eq. (11).
Next we calculate η1 in Eq. (25) as∞nζ dr ry2 e−βU (r,1)
nζ 0 dr r 4 e−βU (r,1)


=
.
(C3)
η1 =
4
12 0∞ dr r 2 e−βU (r,1)
dr e−βU (r,1)
For a harmonic potential, the integration of r m can be generally expressed as
 ∞
 ∞
m

2
2
dr r m e−p(r−q) =du (u + q)m e−pu =gm (p, q) =
0

m!
qm−n
(m
−
n)!n!
n=0

−q

 ∞
−q

du un e−pu .
2The last integral can be further performed as follows.


 ∞
 ∞
 0
 pq2
p−(n+1)/2
n −pu2n −pu2n(n+1)/2−1 −t(n+1)/2−1 −t
du u e+du u e
=
dt t
e +
dt t
e ,
(−1)
2
−q
0
0
0


 ∞
 ∞
p−(n+1)/2[1 + (−1)n ]=dt t (n+1)/2−1 e−t − (−1)ndt t (n+1)/2−1 e−t .
20pq2
Finally, gm (p, q) becomesgm (p, q) =






m
m!
1
n+1
n+1
p−(n+1)/2 qm−n [1 + (−1)n ]
− (−1)n , pq2 ,
2 n=0 (m − n)!n!
2
2

∞
∞
where (x) = 0 dt t x−1 e−t and (x, α) = α dt t x−1 e−t are
the gamma function and the incomplete gamma function of
the second kind, respectively [31].

APPENDIX D: HYDRODYNAMIC INTERACTIONSBETWEEN TWO SPHERES

(C4)

(C5)

(C6)

By taking m = 1, 2 in gm (p, q), Eq. (C6), the dimensionlesscombination ζ h is obtained asa g1 ( , 1)
ζ h( ) =0 g2 ( , 1)√√e− + π [1 + erf ( )]a=√√ . (D2)− +π [1 + 1/(2 )][1 + erf ( )]
0 e

(D1)

For large dimers, a/ 0  1, the hydrodynamic effects become negligible. The limiting behavior of h for the Hookean, 1, and stiff Fraenkel dimers,  1, is given, respectively, by
⎧
2a
⎪
⎪
⎨
 1,
(D3)
ζ h( ) = a0 π 1
⎪
⎪
 1.
⎩
0 1 + 1/(2 )[1] B. Alberts, A. Johnson, P. Walter, J. Lewis, and M. Raff,
Molecular Biology of the Cell (Garland Science, New York,
2008).
[2] M. Gerstein, A. M. Lesk, and C. Chothia, Biochemistry 33,
6739 (1994).

[3] Y.
TogashiandA.
S.
Mikhailov,
Proc.
Natl.
Acad.
Sci.
USA
104,
8697
(2007).
[4] T. Sakaue, R. Kapral, and A. S. Mikhailov, Eur. Phys. J. B 75,
381 (2010).

In this Appendix, we present the calculation of Eq. (35)for the Fraenkel dimer model. With the assumption that thefluid is isotropic, the Oseen tensor becomes δαβ /(6π ηs r).
Substituting it into Eq. (35) yields
dr ψ0 (r)/r
1

.
h=
6π ηs
dr ψ0 (r)012610-10

SHEAR VISCOSITY OF TWO-STATE ENZYME SOLUTIONS[5] C. Echeverria, Y. Togashi, A. S. Mikhailov, and R. Kapral,
Phys. Chem. Chem. Phys. 13, 10527 (2011).
[6] H. Y. Aviram, M. Pirchi, H. Mazal, Y. Barak, I. Riven, and G.
Haran, Proc. Natl. Acad. Sci. USA 115, 3243 (2018).
[7] Y. Zhang and H. Hess, ACS Cent. Sci. 5, 939 (2019).
[8] H. S. Muddana, S. Sengupta, T. E. Mallouk, A. Sen, and P. J.
Butler, J. Am. Chem. Soc. 132, 2110 (2010).
[9] C. Riedel, R. Gabizon, C. A. M. Wilson, K. Hamadani, K.
Tsekouras, S. Marqusee, S. Pressé, and C. Bustamante, Nature(London) 517, 227 (2015).
[10] P. Illien, X. Zhao, K. K. Dey, P. J. Butler, A. Sen, and R.
Golestanian, Nano Lett. 17, 4415 (2017).
[11] S. Sengupta, K. K. Dey, H. S. Muddana, T. Tabouillot, M. E.
Ibele, P. J. Butler, and A. Sen, J. Am. Chem. Soc. 135, 1406
(2013).
[12] A.-Y. Jee, S. Dutta, Y.-K. Cho, T. Tlusty, and S. Granick, Proc.
Natl. Acad. Sci. USA 115, 14 (2018).
[13] K. K. Dey, F. Y. Pong, J. Breffke, R. Pavlick, E. Hatzakis,
C. Pacheco, and A. Sen, Angew. Chem. Int. Ed. 55, 1113
(2016).
[14] X. Zhao, K. K. Dey, S. Jeganathan, P. J. Butler, U. M. CórdovaFigueroa, and A. Sen, Nano Lett. 17, 4807 (2017).
[15] P. Illien, T. Adeleke-Larodo, and R. Golestanian, Europhys.
Lett. 119, 40002 (2017).
[16] T. Adeleke-Larodo, P. Illien, and R. Golestanian, Eur. Phys. J.
E 42, 39 (2019).
[17] R. Golestanian, Phys. Rev. Lett. 115, 108102 (2015).
[18] A. S. Mikhailov and R. Kapral, Proc. Natl. Acad. Sci. USA 112,
E3639 (2015).
[19] R. Kapral and A. S. Mikhailov, Physica D 318-319, 100 (2016).
[20] Y. Hosaka, K. Yasuda, R. Okamoto, and S. Komura, Phys. Rev.
E 95, 052407 (2017).
[21] M. Xu, J. L. Ross, L. Valdez, and A. Sen, Phys. Rev. Lett. 123,
128101 (2019).
[22] Y. Zhang, M. J. Armstrong, N. M. B. Kazeruni, and H. Hess,
Nano Lett. 18, 8025 (2018).
[23] J.-P. Günther, G. Majer, and P. Fischer, J. Chem. Phys. 150,
124201 (2019).
[24] V. Armoškaitė, K. Ramanauskienė, and V. Briedis, Afr. J. Pham.
Pharmacol. 6, 1685 (2012).

PHYSICAL REVIEW E 101, 012610 (2020)[25] H. Flechsig and A. S. Mikhailov, J. R. Soc., Interface 16,
20190244 (2019).
[26] Y. Hosaka, S. Komura, and A. S. Mikhailov (unpublished).
[27] R. B. Bird, R. C. Armstrong, O. Hassager, and C. F. Curtiss,
Dynamics of Polymeric Liquids, Vol. 2 (Wiley, New York,
1987).
[28] M. Doi, Soft Matter Physics (Oxford University Press, Oxford,2013).
[29] G. K. Fraenkel, J. Chem. Phys. 20, 642 (1952).
[30] R. B. Bird, C. F. Curtiss, and K. J. Beers, Rheol. Acta 36, 269(1997).
[31] M. Abramowitz and I. A. Stegun, Handbook of Mathematical
Functions (Dover, New York, 1972).
[32] S. C. Kou, B. J. Cherayil, W. Min, B. P. English, and X. S. Xie,
J. Phys. Chem. 109, 19068 (2005).
[33] B. P. English, W. Min, A. M. van Oijen, K. T. Lee, G. Luo, H.
Sun, B. J. Cherayil, S. C. Kou, and X. S. Xie, Nat. Chem. Bio.
2, 87 (2006).
[34] L. Michaelis and M. L. Menten, Biochem. Z. 49, 333-369(1913).
[35] J. Cao, J. Phys. Chem. B 115, 5493 (2011).
[36] N. G. van Kampen, Stochastic Processes in Physics and Chemistry (Elsevier Science, New York, 1992).
[37] H. P. Lu, L. Xun, and X. S. Xie, Science 282, 1877 (1998).
[38] S. Xie, Single Mol. 2, 229 (2001).
[39] S. Saha, S. Ghose, R. Adhikari, and A. Dua, Phys. Rev. Lett.
107, 218301 (2011).
[40] R. Grima, Phys. Rev. Lett. 102, 218103 (2009).
[41] N. Oppenheimer and H. Diamant, Biophys. J. 96, 3041
(2009).
[42] N. Oppenheimer and H. Diamant, Phys. Rev. E 82, 041912(2010).
[43] R. B. Bird and H. R. Warner, Trans. Soc. Rheol. 15, 741(1971).
[44] M. Doi and S. F. Edwards, The Theory of Polymer Dynamics(Oxford University Press, New York, 1986).
[45] H. R. Warner, Ind. Eng. Chem. Fundam. 11, 379(1972).
[46] P. Hänggi, P. Talkner, and M. Brokovec, Rev. Mod. Phys. 62,
251 (1990).

012610-11