#  Radiative neutron capture reaction rates for r-process nucleosynthesis

Vinay Singh1, Joydev Lahiri2, Malay Kanti Dey3 and D. N. Basu4 Variable Energy
Cyclotron Centre, 1/AF Bidhan Nagar, Kolkata 700064, INDIA<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

About half of the elements beyond iron are synthesized in stars by rapid-
neutron capture process (r-process). The stellar environment provides very
high neutron flux in a short time ($\sim$ seconds) which is conducive for the
creation of progressively neutron-rich nuclei till the waiting point is
reached after which no further neutron capture reactions proceed. At this
point such extremely neutron-rich nuclei become stable via $\beta^{-}$ decay.
A detailed understanding of the r-process remains illusive. In the present
work, we explore the radiative neutron-capture (n,$\gamma$) cross sections and
reaction rates around the r-process peak near mass number eighty. The inherent
uncertainties remain large in some cases, particularly in case of neutron-rich
nuclei. When the low-energy enhancement exists, it results in significant
increase in the reaction rate for neutron-capture. Keywords: Binding energies
$\&$ masses, r-process, (n,$\gamma$) cross sections, Level density,
Nucleosynthesis.

###### pacs:

23.40.-s; 24.10.-i; 26.30.-k; 26.50.+x; 96.10.+i; 97.10.Cv; 98.80.Ft; 21.10.Dr

##  I Introduction

The nucleosynthesis via rapid neutron capture process (r-process) occurring in
astrophysical environment is regarded as carrying the load of producing nearly
half of the heavy elements abundant in nature. This broad-spectrum problem
poses a great challenge for physics. The properties of nuclei like masses,
neutron capture and $\beta$-decay rates along with $\beta$-delayed neutron
emission probabilities are crucial physical quantities which are required for
the calculations of r-process nucleosynthesis. While nuclear properties in the
neighborhood of stability are well known, much remains to be explored in the
domain of neutron-rich nuclei far away from the stability line which may aid
the r-process. The studies of sensitivity ascertain the response of a change
in nuclear physics input(s). It facilitates fixing of critical nuclear
properties that decide observed final elemental abundances prevalent in
nature.

The recognition of locations responsible for the synthesis of heavy elements
DC03 ; DC13  is a matter of great interest in astrophysics. It is well known
for decades that the nuclear abundances in our solar system of elements which
are heavier than iron can be partitioned in two major processes that
synthesize them. In the slow neutron capture process (s-process), substantial
time elapses for the $\beta$-decay to follow before capturing another neutron.
This process produces nuclei around the line of stability. On the contrary,
rapid neutron capture process (r-process) is rather fast in timescale
producing nuclei on the neutron-rich side of the valley of stability. There
also exists another process called the proton capture process (p-process) that
produces primarily elements with low abundances on the proton-rich side Bu57 ;
Ca57 . While the basic physics Wa97  of the s-process Ka11  and heavy
p-process Ar03 ; Th10  is better understood and astrophysical locations well
identified, the same do not hold in case of the r-process Qi07 ; Ar07 ; Th11 .

In principle, by subtracting the p- and s-process contributions An89 ; Ar99
from the solar system abundances, the r-process pattern can be extracted which
consists of three predominant peaks in the abundances at mass numbers $\sim$
80, 130 and 195, linked with the neutron magic numbers at 50, 82 and 126. Of
the order of 100 neutrons per seed nucleus are necessary to produce the
heaviest r-process elements. Additional constraints arise from observations of
r-process elements in old stars in the galactic halo Sn08 ; Ro14  and
meteoritic data Wa96 . These data point towards specific origins for the
r-process nuclei of light mass number $<$ 120 and heavy mass number $>$ 120\.
The distribution of main r-process elemental abundances is much the same among
r-process enhanced halo stars and is similar to the residuals of the solar
system. This is indicative of the fact that since early galactic times Ma90 ;
Ar04 ; Ko14 ; Ma15  identical process created these elements that operates in
a consistent manner. Though the timescale argument is obeyed by the core-
collapse supernovae, the promise of early studies Me92 ; Wo94  is yet to be
achieved by the simulations Arc07 ; Fi10 ; Hu10 ; Ro12  of recent times. The
neutron star mergers cause low temperature outflows which are highly neutron-
rich La74 ; Me89 ; Fr99 ; Wa14 ; Ju15  and due to fission recycling Go11 ;
Ko12  have consistent abundance pattern. But mergers suffer from time delay
Wa15  uncertainty. Other major astrophysical sites explored comprise hot
accretion disk outflows from the mergers of neutron stars or neutron star-
black hole Wa14 ; Ju15 ; Su08 ; Pe14 , supernova neutron-rich jets Wi12 ; No15
; Ts15 , gamma-ray burst collapsar outflows Pr03 ; Su06 ; Ma12 , shocked
surface layers of O-Ne-Mg cores Wa03 ; Ja08  and neutrino-induced
nucleosynthesis in the helium shell of exploding massive stars Ba11 .

The above mentioned astrophysical objects are governed by definitive
environments such as initial composition, temporal distribution of density and
temperature and neutron-richness which lead to elemental abundances that are
unique and bear the signatures of these conditions. The preciseness and near-
universality Ro10  of the pattern should have led to the main astrophysical
sites for the r-process. But the uncertainties of the nuclear rates used in
the network calculations for predicting the r-process estimates would cause
large scale correlations to be unreliable. Moreover, the complications
involved in simulating astrophysical environments and idiosyncrasies of the
properties of a large number of neutron-rich nuclei participating in the
r-process lead to additional uncertainties. In the present work, the radiative
neutron-capture (n,$\gamma$) cross sections and reaction rates around the
r-process peak near mass number eighty have been calculated and inherent
uncertainties in case of neutron-rich nuclei have been investigated.

##  II Theoretical formalism

The thermonuclear reaction rates can be obtained by convoluting fusion cross
sections with Maxwell-Boltzmann distribution of energies. These cross sections
can vary by several orders of magnitude across the required energy range. The
low energy fusion cross sections $\sigma$, some of which are not sufficiently
well known, can be obtained from laboratory experiments. However, there are
cases, in particular involving the weak interaction such as the basic p$+$ p
fusion to deuterium in the solar p-p chain, where no experimental data are
available and one completely relies on theoretical calculations Ad11 . The
theoretical estimates of the thermonuclear reaction rates depend on the
various approximations used. Several factors influence the measured values of
the cross sections. We need to account for the Maxwellian-averaged
thermonuclear reaction rates in the network calculations used in primordial
and stellar nucleosynthesis.

The reaction rates used in the Big Bang Nucleosynthesis (BBN) reaction network
have temperature dependences except 6Li(n,$\gamma$)7Li, 10B(n,$\gamma$)11B,
12C(n,$\gamma$)13C and 14N(n,$\gamma$)15N which are constant with respect to
temperature. The computer code TALYS Talys  allows a comprehensive
astrophysical reaction rate calculations apart from other nuclear physics
calculations. To a good approximation, in the interior of stars the assumption
of a thermodynamic equilibrium holds and nuclei exist both in the ground and
excited states. This assumption along with cross sections calculated from
compound nucleus model for various excited states facilitates Maxwellian-
averaged reaction rates. For stellar evolution models this is quite an
important input. The nuclear reaction rates are generally evaluated using the
statistical model Rauscher97 ; Rauscher10  and astrophysical calculations
mostly use these reaction rates. Stellar reaction rate calculations have been
routinely done in past Rauscher00 ; Rauscher01 . However, TALYS has extended
these Hauser-Feshbach statistical model (HF) Ha52  calculations by adding some
new and important features. Apart from coherent inclusion of fission channel
it also includes reaction mechanism that occurs before equilibrium is reached,
multi-particle emission, competition among all open channels, width
fluctuation corrections in detail, coupled channel description in case of
deformed nuclei and level densities that are parity-dependent. The nuclear
models are also normalized for available experimental data using separate
approaches such as on photo-absorption data, the E1 resonance strength or on
s-wave spacings, the level densities.

###  II.1  Cross sections for radiative capture

In the low energy domain, compound nucleus is formed by the fusion of the
projectile and the target nuclei. While the total energy $E^{tot}$ is fixed
from energy conservation, the total spin $J$ and parity $\Pi$ can have a range
of values. The reaction obeys the following conservation laws,

$E_{a}+S_{a}\ =\ E_{a^{\prime}}+E_{x}+S_{a^{\prime}}=E^{tot},~{}~{}~{}~{}{\rm
energy~{}conservation},$ $s+I+l\ =\
s^{\prime}+I^{\prime}+l^{\prime}=J,~{}~{}~{}~{}{\rm
angular~{}momentum~{}conservation},$

$\pi_{0}\Pi_{0}(-1)^{l}\ =\pi_{f}\Pi_{f}(-1)^{l^{\prime}}=\Pi,~{}~{}~{}~{}{\rm
parity~{}conservation}.$

The formula for binary cross section, assuming the compound nucleus model, is
given by

$\displaystyle\sigma_{\alpha\alpha^{\prime}}^{comp}\ $ $\displaystyle=\
D^{comp}\frac{\pi}{k^{2}}\sum_{J=mod(I+s,1)}^{l_{\max}+I+s}\sum_{\Pi=-1}^{1}\frac{2J+1}{(2I+1)(2s+1)}$
(1)
$\displaystyle\sum_{j=|J-I|}^{J+I}\sum_{l=|j-s|}^{j+s}\sum_{j^{\prime}=|J-I^{\prime}|}^{J+I^{\prime}}\sum_{l^{\prime}=|j^{\prime}-s^{\prime}|}^{j^{\prime}+s^{\prime}}\delta_{\pi}(\alpha)\delta_{\pi}(\alpha^{\prime})$

$\displaystyle\times\ \frac{T_{\alpha lj}^{J}(E_{a})\langle
T_{\alpha^{\prime}l^{\prime}j^{\prime}}^{J}(E_{a^{\prime}})\rangle}{\sum_{\alpha^{\prime\prime},l^{\prime\prime},j^{\prime\prime}}\delta_{\pi}(\alpha^{\prime\prime})\langle
T_{\alpha^{\prime\prime}l^{\prime\prime}j^{\prime\prime}}^{J}(E_{a^{\prime\prime}})\rangle}W_{\alpha
lj\alpha^{\prime}l^{\prime}j^{\prime}}^{J}$

where

$E_{a}=$ the energy of the projectile

$l=$ the orbital angular momentum of the projectile

$s=$ the spin of the projectile

$j=$ the total angular momentum of the projectile

$\pi_{0}=$ the parity of the projectile

$\delta_{\pi}(\alpha)=\left\\{\begin{array}[]{ll}1{}{}&{\rm
if}~{}~{}(-1)^{l}\pi_{0}\Pi_{0}=\Pi\\\ 0{}{}&{\rm
otherwise}\end{array}\right.$

$\alpha=$ the designation of the channel for the initial projectile-target
system:

$\alpha=\\{a,\ s,\ E_{a},\ E_{x}^{0},\ I,\ \Pi_{0}\\}$, where $a$ and
$E_{x}^{0}$ are the type of the projectile and the excitation energy (which is
zero usually) of the target nucleus, respectively

$l_{\max}=$ the maximum l-value of the projectile

$S_{a}=$ the separation energy

$E_{a^{\prime}}=$ the energy of the ejectile

$l^{\prime}=$ the orbital angular momentum of the ejectile

$s^{\prime}=$ the spin of the ejectile

$j^{\prime}=$ the total angular momentum of the ejectile

$\pi_{f}=$ the parity of the ejectile

$\delta_{\pi}(\alpha^{\prime})=\left\\{\begin{array}[]{ll}1{}{}&{\rm
if}~{}~{}(-1)^{l^{\prime}}\pi_{f}\Pi_{f}=\Pi\\\ 0{}{}&{\rm
otherwise}\end{array}\right.$

$\alpha^{\prime}=$ the designation of channel for the ejectile-residual
nucleus final system:

$\alpha^{\prime}=\\{a^{\prime},\ s^{\prime},\ E_{a^{\prime}},\ E_{x},\
I^{\prime},\ \Pi_{f}\\}$, where $a^{\prime}$ and $E_{x}$ are the type of the
ejectile and the residual nucleus excitation energy, respectively

$I=$ the spin of target nucleus

$\Pi_{0}=$ the parity of target nucleus

$I^{\prime}=$ the spin of residual nucleus

$\Pi_{f}=$ the parity of residual nucleus

$J=$ the total angular momentum of the compound system

$\Pi=$ the parity of the compound system

$D^{comp}=$ the depletion factor so as to take into account for pre-
equilibrium and direct effects

$k=$ the wave number of the relative motion

$T=$ the transmission coefficient

$W=$ the correction factor for width fluctuation (WFC).

Figure 1: Plots of the Hauser-Feshbach estimates of (n,$\gamma$) cross
sections as functions of incident neutron energy for Fe isotopes.

###  II.2  Radiative neutron capture reaction rates

The velocities of both the targets and projectiles obey Maxwell- Boltzmann
distributions corresponding to ionic plasma temperature $T$ at the site. The
astrophysical nuclear reaction rate can be calculated by folding the Maxwell-
Boltzmann energy distribution for energies $E$ at the given temperature $T$
with the cross section given by Eq.(2). Additionally, target nuclei exist both
in ground and excited states. The relative populations of various energy
states of nuclei with excitation energies $E_{x}^{\mu}$ and spins $I^{\mu}$ in
thermodynamic equilibrium follows the Maxwell-Boltzmann distribution. In order
to distinguish between different excited states the superscript $\mu$ is used
along with the incident $\alpha$ channel in the formulas that follow. Taking
due account of various target nuclei excited state contributions, the
effective nuclear reaction rate in the entrance channel
$\alpha\rightarrow\alpha^{\prime}$ can be finally expressed as

$N_{A}\langle\sigma v\rangle_{\alpha\alpha^{\prime}}^{*}(T)=\left(\frac{8}{\pi
m}\right)^{1/2}\frac{N_{A}}{(kT)^{3/2}G(T)}\times\ $ (2)

$\displaystyle\int_{0}^{\infty}\sum_{\mu}\frac{(2I^{\mu}+1)}{(2I^{0}+1)}\sigma_{\alpha\alpha^{\prime}}^{\mu}(E)E\exp\left(-\frac{E+E_{x}^{\mu}}{kT}\right)dE,$

where $N_{A}$ is the Avogadro number which is equal to 6.023$\times 10^{23}$,
$k$ and $m$ are the Boltzmann constant and the reduced mass in the $\alpha$
channel, respectively, and

$G(T)=\displaystyle\sum_{\mu}(2I^{\mu}+1)/(2I^{0}+1)\exp(-E_{x}^{\mu}/kT)$

is the temperature dependent normalized partition function. By making use of
the reciprocity theorem Ho76 , the reverse reaction cross sections or rates
can also be estimated.

Figure 2: Plots of the Hauser-Feshbach estimates of (n,$\gamma$) cross
sections as functions of incident neutron energy for 59Co, 74Ge and 75As.

##  III  The nuclear physics of r-process

In nuclear astrophysics, the r-process is a series of nuclear reactions which
is responsible for the synthesis of approximately half of heavy nuclei
resulting origination of elements beyond iron. The neutron star mergers, Type
II supernovae and low-mass supernovae Ba06  are thought to be the three
probable r-process candidate sites where the appropriate conditions necessary
for nucleosynthesis are expected to exist.

Figure 3: Plots of the Hauser-Feshbach estimates of (n,$\gamma$) cross
sections as functions of incident neutron energy for Ni isotopes.

In a Type II supernova, shortly after the extreme compression of electrons
$\beta^{-}$ decay is forbidden. This happens due to the high electron density
which occupies all the states available to free electrons up to a Fermi energy
that is larger than the $\beta^{-}$ decay energy. However, these free
electrons capture by nuclei do occur causing increasing neutronization of
nuclear matter. This fact causes high temperatures and extremely high density
$\sim$ 1024 of free neutrons per cm3 Bu57  that can not decay. At this stage
it bloats and the expansion cools it and neutron capture by heavy nuclei
present in environment advances much rapidly than the $\beta^{-}$ decay.
Consequently, the r-process progresses along the path of the neutron drip line
and neutron-rich unstable nuclei are synthesized.

The dampening process of photo-disintegration, the rapid decrease in the
neutron-capture cross section as it approaches closed neutron shells and the
magnitude of nuclear stability in the domain of heavy-isotopes are the three
major processes which influence the ascent of neutron drip line. Thus the
neutron-rich, weakly bound nuclei are formed with neutron separation energies
as low as 2 MeV Bu57 ; Th04  due to neutron captures during r-process
nucleosynthesis. At this stage, the neutron capture is temporarily paused as
it reaches the neutron shell closures at N = 50, 82 and 126. These pauses or
the so-called waiting points are characterized by increased binding energy
compared to heavier isotopes. As a consequence, it results in low neutron
capture cross sections and a growth of semi-magic nuclei which are more stable
towards $\beta$ decay. Additionally, nuclei beyond the shell closures, owing
to their proximity to the drip line, are inclined to $\beta$ decay briskly.
The $\beta$ decay happens before further neutron capture Wang15  for these
nuclei. The waiting point nuclei then prefer $\beta$ decay to move toward
stability ahead of any more neutron capture Bu57 , causing deceleration or
freeze-out of the reaction.

The decrease in the stability of nuclei puts an end to the r-process when its
heaviest nuclei develop instability towards spontaneous fission. At this point
the total number of nucleons reaches two hundred seventy but before that
sufficiently low fission barrier might induce fission at neutron capture and
continuing towards neutron drip line Bo72  terminates. After the decrease in
neutron flux, these radioactive nuclei which are highly unstable undergo
successive $\beta$ decays rapidly until they approach neutron-rich nuclei
which are more stable Cl68 . In neutron-rich predecessor nuclei, the r-process
produces an abundance pattern of radioactive nuclei $\sim$ 10 amu below the
peaks of slow neutron-capture process after their decays back to stability
whereas the slow neutron-capture process produces an abundance of closed
neutron shell stable nuclei.

Figure 4: Plots of the Hauser-Feshbach estimates of (n,$\gamma$) cross
sections as functions of incident neutron energy for Cu and Ga isotopes.

Figure 5: Plots of the Hauser-Feshbach estimates of (n,$\gamma$) cross
sections as functions of incident neutron energy for Zn isotopes.

Figure 6: Plots of the Hauser-Feshbach estimates of (n,$\gamma$) cross
sections as functions of incident neutron energy for Se isotopes.

Figure 7: Plots of the (n,$\gamma$) reaction rates as functions of T9 for Fe
isotopes.

Figure 8: Plots of the (n,$\gamma$) reaction rates as functions of T9 for
59Co, 74Ge and 75As.

Figure 9: Plots of the (n,$\gamma$) reaction rates as functions of T9 for Ni
isotopes.

Figure 10: Plots of the (n,$\gamma$) reaction rates as functions of T9 for Cu
and Ga isotopes.

Figure 11: Plots of the (n,$\gamma$) reaction rates as functions of T9 for Zn
isotopes.

Figure 12: Plots of the (n,$\gamma$) reaction rates as functions of T9 for Se
isotopes.

Figure 13: Comparison of the predictions of 56Fe(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

Figure 14: Comparison of the predictions of 59Co(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

Figure 15: Comparison of the predictions of 62Ni(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

Figure 16: Comparison of the predictions of 65Cu(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

Figure 17: Comparison of the predictions of 66Zn(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

Figure 18: Comparison of the predictions of 69Ga(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

Figure 19: Comparison of the predictions of 74Ge(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

Figure 20: Comparison of the predictions of 75As(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

Figure 21: Comparison of the predictions of 78Se(n,$\gamma$) reaction rates
(HF) with the data from JINA REACLIB kakd ; ks .

##  IV  Results and discussion

###  IV.1  Radiative neutron capture cross section calculations

We have calculated theoretically the radiative neutron capture (n,$\gamma$)
cross sections of astrophysical importance for Fe, Co, Ni, Cu, Zn, Ga, Ge, As
and Se isotopes using the Hauser-Feshbach statistical model reaction
calculations. The calculations have been performed using the most recent level
density based on temperature-dependent Hartree-Fock-Bogolyubov calculations
using the Gogny force and for the gamma-ray strength function Brink-Axel
Lorentzian has been used Talys . In Figs.1-6, the plots of these estimates for
the radiative neutron capture cross sections as functions of incident neutron
energy for the above mentioned elements and their different isotopes have been
illustrated.

###  IV.2  Reaction rate calculations for neutron rich nuclei

The cross sections for nuclear reaction and its convolution with Maxwell-
Boltzmann distribution of energies are important for the explanation of
various processes occurring under extreme conditions Bu57 ; Fo64 ; Cl83 . In
the main-sequence stars and compact stars which are in their ultimate stages
of evolution, such environments of very high density and temperature prevail.
The exothermic fusion reactions causes nuclear explosions in the surface
layers of the accreting white dwarfs (nova events), in the cores of massive
accreting white dwarfs (type Ia supernovae) Ni97 ; Ho06  and in the surface
layers of accreting neutron stars (type I X-ray bursts and superbursts St06 ;
Sc03 ; Cu06 ; Gu07 ). Precise knowledge of the rates of thermonuclear
reactions obtained by folding Maxwell-Boltzmann distribution of energies with
energy dependent cross sections becomes necessary for describing these
astrophysical phenomena. The Maxwellian-averaged thermonuclear reaction rate
per particle pair $<\sigma v>$ at temperature $T$, can be represented by the
integral Ad11 ; Fo67 ; Bo08  described below:

$<\sigma v>=\Big{[}\frac{8}{\pi
m(kT)^{3}}\Big{]}^{1/2}\int\sigma(E)E\exp(-E/kT)dE,$ (3)

where $v$ is the relative velocity, $E$ is the energy in the centre-of-mass
system, $k$ and $m$ are the Boltzmann constant and the reduced mass of the
reacting nuclei, respectively. Thus, the reaction rate between two nuclei can
be written as $r_{12}=\frac{n_{1}n_{2}}{1+\delta_{12}}<\sigma v>$ where
$n_{1}$ and $n_{2}$ are the number densities of nuclei of types 1 and 2. The
Kronecker delta $\delta_{12}$ prevents double counting in case of identical
particles. In Figs.7-12, the plots of the (n,$\gamma$) reaction rates as
functions of T9 for Fe, Co, Ni, Cu, Zn, Ga, Ge, As and Se isotopes have been
shown.

Some key reactions Su14  which may have particularly large impact on the final
abundances in the A $\sim$ 80 region have been investigated. Several works
Fo88 ; Ma89 ; Sm93 ; An99 ; An04  regarding reaction rates have been done in
past. Except a few neutron induced reactions, all other reaction rates have
temperature dependencies. The neutron capture cross sections show $\approx
1/v$ behavior at very low energies in the thermal domain. Therefore in Eq.(3),
using $\sigma(E)\propto E^{-1/2}$ one finds immediately that at thermal
energies the reaction rates are more or less constant with respect to
temperature. However, for very low energy neutrons ($\sim$ 0.025 eV) only this
fact is true where energies are $\sim$ eV and below. On the other hand, at
energies in the region of astrophysical interest, the capture cross sections
for the neutron induced reactions can be best described by
$\sigma(E)=\frac{R(E)}{v}$ Bl55 , where $R(E)$ varies gently as a function of
energy Mu10 . It is, therefore, similar to the astrophysical S-factor and one
expects $\langle\sigma v\rangle$ to be temperature dependent. Due to the above
reason, since Bao-Käppeler fit baka  to experimental has been independent of
temperature we have excluded it and retained I. Dillman et al. (ka,kd) kakd
and KADoNiS (ks) ks  for comparison with experimental results. The comparison
of the present (n,$\gamma$) reaction rate calculations with the JINA REACLIB
reaction rates baka ; kakd ; ks  has been presented in Figs.13-21.

##  V  Summary and conclusion

To summarize, in the present work the theoretical predictions of radiative
neutron capture (n,$\gamma$) cross sections of astrophysical importance and
the reaction rates for Fe, Co, Ni, Cu, Zn, Ga, Ge, As and Se isotopes using
the Hauser-Feshbach statistical model reaction calculations have been
investigated. It is observed that the experimental results are uncertain by a
few orders of magnitude for nuclei even in the vicinity of the valley of
stability. Some key reactions Su14  which may have significantly large impact
on the final abundances in the region of mass number around eighty have been
explored. The calculations of the (n,$\gamma$) reaction rates have been
compared with the JINA REACLIB reaction rates. Since in several cases large
deviations among fits to experimental data of ka, kd and ks do exist,
estimates of present calculations can be termed as good. In addition, it is
recognized that the uncertainties due to the factors such as level densities
and mass models may have substantial effects on the rates while the low-energy
upbend in the $\gamma$-strength function has a little (though non-negligible)
effect on the rates.

To conclude, it is envisaged that to constrain the (n,$\gamma$) reaction rates
near the mass region eighty there is an acute need of more data. In order to
exclude or establish certain model inputs, new experimental techniques,
namely, the surrogate method for neutron rich nuclei Koz12  and the beta-Oslo
method Sp14  may contribute some crucial information of paramount importance.

###### Acknowledgements.

One of the authors (DNB) acknowledges support from Science and Engineering
Research Board, Department of Science and Technology, Government of India,
through Grant No. CRG/2021/007333.

## References

  * (1) Connecting Quarks with the Cosmos: Eleven Science Questions for the New Century, The National Academies Press, Washington, DC (2003).
  * (2) National Research Council, Nuclear Physics: Exploring the Heart of Matter, The National Academies Press, Washington, DC (2013).
  * (3) E. M. Burbidge, G. R. Burbidge, W. A. Fowler, F. Hoyle, Rev. Mod. Phys. 29, 547 (1957).
  * (4) A. G. W. Cameron, Chalk River Rep. (1957).
  * (5) G. Wallerstein, I. Iben Jr., P. Parker, A. M. Boesgaard, G. M. Hale, A. E. Champagne, C. A. Barnes, F. Käppeler, V. V. Smith, R. D. Hoffman, F.X. Timmes, C. Sneden, R. N. Boyd, B. S. Meyer, D. L. Lambert, Rev. Mod. Phys. 69, 995 (1997).
  * (6) F. Käppeler, R. Gallino, S. Bisterzo, W. Aoki, Rev. Mod. Phys. 83, 157 (2011).
  * (7) M. Arnould, S. Goriely, Phys. Rep. 384, 1 (2003).
  * (8) F.-K. Thielemann, I. Dillmann, K. Farouqi, T. Fischer, C. Fröhlich, A. Kelic-Heil, I. Korneev, K.-L. Kratz, K. Langanke, M. Liebendörfer, I. V. Panov, G. Martinez-Pinedo, T. Rauscher, J. Phys. Conf. Ser. 202 (1), 012006 (2010).
  * (9) Y.-Z. Qian, G. J. Wasserburg, Phys. Rep. 442, 237 (2007).
  * (10) M. Arnould, S. Goriely, K. Takahashi, Phys. Rep. 450, 97 (2007).
  * (11) F.-K. Thielemann, A. Arcones, R. Käppeli, M. Liebendörfer, T. Rauscher, C. Winteler, C. Fröhlich, I. Dillmann, T. Fischer, G. Martinez-Pinedo, K. Langanke, K. Farouqi, K.-L. Kratz, I. Panov, I. K. Korneev, Prog. Part. Nucl. Phys. 66, 346 (2011).
  * (12) E. Anders, N. Grevesse, Geochim. Cosmochim. Acta 53, 197 (1989); M. R. Mumpower et al., Progress in Particle and Nuclear Physics 86, 86 (2016).
  * (13) C. Arlandini, F. Käppeler, K. Wisshak, R. Gallino, M. Lugaro, M. Busso, O. Straniero, Astrophys. J. 525, 886 (1999).
  * (14) C. Sneden, J. J. Cowan, R. Gallino, Ann. Rev. Astron. Astrophys. 46, 241 (2008).
  * (15) I. U. Roederer, G. W. Preston, I. B. Thompson, S. A. Shectman, C. Sneden, Astrophys. J. 784, 158 (2014).
  * (16) G. J. Wasserburg, M. Busso, R. Gallino, Astrophys. J. 466, L109 (1996).
  * (17) G. J. Mathews, J. J. Cowan, Nature 345, 491 (1990).
  * (18) D. Argast, M. Samland, F.-K. Thielemann, Y.-Z. Qian, Astronom. Astrophys. 416, 997 (2004).
  * (19) Y. Komiya, S. Yamada, T. Suda, M. Y. Fujimoto, Astrophys. J. 783, 132 (2014).
  * (20) F. Matteucci, D. Romano, A. Arcones, O. Korobkin, S. Rosswog, Mon. Not. R. Astron. Soc. 447, 326 (2015).
  * (21) B. S. Meyer, G. J. Mathews, W. M. Howard, S. E. Woosley, R. D. Hoffman, Astrophys. J. 399, 656 (1992).
  * (22) S. E. Woosley, J. R. Wilson, G. J. Mathews, R. D. Hoffman, B. S. Meyer, Astrophys. J. 433, 229 (1994).
  * (23) A. Arcones, H.-T. Janka, L. Scheck, Astronom. Astrophys. 467, 1227 (2007).
  * (24) T. Fischer, S. C. Whitehouse, A. Mezzacappa, F.-K. Thielemann, M. Liebendörfer, Astronom. Astrophys. 517, A80 (2010).
  * (25) L. Hüdepohl, B. Müller, H.-T. Janka, A. Marek, G. G. Raffelt, Phys. Rev. Lett. 104, 251101 (2010).
  * (26) L. F. Roberts, S. Reddy, G. Shen, Phys. Rev. C 86, 065803 (2012).
  * (27) J. M. Lattimer, D.N. Schramm, Astrophys. J. 192 L145 (1974).
  * (28) B. S. Meyer, Astrophys. J. 343, 254 (1989).
  * (29) C. Freiburghaus, S. Rosswog, F.-K. Thielemann, Astrophys. J. 525, L121 (1999).
  * (30) S. Wanajo, Y. Sekiguchi, N. Nishimura, K. Kiuchi, K. Kyutoku, M. Shibata, Astrophys. J. 789, L39 (2014).
  * (31) O. Just, A. Bauswein, R. A. Pulpillo, S. Goriely, H.-T. Janka, Mon. Not. R. Astron. Soc. 448, 541 (2015).
  * (32) S. Goriely, A. Bauswein, H.-T. Janka, Astrophys. J. 738, L32 (2011).
  * (33) O. Korobkin, S. Rosswog, A. Arcones, C. Winteler, Mon. Not. R. Astron. Soc. 426, 1940 (2012).
  * (34) D. Wanderman, T. Piran, Mon. Not. R. Astron. Soc. 448, 3026 (2015).
  * (35) R. Surman, G. C. McLaughlin, M. Ruffert, H.-T. Janka, W. R. Hix, Astrophys. J. 679, L117 (2008).
  * (36) A. Perego, S. Rosswog, R. M. Cabezón, O. Korobkin, R. Käppeli, A. Arcones, M. Liebendörfer, Mon. Not. R. Astron. Soc. 443, 3134 (2014).
  * (37) C. Winteler, R. Käppeli, A. Perego, A. Arcones, N. Vasset, N. Nishimura, M. Liebendörfer, F.-K. Thielemann, Astrophys. J. 750 L22 (2012).
  * (38) Nobuya Nishimura, Tomoya Takiwaki, Friedrich-Karl Thielemann, Astrophys. J. 810, 109 (2015).
  * (39) T. Tsujimoto, N. Nishimura, Astrophys. J. Lett. 811, L10 (2015).
  * (40) J. Pruet, S. E. Woosley, R. D. Hoffman, Astrophys. J. 586, 1254 (2003).
  * (41) R. Surman, G. C. McLaughlin, W. R. Hix, Astrophys. J. 643, 1057 (2006).
  * (42) A. Malkus, J. P. Kneller, G. C. McLaughlin, R. Surman, Phys. Rev. D 86, 085015 (2012).
  * (43) S. Wanajo, M. Tamamura, N. Itoh, K. Nomoto, Y. Ishimaru, T. C. Beers, S. Nozawa, Astrophys. J. 593, 968 (2003).
  * (44) H.-T. Janka, B. Müller, F. S. Kitaura, R. Buras, Astronom. Astrophys. 485, 199 (2008).
  * (45) P. Banerjee, W. C. Haxton, Y.-Z. Qian, Phys. Rev. Lett. 106, 201104 (2011).
  * (46) I. U. Roederer, J. J. Cowan, A. I. Karakas, K.-L. Kratz, M. Lugaro, J. Simmerer, K. Farouqi, C. Sneden, Astrophys. J. 724, 975 (2010).
  * (47) E. G. Adelberger et al., Rev. Mod. Phys. 83, 195 (2011).
  * (48) A. J. Koning, S. Hilaire, and M. C. Duijvestijn. In Proceedings of the International Conference on Nuclear Data for Science and Technology, 22-27 April 2007, Nice, France. 2008. pp. 211-214. Arjan Koning, Stephane Hilaire and Stephane Goriely, TALYS-1.8 A nuclear reaction program, December 26, (2015).
  * (49) T. Rauscher, F.-K. Thielemann and K.-L. Kratz, Phys. Rev. C 56, 1613 (1997).
  * (50) T. Rauscher, Phys. Rev. C 81, 045807 (2010).
  * (51) T. Rauscher and F.-K. Thielemann, Atom. Data Nucl. Data Tables 75, 1 (2000).
  * (52) T. Rauscher and F.-K. Thielemann, Atom. Data Nucl. Data Tables 79, 47 (2001).
  * (53) W. Hauser and H. Feshbach, Phys. Rev. 87, 366 (1952).
  * (54) J. A. Holmes, S. E. Woosley, W. A. Fowler and B. A. Zimmerman, Atom. Data Nucl. Data Tables 18, 306 (1976).
  * (55) Amy Bartlett, Joachim Görres, Grant J. Mathews, Kaori Otsuki, Michael Wiescher, Dieter Frekers, Alberto Mengoni and Jeffrey Tostevin, Phys. Rev. C 74, 015802 (2006).
  * (56) M. Thoennessen, Reports on Progress in Physics 67, 1187 (2004).
  * (57) Rui Wang and Lie-Wen Chen, Phys. Rev. C 92, 031303(R) (2015).
  * (58) R. Boleu, S. G. Nilsson and R. K. Sheline, Physics Letters B 40, 517 (1972).
  * (59) D. D. Clayton, Principles of Stellar Evolution and Nucleosynthesis, Mc-Graw-Hill, pp. 577–91, ISBN 978-0226109534 (1968).
  * (60) W. A. Fowler and F. Hoyle, Astrophys. J. Suppl. 9, 201 (1964); Appendix C.
  * (61) D. D. Clayton, Principles of Stellar Evolution and Nucleosynthesis (University of Chicago Press, Chicago, 1983).
  * (62) J. C. Niemeyer, S. E. Woosley, Astrophys. J. 475, 740 (1997).
  * (63) P. Höflich, Nucl. Phys. A 777, 579 (2006).
  * (64) T. Strohmayer, L. Bildsten, New views of thermonuclear bursts in: W. H. G. Lewin, M. Van der Klis (Eds.), Compact Stellar X-ray Sources, (Cambridge University Press, Cambridge, London.), 113 (2006).
  * (65) H. Schatz, L. Bildsten, A. Cumming, Astrophys. J. 583, L87 (2003).
  * (66) A. Cumming, J. Macbeth, J.J.M. in t Zand, D. Page, Astrophys. J. 646, 429 (2006).
  * (67) S. Gupta, E. F. Brown, H. Schatz, P. Möller, K.-L. Kratz, Astrophys. J. 662, 1188 (2007).
  * (68) W. A. Fowler, G. R. Caughlan and B. A. Zimmerman, Ann. Rev. Astron. Astrophys. 5 (1967) 525.
  * (69) R. N. Boyd, An Introduction to Nuclear Astrophysics (University of Chicago, Chicago, 2008), 1st ed.
  * (70) R. Surman, M. Mumpower, R. Sinclair, K. L. Jones, W. R. Hix and G. C. McLaughlin, AIP Advances 4, 041008 (2014).
  * (71) G. R. Caughlan and W. A. Fowler, Atom. Data Nucl. Data Tables 40, 283 (1988).
  * (72) R. A. Malaney and W. A. Fowler, Astrophys. J. 345, L5 (1989).
  * (73) M. S. Smith, L. H. Kawano and R. A. Malaney, Astrophys. J. Suppl. 85, 219 (1993).
  * (74) C. Angulo et al., Nucl. Phys. A 656, 3 (1999).
  * (75) P. Descouvemont, A. Adahchour, C. Angulo, A. Coc and E. Vangioni-Flam, Atom. Data Nucl. Data Tables 88, 203 (2004).
  * (76) J. M. Blatt and V. F. Weisskopf, Theoretical Nuclear Physics (John Wiley $\&$ Sons, New York; Chapman $\&$ Hall Limited, London.)
  * (77) Tapan Mukhopadhyay, Joydev Lahiri and D. N. Basu, Phys. Rev. C 82, 044613 (2010); ibid Phys. Rev. C 83, 067603 (2011).
  * (78) Z. Y. Bao, F. Käppeler, Atom. Data Nucl. Data Tables 36, 411 (1987).
  * (79) I. Dillman, M. Heil, R. Käppeler, R. Plag, T. Rauscher, F.-K. Thielemann, AIP Conf. Proc. 819, 123 (2007).
  * (80) The third update of the ‘Karlsruhe Astrophysical Database of Nucleosynthesis in Stars’, I. Dillmann, R. Plag, F. Käppeler, T. Rauscher, Proceeding of the workshop ”EFNUDAT Fast Neutrons - scientific workshop on neutron measurements, theory $\&$ applications” held on April 28-30, 2009 at Geel, Belgium.
  * (81) Kozub et al., Phys. Rev. Lett. 109, 172501 (2012).
  * (82) A. Spyrou et al., Phys. Rev. Lett. 113, 232502 (2014).


# Field theory of enzyme-substrate systems with restricted long-range
interactions

Fabrizio Olmeda  Max Planck Institute for the Physics of Complex Systems,
Nöthnitzer Strasse 38, Dresden, D-01138, Germany  Institute of Science and
Technology Austria, Am Campus 1, 3400 Klosterneuberg, Austria  Steffen Rulands
<EMAIL_ADDRESS>Arnold Sommerfeld Center for Theoretical Physics and Center for
NanoScience, Department of Physics, Ludwig-Maximilians-Universität München,
Theresienstrasse 37, D-80333 Munich, Germany  Max Planck Institute for the
Physics of Complex Systems, Nöthnitzer Strasse 38, Dresden, D-01138, Germany

###### Abstract

Enzyme-substrate kinetics form the basis of many biomolecular processes. The
interplay between substrate binding and substrate geometry can give rise to
long-range interactions between enzyme binding events. Here, we study a
general model of enzyme-substrate kinetics with restricted long-range
interactions described by an exponent $-\lambda$. We employ a coherent-state
path integral and renormalization group approach to calculate the first moment
and two-point correlation function of the enzyme-binding profile. We show that
starting from an empty substrate the average occupancy follows a power law
with an exponent $1/(1-\lambda)$ over time. The correlation function decays
algebraically with two distinct spatial regimes characterized by exponents
$-\lambda$ on short distances and $-(2/3)(2-\lambda)$ on long distances. The
crossover between both regimes scales inversely with the average substrate
occupancy. Our work allows to associate experimental measurements of bound
enzyme locations with their binding kinetics and the spatial confirmation of
the substrate.

The binding of enzymes to substrates can catalyze chemical reactions. Enzyme-
substrate kinetics therefore form the basis of many biochemical processes.
Such kinetics include, for example, the binding of polymerase molecules to the
DNA for the transcription of genes to mRNA molecules [1], the binding of
oxidase enzymes to the outer membrane of mitochondria for the oxidation of
monoamines, or the multi-site phosphorylation of proteins [2].

Enzymes can interact in many ways when binding to the substrate. Transcription
factors, which control the transcription of genes, often need to bind
combinatorically together with other transcription factors in order to
initiate transcription [3]. Some enzymes seem to slide diffusively along the
substrate, such as the methyl transferase DNMT1 [4]. In many biologically
relevant scenarios the substrate undergoes confirmational changes when bound
by an enzyme [5]. This is, for example, used to increase specificity of enzyme
binding in a kinetic proof-reading scheme [6]. Many enzymes that bind to the
DNA, including chromatin modifiers, cause confirmation changes thereof, most
notably the compaction and decompaction of the DNA [7]. Positions that are far
apart measured along the DNA sequence might then be close in physical space.
Due to the feedback of enzyme binding with chromatin confirmation this give
rise to effective long-range interactions between binding events.

Experimentally, novel technologies in molecular biology allow quantifying the
binding locations of enzymes to substrates like membranes or the DNA in a
static manner. These technologies include super-resolution microscopy such as
DNA Paint [8], which can image a wide range of enzyme binding events with the
resolution of individual enzymes [9]. For one-dimensional substrates, single-
cell sequencing technologies allow measuring the consequences of enzyme
binding, such as chemical modifications of the DNA or histone tails, with the
resolution of single base pairs and in individual cells [10, 11]. A
theoretical prediction about the relation between the statistics of bound
positions on the substrate and the kinetics of enzyme binding would set the
basis for concluding the enzyme kinetics and substrate confirmation from such
experiments.

In thermal equilibrium, the binding and unbinding of enzymes is strictly
constrained by the condition of detailed balance. In this case, the
equilibrium enzyme binding profile is predicted by the Gibbs free energy [12].
In the context of cell biology, enzyme-substrate kinetics is often out-of-
equilibrium. Prominent examples of such kinetics are chromatin modifies, which
catalyze the irreversible deposition of epigenetic marks on the DNA or histone
tails and may lead to confirmational changes of the DNA. In these cases, a
description in terms of thermodynamic potentials is not feasible and the
theory of stochastic processes provides a general framework for describing
non-equilibrium enzyme-substrate kinetics [13, 14, 15, 16].

Biological substrates and enzymes are often highly nonlinear and the spatial
dimensions cannot be neglected [17, 16, 18, 19, 20, 21]. In these cases,
standard approximation schemes of master equations and Langevin equations do
not exist or are difficult to use. In these cases, path-integral
representations of stochastic processes provide a versatile and and powerful
framework for studying stochastic processes. These include the Martin-Siggia-
Rose-Janssen-de Dominicis functional integral representation of stochastic
differential equations and the coherent-state path integral representation of
master equations [22]. In both cases, expectation values of observables can be
expressed in terms of deterministic path integrals over a pair of conjugated
fields. Path-integral representations allow for the application of powerful
theoretical tools, such as perturbation theory and renormalization group
theory [23], in the context of stochastic processes.

In this work, we develop a theoretical framework of the out-of-equilibrium
enzyme-binding kinetics on substrates. Specifically, we use a coherent-state
path-integral approach in combination with renormalization-group calculations
to obtain the moments of the distribution of enzymes bound on the substrate.
We show that the first moment, the fraction of bound sites, increases
according to a power-law with an exponent depending on the decay exponent of
long-range interactions. The two-point correlation function decays in two
temporal regimes characterized by different exponents. For short distances, it
is dominated by active feedback between enzyme binding events and for long
distances by conservative noise.

Figure 1: Schematic depicting the stochastic enzyme-substrate model. In this
model, enzymes bind to a substrate with a binding rate that depends on the
distance to the nearest bound sites. In the context of cell biology,
substrates like DNA have dynamic geometries, such that local interactions in
space lead to long-range interactions along the linear lattice.

##  I Model definition

We consider a general model for enzyme-substrate kinetics. In this model, the
substrate is represented by a one-dimensional lattice of size $N$. On this
lattice, each position is a potential binding site indexed by
$i\in\\{1,\ldots,N\\}$. Enzymes can bind to and unbind from lattice sites with
rates that depend on the positions of other bound enzymes. This dependence of
the binding rates is encoded in an interaction kernel $J_{ij}$ which gives the
contribution of a bound enzyme at position $j$ to the enzyme binding rate at
position $i$.

Since in biological systems infinite-range interactions are implausible we
here restrict these interactions to the nearest bound sites. Specifically, if
the nearest bound sites are at positions $j_{1}$ and $j_{2}$ the interaction
kernel is $J_{i,j_{1}}+J_{i,j_{2}}$, where the two terms correspond to the
contribution from the left and right bound neighbour, respectively. Such
interactions are in the literature known as restricted long-range interactions
[24, 25, 26]. We further only consider translationally invariant kernels of
the form $J(|i-j|)$. The unbinding rate of the enzyme is constant and we
denote it by $u$. In equilibrium, the binding and unbinding rates would be
related by the condition of detailed balance. Here, since many enzymatic
processes are microscopically irreversible through the conversion of ATP to
ADP, we do not make this assumption and take the unbinding rate to be
independent. We further allow enzymes to perform random walks along the
substrate with a rate $\epsilon$.

With this, the state of the system is described by a random vector
$\boldsymbol{\sigma}$ with entries representing enzyme occupancy:
$\sigma_{i}=1$ if site $i$ is bound and $\sigma_{i}=0$ otherwise. The time
evolution of the probability of finding a given binding profile
$\boldsymbol{\sigma}$ at a time $t$ then follows a master equation of the form

$\displaystyle\frac{\partial P(\boldsymbol{\sigma})}{\partial t}=$
$\displaystyle\sum_{i=1}^{N}\sum_{l=i+1}^{N}J_{i,i+l}\left(\prod_{j=1}^{l-1}\bar{\sigma}_{i+j}\right)\sigma_{l}[\sigma_{i}P(\boldsymbol{\bar{\sigma}}_{i})-\bar{\sigma_{i}}P(\boldsymbol{\sigma})]+u\sum_{i=1}^{N}\sigma_{i}\left[P(\boldsymbol{\bar{\sigma}}_{i})-P(\boldsymbol{\sigma})\right]$
(1)
$\displaystyle+\epsilon\sum_{i=1}^{N}\left[\bar{\sigma}_{i}\sigma_{i+1}P(\boldsymbol{\bar{\sigma}}_{i})-\sigma_{i}\bar{\sigma}_{i+1}P(\boldsymbol{\bar{\sigma}}_{i+1})\right]+\text{l.
n. n.}+\text{l. h.}\,.$ (2)

Here, we introduced a notation for the vector $\boldsymbol{\bar{\sigma}}_{i}$
which has the same elements as $\boldsymbol{\sigma}$, but the value of
$\sigma_{i}$ at position $i$ is replaced by $1-\sigma_{i}$. Similarly, we
define $\bar{\sigma}_{i}=1-\sigma_{i}$. In Eq. (1) the first term describes
the contribution of interactions with the right nearest neighbor to the
enzyme-binding rate. The contribution from the left nearest neighbor takes the
same form and is here abbreviated by l. n. n.. The second terms describes
unbinding events and the third term enzyme hopping to to the right. Hopping
events to the left again take the same form and are abbreviated by l. h..

##  II Coherent-state path integral formulation of the master equation

Having defined a general model for the kinetics of enzyme-substrate binding in
terms of Eq. (1), we now investigate the non-stationary distribution
$P(\boldsymbol{\sigma},t)$ in terms of its moments. To this end, we employ a
coherent-state path integral formulation of the master equation, Eq. (1) [22].
Taking the semi-classical limit we will compare the prediction for the first
and higher-order moments to numerical simulations.

To define the path-integral representation we first introduce a Fock space, in
which we represent lattice configurations in bra-ket notation. We define
creation operators $a_{i}^{\dagger}$which formally represent the binding of
enzymes to site $i$. With these, we can express a given state
$\ket{\sigma}=\ket{\sigma_{1},\ldots,\sigma_{N}}$ as the repeated application
of creation operators on the empty lattice, $\ket{0}$,

$\ket{\sigma}=a_{1}^{\dagger\sigma_{1}}...a_{N}^{\dagger\sigma_{N}}\ket{0}\,.$
(3)

We also define annihilation operators, $a_{i}$, which represent the unbinding
of enzymes from a given site $i$. The annihilation operators are conjugate to
the creation operators, such that both act on a given lattice configuration as
follows,

$\begin{split}a_{i}^{\dagger}\ket{\sigma}&=\ket{\sigma_{i}+1}\,,\\\
a_{i}\ket{\sigma}&=\sigma_{i}\ket{\sigma_{i}-1}\,.\end{split}$ (4)

With these definitions, the creation and annihilation operators follow
standard commutation rules,

$[a_{i},a^{\dagger}_{i}]=1.$ (5)

We can now write the probability distribution in Fock space formally as [27],

$|P(t)\rangle=\sum_{\boldsymbol{\sigma}}P(\boldsymbol{\sigma},t)a_{1}^{\dagger\sigma_{1}}...a_{N}^{\dagger\sigma_{N}}|0\rangle\,.$
(6)

Using this notation, we can formally rewrite the master equation in terms of
the creation and annihilation operators,

$\partial_{t}|P(t)\rangle=-H|P(t)\rangle\,,$ (7)

with a ”Hamiltonian” $H$ defined as

$\displaystyle H=$
$\displaystyle\sum_{i=1}^{N}\sum_{l=1}^{N-i}J_{i,i+l}\prod_{j=1}^{l-1}\Big{[}(1-a^{\dagger}_{i+j}a_{i+j})$
(8) $\displaystyle\quad\quad\cdot
a^{\dagger}_{i+l}a_{i+l}(a^{\dagger}_{i}\hat{\delta}_{\sigma_{i,0}}-\hat{\delta}_{\sigma_{i,0}})\Big{]}$
(9) $\displaystyle+$
$\displaystyle\sum_{i=1}^{N-1}\epsilon(1-a_{i}a^{\dagger}_{i+1})\hat{\delta}_{\sigma_{i},1}\delta_{\sigma_{i+1},0}$
(10) $\displaystyle+\text{l. n. n.}+\text{l. h.}\,.$ (11)

The operators, $\hat{\delta}_{\sigma_{i},k}$ are equal to the identity if a
site $i$ is occupied by $k$ enzymes and 0 otherwise. These operators enforce
that only a single enzyme can be bound at a given site. Although we will refer
to $H$ as a Hamiltonian, $H$ does not represent an energy as it is not
necessarily hermitian unless detailed balance is fulfilled.

With this, we can write the expectation value of any observable
$A(\boldsymbol{\sigma},t)$ as

$\langle
A(\boldsymbol{\sigma},t)\rangle=\sum_{\boldsymbol{\sigma}}A(\boldsymbol{\sigma})P(\boldsymbol{\sigma},t)\,.$
(12)

We can rewrite this equation as

$\langle
A(\boldsymbol{\sigma},t)\rangle=\bra{0}\prod_{i}\sigma_{i}e^{a_{i}}A(\boldsymbol{\sigma})\ket{P(t)}\,,$
(13)

where we introduce a coherent state basis $\bra{0}e^{a}$, which has the
property to be the left eigenstate of the creation operator $a^{\dagger}$,

$\bra{0}e^{a}a^{\dagger}=\sum_{n=1}^{\infty}\frac{\bra{0}}{n!}a^{n}a^{\dagger}=\bra{0}e^{a}.$
(14)

Having defined all the rules of the operators (cf. Appendix A) we now proceed
to derive the path-integral representation of the master equation. To this
end, we define a continuous field, $\phi_{i}(t)$, giving the local density of
bound sites at position $i$. We then use the decomposition of 1,

$1=\int\text{d}\phi_{i}\,\text{d}\hat{\phi}_{i}\,e^{-\hat{\phi}_{i}\phi_{i}}e^{\phi_{i}a_{i}^{\dagger}}\ket{0}\bra{0}e^{\hat{\phi}_{i}a_{i}}\,,$
(15)

to derive a path-integral representation of the master equation following
standard steps [27, 28]. This decomposition introduces a second field,
$\hat{\phi}$, which is called response field or conjugated field. In this
path-integral representation the moments of the distribution
$P(\boldsymbol{\phi},t)$ can be expressed in terms derivatives of the
generating functional [29, 30],

$Z[\boldsymbol{h},\phi,\hat{\phi}]=\int\mathcal{D}[\boldsymbol{h},\phi,\hat{\phi}]e^{-S[\phi,\hat{\phi}]+\int\text{d}t\sum_{i}\left[h_{i}\phi_{i}+\hat{h}_{i}\hat{\phi_{i}}\right]}\,,$
(16)

where we used the notation $\mathcal{D}[\boldsymbol{h},\phi,\hat{\phi}]$ to
denote the integrals over all positions
$\prod_{i}\text{d}\phi_{i}\text{d}\hat{\phi}_{i}\text{d}h_{i}$. The argument
in the exponential is an action of the form,

$S[\hat{\phi},\phi]=-\sum_{i}\phi_{i}(t_{f})+\int_{0}^{t_{f}}\mathrm{d}t\sum_{i}\left(\hat{\phi}_{i}\partial_{t}\phi_{i}+H_{i}[\hat{\phi},\phi]\right)\,.$
(17)

with Hamiltonian,

$H_{i}=(1-\hat{\phi}_{i})e^{-\hat{\phi}_{i}\phi_{i}}\sum_{l=1}^{N-i}\prod_{j=1}^{l-1}J_{i,i+l}(l)\hat{\phi}_{i+l}\phi_{i+l}(1-\phi_{i+j})+\text{o.
t.}\,,$ (18)

where we summarized the hopping and left-nearest-neighbor interactions by o.
t..

In order to evaluate the path integral, we define how the $\hat{\delta}$
operators in $H$ act on the coherent-state basis by following the rules in
Ref. [31],

$\begin{split}&\bra{\phi}a_{i}^{\dagger}\hat{\delta}_{{n_{i}},m_{i}}\ket{\phi_{i}}=\frac{1}{m_{i}!}\hat{\phi_{i}}\left(\hat{\phi_{i}}\phi_{i}\right)^{m_{i}}e^{-\phi_{i}\hat{\phi_{i}}}\,,\\\
&\bra{\phi}a_{i}\hat{\delta}_{n_{i},m_{i}}\ket{\phi}=\frac{1}{(m_{i}-1)!}\phi_{i}\left(\hat{\phi_{i}}\phi_{i}\right)^{m-1}e^{-\phi_{i}\hat{\phi_{i}}}\,,\\\
&\bra{\phi}\hat{\delta}_{n_{i},m_{i}}\ket{\phi}=\frac{1}{m_{i}!}\left(\hat{\phi_{i}}\phi_{i}\right)^{m_{i}}e^{-\phi_{i}\hat{\phi_{i}}}\,.\end{split}$
(19)

Many different approaches have been proposed to deal with the particle
exclusion enforced by the $\hat{\delta}$ operator. Examples for this are
fermionic field theories [32] or more recent approaches using negative rates
[33]. Here, we will use the hard-bosonic path-integral approach [31]. The
benefit of using this formalism is that it gives a characteristic length
scale, which distinguishes different spatial regimes of the correlation
function and which will be important for its calculation.

###  II.1 Semiclassical limit of the field theory

In order to obtain the moments of $P(\boldsymbol{\sigma},t)$, we now derive
the semiclassical solution of the field theory. To this end, we will evaluate
the path integral in Eq. (16) for fields that make the action $S$ extremal.
This limit will be shown to be valid if the average enzyme occupancy is small.

Solutions to this model with _local_ interactions are known or can easily be
obtained using the methods we describe below. Here, we study interactions that
are motivated by the interplay between enzyme binding and geometric changes of
the substrate. For example, enzyme binding could compact the substrate, such
that sites far away on the lattice are in close proximity in real space. Such
interactions are therefore long-ranged in nature and we consider a general
class of non-local interaction kernels of the form $J_{i,i+l}=1/l^{\lambda}$.

In a first step, we rewrite the Hamiltonian, Eq. (18), in continuous space
with coordinate $s$. Upon introducing a spatial discretization $\sum_{i}\Delta
s\rightarrow\int ds$, where $\Delta s$ is the lattice spacing and $s$ a
continuous coordinate, the Hamiltonian in the action, Eq.(18) becomes,

$\begin{split}&H[\hat{\phi},\phi]=J(1-\hat{\phi}(s))e^{-\phi\hat{\phi}}\cdot\\\
&\cdot\Big{[}\int_{0}^{N-s}\mathrm{d}y\frac{\hat{\phi}(s+y)\phi(s+y)}{y^{\lambda}}e^{-\int_{z=0}^{y}\mathrm{d}z\hat{\phi}(s+z)\phi(s+z)}\\\
&+\int_{0}^{s}\mathrm{d}y\frac{\hat{\phi}(s-y)\phi(s-y)}{y^{\lambda}}e^{-\int_{z=0}^{y}\mathrm{d}z\hat{\phi}(s-z)\phi(s-z)}\Big{]}\,,\end{split}$
(20)

Here we choose the binding rates $J_{ij}$ to be equal throughout the lattice,
i.e. we do not consider disorder. As the integrals in Eq. (20) diverge for
$\lambda\geq 1$, we henceforth require that $\lambda<1$.

Expanding Eq. (20) to first order in the exponentials, to the second order in
the terms $\left.\hat{\phi}(s\pm y)\phi(s\pm y)\middle/y^{\lambda}\right.$ and
extending the upper limit of integration to infinity the terms containing the
integrals simplify to

$\left[2\left(\hat{\phi}\phi\right)^{\lambda}+\left(\phi\hat{\phi}\right)^{\lambda-3}\left(2-3\lambda+\lambda^{2}\right)\frac{\partial^{2}\left(\phi\hat{\phi}\right)}{\partial
s^{2}}\right].$ (21)

We now perform a semi-classical approximation of the generating functional by
employing a saddle-note bifurcation of the path integral. To this end, we
minimize the action, $\left.\delta
S\middle/\delta\hat{\phi(s)}\right.|_{\hat{\phi}(s)=1}=0$, while setting
$\hat{\phi}(s)=1$ to ensure probability conservation (cf. Ref. [28]). With
this, we obtain a partial differential equation describing the time evolution
of the enzyme-binding profile along the substrate, $\phi(s)$, which for small
values of $\phi(s)$ is,

$\begin{split}\frac{\partial\phi(s)}{\partial\tilde{t}}=&e^{-\phi(s)}\left\\{\phi(s)^{\lambda}+\left[\mu\phi(s)^{\lambda-3}+\epsilon\right]\partial^{2}_{s}\phi(s)\right\\}\\\
&+e^{-\phi(s)}\eta(s,t)\,.\end{split}$ (22)

Here, we rescaled time according to $\tilde{t}=t2J\Gamma(1-\lambda)$ and
defined $\mu=\left.(2-3\lambda+\lambda^{2})\middle/2\right.>0$. We reuse the
symbol $\epsilon$ for the diffusion constant. The factor $e^{-\phi(s)}$ is a
direct consequence of the site restriction, which intuitively has the role of
preventing unbounded growth. e In the hard-boson path-integral representation,
Eq. (20), the term
$\exp\left[-\int_{z=0}^{y}\mathrm{d}z\hat{\phi}(s+z)\phi(s+z)\right]$ in the
Hamiltonian is approximated by an exponential of the space-averaged fields,
$\exp\left[-y\langle\hat{\phi}(s)\phi(s)\rangle\right]$ if the fields are
slowly varying. This defines an effective exponential cutoff of the
interactions at a characteristic length $\left.1\middle/y^{\lambda}\right.$.
Because contributions from sites separated by a long distance are effectively
suppressed, this justifies the expansion of the limit in Eq. (20) to $\infty$.

The first spatially homogeneous moment, $\phi_{0}(t)$, i.e. the average number
of bound enzymes, $\langle\sigma\rangle$, is the solution of Eq. (22) when
neglecting diffusion. With an almost empty substrate as an initial condition,
the term $e^{-\phi}\phi^{\lambda}$ on the right hand side of Eq. (22) is to
first order well approximated by $\phi^{\lambda}$. Therefore, at short times,
the first moment follows a power-law,

$\langle\sigma\rangle\sim t^{\frac{1}{1-\lambda}}\,.$ (23)

This is the first result of this manuscript. In Fig. 2 we compare the
analytical solution to stochastic simulations using Gillespie’s algorithm
[34]. See Appendix B for details of the numerical implementation.

At large times, site exclusion becomes dominant and the solution approaches a
steady-state solution, $\phi^{*}$, logarithmically. The steady-state
concentration of bound enzymes is determined by a balance between binding and
unbinding events and it is given by

$\phi^{*}=(1-\lambda)W\left(\frac{\bar{u}^{\frac{1}{\lambda-1}}}{1-\lambda}\right)\,,$
(24)

where $W$ is the Lambert $W$ function and $\bar{u}=u/[2J\Gamma(1-\lambda)]$ is
the relative strength of unbinding and interactions. Therefore, if the
coupling, $J$, is much larger then $u$, the steady-state concentration of
bound enzymes decreases with $\bar{u}^{-1/(1-\lambda)}$. In the limit that
interactions are much weaker compared to degradation the steady-state
concentration in creases logaithmically,
$\phi^{*}\sim\ln\left[\bar{u}^{-1/(1-\lambda)}\right]$.

Taken together, starting from an almost empty substrate, the average number of
occupied binding sites increases following a power law up to a point where
higher-order interactions due to site-exclusion become dominant. As an
example, for enzymes binding to a fully compacted, space-filling polymer,
$\lambda$ takes a value of $1/3$ such that the number of bound sites increases
with an exponent of $3/2$.

As a remark, at the level of a mean-field description, we can map the enzyme-
substrate model to a Smoluchowski equation [35] describing the size of domains
where all sites are bound. In this view, a binding event between two already
occupied sites is equivalent to a fragmentation event splitting the domain
into two smaller intervals. Restricting our attention to the first moment, we
can rewrite our model in the mean-field limit as a fragmentation equation,

$\frac{\partial c(x)}{\partial
t}=\int_{0}^{\infty}\mathrm{d}y\,J(x,y-x)c(y)-c(x)\int_{0}^{x}\mathrm{d}y\,J(y,x-y)\,,$
(25)

where $c(x,t)$ is the number of domains of size $x$ at time $t$ and we omitted
the time-dependency for a shorter notation. For the chosen kernel,
$J(x,y)=\left.1\middle/x^{\lambda}\right.+\left.1\middle/y^{\lambda}\right.$,
and after defining the moments as
$M_{\alpha}=\int\mathrm{d}x\,c(x)x^{\alpha}$, we obtain
$M_{\alpha+\lambda+1}\sim
t^{-\left.(\alpha+\lambda)\middle/(\lambda+1)\right.}$. The average occupancy,
which is the zeroth moment, scales as $\langle\sigma\rangle\sim
t^{\left.1\middle/(\lambda-1)\right.}$. This equation gives the same time
evolution of the number of enzymes bound to the substrate as the path-integral
methods.

Figure 2: Time evolution of the average enzyme occupancy obtained from
stochastic simulations of the master equation (1) for $\lambda=1/3$ and
$10^{7}$ lattice sites (solid line). The simulation was initialized with For
early times, the theory predicts and increase following a power law with an
exponent of $3/2$ which is indicated by the dashed line.

###  II.2 Correlation functions

To compute correlation functions and scaling exponents we need to take the
derivative of the generating function (16). However, such an approach is not
feasible for three reasons.

The first reason is that we would need to expand the action and this attempt
fails due to the presence of a mass term in the field theory. Therefore, the
field theory is not scale-invariant. In order to overcome this problem, we
could try to do a change of variables in the action Eq.(20), and make a Hopf-
Cole transformation around a dynamical mean field solution, such that in that
reference frame the theory is massless. Even though the theory will be
renormalizable, the exponent we get from such a calculation does not agree
with the result of numerical simulations. The second reason why such an
approach fails is that we would need to consider perturbations of any order in
the field theory and a one-loop calculation would not be sufficient [36]. The
third and main reason for the failure of this approach is that, after
expanding around a base state, integrals of the form in Eq. (20) cannot be
approximated without losing the length scale, $1/\langle\phi\rangle$,
associated with the effective cutoff of the interactions.

The key insight to calculate the correlation function is that Eq. (20) gives
rise to two spatial regimes: at short distances, interactions are long-ranged
following a power law decay with exponent $1/\lambda$ while for distances much
larger than $1/\langle\phi\rangle$ interactions are effectively screened. In
the following, we will therefore derive the correlation function separately
for these two regimes using renormalization group methods and perturbation
theory. We will then confirm these results with numerical simulations.

####  II.2.1 Short-distance regime

To describe the short-distance behavior of correlation functions, we start
from the action, Eq. (20), and, after taking the semi-classical approximation,
we expand it to first order in the fields and their derivatives,

$\begin{split}&\partial_{t}\phi(s,t)=\int_{0}^{s}dy\,\phi(y)|s-y|^{-\lambda}e^{-\int_{z=0}^{s-y}dz\,\phi(z)}+\\\
&\int_{0}^{s}dy\,\partial_{y}\phi(y)|s-y|^{1-\lambda}e^{-\int_{z=0}^{s-y}dz\,\phi(z)}\,,\end{split}$
(26)

where $s$ is again a continuous coordinate giving the position on the
substrate. Here, for notation clarity, we omitted the noise terms, hopping
terms, and integrals of the same form describing interactions with the right
nearest bound site but we will consider them at the end of the derivation.

The interaction kernel Eq. (20) has the form
$|s-y|^{-\lambda}e^{-\int_{z=0}^{s-y}dz\,\phi(z)}$, with $\lambda<1$. By
considering a perturbation $h(s,t)$ around the mean-field solution,
$\phi_{0}(t)$, i.e. $\phi(s,t)=\phi_{0}(t)+h(s,t)$, Eq.
(LABEL:perturbationlangevin) can be expressed to first order as

$\displaystyle\partial_{t}\phi_{0}+\partial_{t}h$
$\displaystyle=e^{-\phi_{0}}\int_{0}^{s}\text{d}y\,\phi_{0}|s-y|^{-\lambda}\left[1-\int_{0}^{s-y}\text{d}z\,h(z)\right]$
$\displaystyle+e^{-\phi_{0}}\int_{0}^{s}\text{d}y\,h(y)|s-y|^{-\lambda}\left[1-\int_{0}^{s-y}\text{d}z\,h(z)\right]$
$\displaystyle+\text{h.o.t}\,.$ (27)

The first two terms of the integral on the right-hand side cancel with the
first one on the left-hand side, which follows from the dynamical mean field
solution, Eq. (22). Making a change of variables, $w=z+y$, we then obtain

$\begin{split}\partial_{t}h=&e^{-\phi_{0}}\int_{0}^{s}\text{d}y\,h(y)|s-y|^{-\lambda}\\\
-&e^{-\phi_{0}}\int_{0}^{s}\text{d}y\int_{y}^{s}\text{d}w\,h(y)|s-y|^{-\lambda}h(w-y)+\xi(s,t)\,.\end{split}{}$
(28)

The second term on the right-hand side of Eq. (28) is the convolution of a
fractional integral of a function and the function itself.

The noise term, $\xi(s,t)$, is Gaussian white noise. Its correlations can be
derived from the field theory, Eq. (20), by identifying terms that are
proportional to $\hat{\phi}^{2}$. In the expansion of the action these terms
comprise both non-conservative noise, which means that they are proportional
to $\hat{\phi}^{2}$, and conservative noise, which means that they are
proportional to $\hat{\phi}^{2}\partial_{s}^{2}\phi$. The perturbation
$h(s,t)$ thus has both conservative and non-conservative components, such that
we make the following ansatz for the noise correlations:

$\langle\xi(s,t)\xi(s^{\prime},t^{\prime})\rangle=\delta(t-t^{\prime})(2\Gamma_{NC}-2\Gamma_{C}\partial_{s}^{2})\delta(s-s^{\prime}).$
(29)

$\Gamma_{C}$ and $\Gamma_{NC}$ are the noise strengths for conservative and
non-conservative noise, respectively, which are functions of $J$ only. We omit
their dependency on the model parameters as this dependence does not influence
the exponents of correlation functions.

Taken together, as the fractional integral scales in Fourier space as
$q^{\lambda-1}$, a spatio-temporal perturbation in Fourier space follows a
nonlinear Langevin equation of the form

$\begin{split}\partial_{t}h(q,t)=&\left(e^{-\phi_{0}(t)}q^{\lambda-1}-q^{2}\right)h(q,t)\\\
&-q^{\lambda-1}e^{-\phi_{0}(t)}h(q,t)^{2}+\xi(q,t)\,.\end{split}$ (30)

From the dynamical mean-field solution of the first moment, Eq. (22), we
obtain, up to a prefactor, that
$e^{-\phi_{0}}=e^{\left[{-t^{\left.1\middle/(1-\lambda)\right.}}\right]}$. In
the frequency domain, we hence obtain at small times, $t\rightarrow 0$ or
$\omega\to\infty$, that

$G_{0}^{-1}h(q,\omega)=-q^{\lambda-1}h(q,\omega)^{2}+\xi(q,\omega)\,,$ (31)

where we have defined the inverse bare propagator as
$G_{0}^{-1}=i\omega+q^{2}-q^{\lambda-1}$. From the above equation, the bare
correlations of $h$ are then

$C_{0}=\left(2\Gamma_{NC}+2\Gamma_{C}q^{2}\right)|G_{0}|^{2}\,.$ (32)

Eq. (31) is of second order and so admits an exact solution from which we
calculate the two-point correlation function, $\langle
h(q,\omega)h(q^{\prime},\omega^{\prime})\rangle$, where the average
$\langle\ldots\rangle$ is performed over the noise.

As we are considering the short-distance regime, we keep only leading-order in
$q$ and perform the inverse Fourier transform. With this, we find that
correlation functions asymptotically approach a power law,

$\langle
h(s,t)h(s^{\prime},t)\rangle=\frac{\text{cos}(\pi\lambda)\left(\absolutevalue{s-s^{\prime}}^{2}\Gamma_{C}+\Gamma_{NC}\lambda(1+\lambda)\right)}{\Gamma(\lambda)\absolutevalue{s-s^{\prime}}^{2+\lambda}}\,.$
(33)

We therefore obtain for the scaling of the correlation function that

$\langle
h(s,t)h(s^{\prime},t)\rangle\sim\absolutevalue{s-s^{\prime}}^{-\lambda}\,,$
(34)

as the conservative noise ($\sim q^{2}$) is leading with respect to non-
conservative noise ($\sim q^{0}$) for short-distances. This result is
intuitive as the spatial correlations are proportional to the exponent of the
scale-free kernel $J$. This derivation is valid for low values of the average
occupancy $\phi_{0}$. We expect that for larger values of $\phi_{0}$ higher-
order corrections become relevant.

####  II.2.2 Long-distance regime

From the bare propagator and correlator in Eq (32) we notice that, at large
distances, second-order spatial derivatives dominate correlations, and due to
non-conservative noise we expect correlation functions to be described by
different exponents compared to the short-distance regime, where conservative
noise was the leading contribution in the scaling, Eq. (33). The mean-field
exponents can be estimated by dimensional analysis [28], where we compute the
scaling of the parameters in the action in order to make it dimensionless.
Following these arguments, we find from Eq. (20) that the value of the
critical exponent of the correlation in the long-distance regime should be
$\chi=\left.-(1+d+\lambda)\middle/3\right.$. This follows from considering
higher-order non-linearities and it is the correct exponent for the long tail
of the correlation function. In the following, we will derive the results of
this scaling argument rigorously with the use of dynamical renormalization
group methods.

We begin with Eq. (LABEL:perturbationlangevin). After linearization it
becomes,

$\displaystyle{}\partial_{t}h(s)$
$\displaystyle=\partial_{s}^{2}h(s)+\int_{0}^{s}\text{d}yh(y)|s-y|^{-\lambda}$
(35)
$\displaystyle-\int_{0}^{s}\text{d}yh(y)|s-y|^{1-\lambda}\left[h(s)-(s-y)\frac{1}{2}\partial_{s}h(s)\right]$
(36) $\displaystyle+\xi(s,t)\,,$ (37)

where we used that $\int_{a}^{b}\text{d}s\,f(s)\approx(b-a)(f(a)+f(b))/2$.
Before proceeding with a renormalization procedure, we note that the non-
linearities of order $h^{2}$ in Eq. (35) can be relevant in this regime.

Figure 3: (a) Overview of diagrammatical elements and diagrams contributing to
the renormalization of the (b) propagator, (c) correlator and (d) vertex
functions. Every intersection corresponds to an integration over the wave-
vectors.

Because we approximated to linear order when obtaining Eq. (35) we must keep
all the other quadratic terms in the field theory, Eq. (20). There is only one
quadratic term which is,

$\left(\hat{\phi}-1\right)\phi\hat{\phi}\int_{0}^{s}\mathrm{d}y\frac{\hat{\phi}(s-y)\phi(s-y)}{y^{\lambda}}e^{-\int_{0}^{y}\mathrm{d}z\hat{\phi}(s-z)\phi(s-z)}.$
(38)

However, after functional minimization of the action it cancels out with the
term involving interactions with the right nearest neighbor. This is not
surprising, because the model does not allow for field theories involving
terms that break the space-reversal symmetry $s\rightarrow-s$. Taken together,
by repeating the same calculations as before and including the additional term
in Eq. (38) we obtain

$\displaystyle\partial_{t}h(s)$
$\displaystyle=\partial_{s}^{2}h(s)+\int_{0}^{s}\text{d}y\,h(y)|s-y|^{-\lambda}$
(39)
$\displaystyle+\frac{1}{2}\int_{0}^{s}\text{d}y\,h(y)|s-y|^{2-\lambda}\partial_{s}h(s)+\xi(s,t)\,.$
(40)

Considering both right and left nearest-neighbor interactions, the advective
terms cancel out, again by the necessity to obey left-right symmetry. We then
include the next highest-order term,

$\begin{split}\partial_{t}h(s)=&\partial_{s}^{2}h(s)+\int_{0}^{x}\text{d}y\,h(y)|s-y|^{-\lambda}\\\
&+\frac{1}{2}\partial_{s}^{2}h(s)\int_{0}^{s}\text{d}y\,h(y)|s-y|^{2-\lambda}+\xi(s)\,.\end{split}$
(41)

The following results are applicable to multi-dimensional systems. We
therefore now generalize to any spatial dimension by writing spatial
coordinate in vector form, $\boldsymbol{s}$. In Fourier space, Eq. (41) can
then be written in compact form,

$\displaystyle G_{0}(\boldsymbol{q})^{-1}h(\boldsymbol{q},\omega)=$
$\displaystyle-\mu\int_{\boldsymbol{k}}W(\boldsymbol{q},\boldsymbol{k})h(\boldsymbol{q},\omega)h(\boldsymbol{k-q},\omega)$
(42) $\displaystyle+\xi(\boldsymbol{q},\omega)\,,$ (43)

where $h(\boldsymbol{q},t)$ is the Fourier transform of $h(s,t)$ and
$\mu=\left.(2-3\lambda+\lambda^{2})\middle/2\right.$,

$h(\boldsymbol{q},\omega)=\int\text{d}\boldsymbol{s}\int\text{d}t\,h\left(\boldsymbol{s},t\right)e^{i\boldsymbol{q}s}e^{i\omega
t}\,.$ (44)

We also defined the free propagator,

$G_{0}^{-1}=\left(i\omega+D\boldsymbol{q}^{2}+J|\boldsymbol{q}|^{-\lambda}\right)\,,$
(45)

in which we reintroduced the dimensional parameters from the adimensional Eq.
(42) in order to study their flow under renormalization. Finally, we defined
the vertex, which accounts for the nonlinear terms, as

$W(\boldsymbol{q},\boldsymbol{k})=\frac{1}{2}\left[\frac{\boldsymbol{q}(\boldsymbol{k-q})}{|\boldsymbol{k-q}|^{3-\lambda}}+\frac{(\boldsymbol{k-q})\boldsymbol{q}}{|\boldsymbol{q}|^{3-\lambda}}\right]\,.$
(46)

In the hydrodynamic limit, $\absolutevalue{\boldsymbol{k}}\rightarrow 0$ and
$\absolutevalue{\boldsymbol{q}}\to 0$, the vertex scales linearly with
$\absolutevalue{\boldsymbol{k}}$ and $\absolutevalue{\boldsymbol{q}}$, which
implies non-renormalization of the vertex function.

The rescaling step of the renormalization group procedure then gives a
rescaling of the form

$\begin{split}\partial_{l}\epsilon&=\left[z-2+A_{D}\right]\epsilon\,,\\\
\partial_{l}\mu&=\left[z+\chi-2+(3-\lambda)\right]\mu\,,\\\
\partial_{l}\Gamma_{C}&=\left[z-2\chi-d-2+A_{\Gamma_{C}}\right]\Gamma_{C}\,,\\\
\partial_{l}\Gamma_{NC}&=\left[z-2\chi-d+A_{\Gamma_{NC}}\right]\Gamma_{NC}\,.\end{split}$
(47)

$A_{\Gamma_{NC}}$ and $A_{\Gamma_{C}}$ depend on all of the parameters and
need to be evaluated perturbatively by calculating the diagrams in Fig. 3.

Following standard calculations for the integrals represented by the diagrams
in Fig. 3 [23], we recover the renormalization-group flow for the parameters,

$\displaystyle\partial_{l}\ln(\epsilon)=$ $\displaystyle z-2$ (48)
$\displaystyle-\frac{K_{d}\mu^{2}}{d\epsilon^{3}}\left[\left(d-2\right)\Gamma_{NC}+\left(d-3\right)\Gamma_{C}\right],$
(49) $\displaystyle\partial_{l}\ln(\mu)=$ $\displaystyle
z+\chi-2+(3-\lambda),$ (50) $\displaystyle\partial_{l}\ln(\Gamma_{C})=$
$\displaystyle z-2\chi-d-2$ (51)
$\displaystyle-\frac{K_{d}\nu^{2}(1+d)(\Gamma_{NC}+\Gamma_{C})^{2}}{2d\sigma_{0}^{3}\Gamma_{C}},$
(52) $\displaystyle\partial_{l}\ln(\Gamma_{NC})=$ $\displaystyle z-2\chi-d\,,$
(53)

where $K_{d}=S_{d}/(2\pi)^{d}$ and $S_{d}$ is the surface area of a
$d$-dimensional sphere.

From the non-renormalization of the non-conserved noise and of the couplings
$\mu$ we get the exact exponent identities
$\chi=\left.(-1-d+\lambda)\middle/3\right.$ and
$z=\left.(-2+d+2\lambda)\middle/3\right.$, where $z$ is the dynamical critical
exponent. In $d=1$, two-point spatial correlations functions then decay with
an exponent $\left.2(\lambda-2)\middle/3\right.$ in the long-distance regime.

Taken together, we find that the correlation function along one-dimensional
substrates decays in two algebraic regimes [Fig. 4(a)],

$\displaystyle\langle
h(s)h(s^{\prime})\rangle=\begin{cases}|s-s^{\prime}|^{-\lambda}.&\text{for
}|s-s|\ll 1/\langle\sigma\rangle\,,\\\
|s-s^{\prime}|^{-\frac{2}{3}(2-\lambda)},&\text{for }|s-s|\gg
1/\langle\sigma\rangle\,.\end{cases}$ (54)

The cross-over between these regimes stems from an effective exponential
cutoff of the long-range interactions. The position of the cross-over scales
with the only length scale in the system, the typical distance between
neighboring occupied sites, $1/\langle\sigma\rangle$, and, intuitively,
separates a regime dominated by long-range interactions and a regime
characterized by passive, conservative fluctuations.

Figure 3(b) shows correlation functions obtained from numerical simulation of
Eq. (1) for the special case of a space-filling substrate ($\lambda=1/3$). In
order to separately emphasize the short-distance and the long-distance regime
we computed correlations at two levels of the average substrate occupancy.
While we cannot confidently estimate the numerical exponent for the short-
distance regime, the simulation data confirm the existence of a crossover and
the exponent of $-10/9$ in the long-distance regime.

Figure 4: (a) The two-point correlation function decays in two spatial regimes
with exponents given by Eq. (54). The position of the cross-over between both
regimes increases with the typical length scale associated with the enzyme
density, which scales as $1/\langle\phi\rangle$. (b) Numerical simulation of
the spatial correlation functions for $\lambda=1/3$ from stochastic
simulations using Gillespie’s algorithm with $10^{7}$ lattice sites. We
calculated the correlation function in two stages of the simulation: for an
average occupancy probability of 0.006 (grey, emphasizing the short-distance
regime) and an average occupancy of 0.108 (black, emphasizing the long-
distance regime). The vertical lines indicate the corresponding crossover
positions. Dashed lines indicate the predicted exponents in both regimes. We
binned the data logarithmically and error bars indicate mean $\pm$ standard
error.

##  III Discussion

In summary, we studied a stochastic enzyme-substrate model where binding
events are correlated via long-range interactions. Such interactions mimic the
effect of conformational changes in the substrate, where positions far apart
along the substrate might be close in physical space. We employed a coherent-
state path-integral representation of the master equation and renormalization-
group theory to calculate the exponents describing the time evolution of the
average occupancy and the

In polymer systems, the exponent $\lambda$ relates to the end-to-end distance
of the polymer. For example, for a space-filling polymer, where $\lambda=1/3$,
the average occupancy increases with $t^{3/2}$ while the correlation function
decays with $|x-x^{\prime}|^{-1/3}$ for short distances and
$|x-x^{\prime}|^{-10/9}$ for long distances. Our results were obtained in a
limit, where two-body interactions dominate. If substrate occupation is high,
we expect that these exponents might be modified to incorporate higher-order
interactions.

Novel technologies in super-resolution microscopy can quantify the location of
individual enzymes bound to membranes or the DNA [9]. Such technologies allow
for the quantification of spatial correlation functions of bound enzymes. Our
work allows associating the spatio-temporal statistics of bound enzymes to the
enzyme-substrate kinetics and spatial confirmation of the substrate. Measuring
the moments along the substrate in an experiment would therefore allow drawing
conclusions about the underlying biochemical processes.

##  IV Acknowledgments

We thank F. Piazza , M. Henkel, and F. Jülicher for helpful feedback and the
entire Rulands group for fruitful discussions. We thank W. Reik, S. Clark, T.
Lohoff, and I. Kafetzopoulos for fruitful discussions about the biological
aspects of this work. This project has received funding from the European
Research Council (ERC) under the European Union’s Horizon 2020 research and
innovation program (grant agreement no. 950349). This project has received
funding from the European Union’s Horizon 2020 research and innovation
programme under the Marie Skłodowska-Curie grant agreement No 101034413.

##  Appendix A Path Integral Formulation of the Master Equation

To construct the field theory we first note that we can formally write a
solution of the master equation as

$|P(t)\rangle=e^{-Ht}|P(0)\rangle\,,$ (55)

where $|P(0)\rangle$ is the initial state (i.e. the probability distribution
of enzyme binding profiles at time $t=0$). The exponential factor can be
rewritten as

$e^{-Ht}=(1-\Delta tH)^{\frac{t}{\Delta t}}=(1-\Delta tH)\cdot(1-\Delta
tH)\cdot\ldots\,.$ (56)

By inserting the identity in the coherent state basis between every factor on
the right-hand side of Eq. (55), we find that the solution for any time
$t_{1}$ can be written as

$\begin{split}&|P(t_{1})\rangle=\int\prod\sigma_{i}d\phi_{i}(t_{1}+\Delta
t)d\hat{\phi}_{i}(t_{1}+\Delta t)d\phi(t_{1})d\hat{\phi}_{i}(t_{1})\\\ &\cdot
e^{-\hat{\phi}_{i}(t_{1})\phi(t_{1})}\cdot e^{-\hat{\phi}_{i}(t_{1}+\Delta
t)\phi(t_{1}+\Delta t)_{i}}\cdot e^{\phi_{i}(t_{1}+\Delta t)a^{\dagger}}\\\
&\cdot\ket{0}\bra{0}e^{\hat{\phi}_{i}(t_{1}+\Delta t)a}\left(1-\Delta
tH\right)e^{\phi_{i}(t_{1})a^{\dagger}}\ket{0}\bra{0}e^{\hat{\phi}_{i}(t_{1})a}\,.\end{split}$
(57)

In this equation, we have to evaluate quantities in the coherent state basis
between the bra and the ket,

$\begin{split}&\bra{0}e^{\hat{\phi_{i}}(t_{1}+\Delta t)a}\left(1-\Delta
tH\right)e^{\phi_{i}(t_{1})a^{\dagger}}\ket{0}=\\\
&e^{\hat{\phi_{i}}(t_{1}+\Delta t)}e^{\phi_{i}(t_{1})}-\Delta
t\bra{0}e^{\hat{\phi_{i}}(t_{1}+\Delta
t)a}\left(H\right)e^{\phi_{i}(t_{1})a^{\dagger}}\ket{0}\approx\\\
&e^{\hat{\phi_{i}}(t_{1}+\Delta t)}e^{\phi_{i}(t_{1})}e^{-\Delta
tH(\hat{\phi_{i}}(t_{1}),\phi_{i}(t_{1}))}\,,\end{split}$ (58)

where $H(\hat{\phi}_{i},\phi_{i})$ is obtained by replacing all $a_{i}$ with
$\phi_{i}$ and $a^{\dagger}_{i}$ with $\hat{\phi}_{i}$, using the previously
introduced decomposition and by taking into account the rules for exclusion
processes, Eq.(19). Repeating this procedure $\left.t\middle/\Delta t\right.$
times for each factor $(1-\Delta tH)$ we end up with an integral, over a
product of three terms $P_{2}$, $P_{3}$, $P_{4}$. The integral is,

$\int\prod\sigma_{i}d\hat{\phi}_{i}(t)d\phi_{i}(t)\ldots
d\hat{\phi}_{i}(\Delta t)d\phi_{i}(\Delta t)d\hat{\phi}_{i}(0)d\phi_{i}(0)\,,$
(59)

which can be rewritten compactly as a functional integral
$\int\mathcal{D}[\phi]\mathcal{D}[\hat{\phi}]$. $P_{2}$ is composed of a
product of terms that can be rewritten utilizing Riemann integration,

$P_{2}=\prod_{t_{1}=\Delta t}^{t}e^{\hat{\phi}(t_{1}+\Delta
t)\phi(t_{1})-\hat{\phi}(t_{1})\phi(t_{1})}\approx e^{-\int
dt\partial_{t}\hat{\phi}\phi}.$ (60)

Finally, there are further $\left.t\middle/\Delta t\right.$ terms coming from
the Hamiltonian evaluated at each time step,

$P_{3}=\prod_{t_{1}=\Delta t}^{t}e^{-\Delta
tH(\hat{\phi}(t_{1}),\phi(t_{1}))}\approx e^{-\int
dtH(\hat{\phi}(t),\phi(t))}\,.$ (61)

The final factor, $P_{4}$, represents initial conditions and we refer to Ref.
[27] for a discussion of this term. With this, any observable can then be
expressed as a path integral of the form

$A(\boldsymbol{\sigma})=\int\mathcal{D}[\phi]\mathcal{D}[\hat{\phi}]A(\boldsymbol{\phi},\boldsymbol{\hat{\phi}}=1)e^{-S[\hat{\phi},\phi]}\,,$
(62)

with

$S[\hat{\phi},\phi]=-\sum_{i}\phi_{i}(t_{f})+\int_{0}^{t_{f}}\mathrm{d}t\sum_{i}\left(\hat{\phi}_{i}(t)\partial_{t}\phi_{i}(t)+H_{i}[\hat{\phi},\phi]\right)\,,$
(63)

where we have performed a partial integration in time and,

$\begin{split}&H_{i}[\hat{\phi},\phi]=(1-\hat{\phi}_{i})e^{-\hat{\phi}_{i}\phi_{i}}\\\
&\left(\sum_{l=1}^{N-i}\prod_{j=1}^{l-1}J_{i,i+l}(l)\hat{\phi}_{i+l}\phi_{i+l}(1-\phi_{i+j})++\textit{o.t.}\right)\,.\end{split}$
(64)

Defining the generating functional of correlations,
$Z[\boldsymbol{h},\phi,\hat{\phi}]$, as

$Z[\boldsymbol{h},\phi,\hat{\phi}]=\int\mathcal{D}[\boldsymbol{h},\phi,\hat{\phi}]e^{-S[\phi,\hat{\phi}]+\int\int
dsdt\left[h\phi+\hat{h}\hat{\phi}\right]}\,,$ (65)

expectation values of products of observables, such as correlation functions,
can then be expressed as functional derivatives with respect to the auxiliary
external field,

$\langle\phi(s,t)\phi(y,t^{\prime})\rangle=\frac{\delta^{2}}{\delta
h(s,t)\delta h(y,t)}Z[\boldsymbol{h},\phi,\hat{\phi}]|_{\boldsymbol{h}=0}\,.$
(66)

##  Appendix B Stochastic simulations

To test the validity of our analytical results we perform extensive stochastic
simulations by integration of the master equation Eq. (1) using Gillespie’s
algorithm [34] with $1e7$ lattice sites, a random initial distribution of
enzymes with a occupancy fraction equal to $10^{-4}$, and interaction strength
$J=1$. Unbinding and hopping rates are set to zero as they do not affect the
critical dynamics and would lead to significantly slower simulations. This is
because the number of lattice sites needs to be large enough to measure the
exponents of the spatial correlations.

## References

  * Alberts  _et al._ [2002] B. Alberts, D. Bray, J. Lewis, M. Raff, K. Roberts, and J. Watson,  _Molecular Biology of the Cell_ , 4th ed. (Garland, 2002).
  * Salazar and Höfer [2009] C. Salazar and T. Höfer, Multisite protein phosphorylation – from molecular mechanisms to kinetic models, The FEBS Journal 276, 3177 (2009).
  * Reiter  _et al._ [2017] F. Reiter, S. Wienerroither, and A. Stark, Combinatorial function of transcription factors and cofactors., Current Opinion in Genetics and Development 43, 73 (2017).
  * Vilkaitis  _et al._ [2005] G. Vilkaitis, I. Suetake, S. Klimašauskas, and S. Tajima, Processive methylation of hemimethylated cpg sites by mouse dnmt1 dna methyltransferase*, Journal of Biological Chemistry 280, 64 (2005).
  * Koshland [1963] D. E. J. Koshland, Correlation of structure and function in enzyme action, Science 142, 1533 (1963).
  * Savir and Tlusty [2007] Y. Savir and T. Tlusty, Conformational proofreading: The impact of conformational changes on the specificity of molecular recognition, PLOS ONE 2, 1 (2007).
  * Dekker and Misteli [2015] J. Dekker and T. Misteli, Long-range chromatin interactions., Cold Spring Harb Perspect Biol 7, a019356 (2015).
  * Schnitzbauer  _et al._ [2017] J. Schnitzbauer, M. T. Strauss, T. Schlichthaerle, F. Schueder, and R. Jungmann, Super-resolution microscopy with dna-paint, Nature Protocols 12, 1198 (2017).
  * Reinhardt  _et al._ [2023] S. C. M. Reinhardt, L. A. Masullo, I. Baudrexel, P. R. Steen, R. Kowalewski, A. S. Eklund, S. Strauss, E. M. Unterauer, T. Schlichthaerle, M. T. Strauss, C. Klein, and R. Jungmann, Ångström-resolution fluorescence microscopy, Nature 617, 711 (2023).
  * Clark  _et al._ [2018] S. J. Clark, R. Argelaguet, C.-A. Kapourani, T. M. Stubbs, H. J. Lee, C. Alda-Catalinas, F. Krueger, G. Sanguinetti, G. Kelsey, J. C. Marioni, O. Stegle, and W. Reik, scnmt-seq enables joint profiling of chromatin accessibility dna methylation and transcription in single cells, Nature Communications 9, 781 (2018).
  * Gawad  _et al._ [2016] C. Gawad, W. Koh, and S. R. Quake, Single-cell genome sequencing: current state of the science, Nature Reviews Genetics 17, 175 (2016).
  * Langmuir [1918] I. Langmuir, The adsorption of gases on plane surfaces of glass, mica and pt, Journal of the American Chemical Society 40, 1361 (1918).
  * Gardiner [2004] C. W. Gardiner,  _Handbook of stochastic methods for physics, chemistry and the natural sciences_ , 3rd ed., Springer Series in Synergetics, Vol. 13 (Springer-Verlag, Berlin, 2004) pp. xviii+415.
  * Reuveni  _et al._ [2014] S. Reuveni, M. Urbakh, and J. Klafter, The role of substrate unbinding in michaelis-menten enzymatic reactions, Biophysical Journal 106, 677a (2014).
  * Morelli  _et al._ [2011] M. J. Morelli, R. J. Allen, and P. R. t. Wolde, Effects of macromolecular crowding on genetic networks., Biophys J 101, 2882 (2011).
  * Jülicher and Bruinsma [1998] F. Jülicher and R. Bruinsma, Motion of rna polymerase along dna: A stochastic model, Biophysical Journal 74, 1169 (1998).
  * Li  _et al._ [2009] G.-W. Li, O. G. Berg, and J. Elf, Effects of macromolecular crowding and dna looping on gene regulation kinetics, Nature Physics 5, 294 (2009).
  * Coppey  _et al._ [2004] M. Coppey, O. Bénichou, R. Voituriez, and M. Moreau, Kinetics of target site localization of a protein on dna: A stochastic approach, Biophysical Journal 87, 1640 (2004).
  * Mirny  _et al._ [2009] L. Mirny, M. Slutsky, Z. Wunderlich, A. Tafvizi, J. Leith, and A. Kosmrlj, How a protein searches for its site on dna: the mechanism of facilitated diffusion, Journal of Physics A: Mathematical and Theoretical 42, 434013 (2009).
  * Buchner  _et al._ [2013] A. Buchner, F. Tostevin, and U. Gerland, Clustering and optimal arrangement of enzymes in reaction-diffusion systems, Physical Review Letters 110, 208104 (2013).
  * Gupta  _et al._ [2017] A. Gupta, A. Milias-Argeitis, and M. Khammash, Dynamic disorder in simple enzymatic reactions induces stochastic amplification of substrate, Journal of The Royal Society Interface 14, 20170311 (2017), https://royalsocietypublishing.org/doi/pdf/10.1098/rsif.2017.0311 .
  * Altland and Simons [2010] A. Altland and B. D. Simons,  _Condensed matter field theory_ (Cambridge University Press, 2010).
  * Zinn-Justin [2002] J. Zinn-Justin,  _Quantum Field Theory and Critical Phenomena; 4th ed._ , International series of monographs on physics (Clarendon Press, Oxford, 2002).
  * Ginelli  _et al._ [2005] F. Ginelli, H. Hinrichsen, R. Livi, D. Mukamel, and A. Politi, Directed percolation with long-range interactions: Modeling nonequilibrium wetting, Physical Review E 71, 026121 (2005).
  * Ginelli  _et al._ [2006] F. Ginelli, H. Hinrichsen, R. Livi, D. Mukamel, and A. Torcini, Contact processes with long range interactions, Journal of Statistical Mechanics: Theory and Experiment 2006, P08008 (2006).
  * Hinrichsen [2007] H. Hinrichsen, Non-equilibrium phase transitions with long-range interactions, Journal of Statistical Mechanics: Theory and Experiment 2007, P07006 (2007).
  * Cardy  _et al._ [2008] J. Cardy, G. Falkovich, and K. Gawedzki,  _Non-equilibrium Statistical Mechanics and Turbulence_ , edited by S. Nazarenko and O. V. Zaboronski, London Mathematical Society Lecture Note Series (Cambridge University Press, 2008).
  * Täuber [2014] U. C. Täuber,  _Critical Dynamics: A Field Theory Approach to Equilibrium and Non-Equilibrium Scaling Behavior_ (Cambridge University Press, 2014).
  * Doi [1976] M. Doi, Second quantization representation for classical many-particle system, Journal of Physics A: Mathematical and General 9, 1465 (1976).
  * Peliti [1985] L. Peliti, Path integral approach to birth-death processes on a lattice, Journal de Physique 46, 1469 (1985).
  * Wijland [2001] F. v. Wijland, Field theory for reaction-diffusion processes with hard-core particles, Physical Review E 63, 022101 (2001).
  * Brunel  _et al._ [2000] V. Brunel, K. Oerding, and F. van Wijland, Fermionic field theory for directed percolation in (1+ 1)-dimensions, Journal of Physics A: Mathematical and General 33, 1085 (2000).
  * Nekovar and Pruessner [2016] S. Nekovar and G. Pruessner, A field-theoretic approach to the wiener sausage, Journal of Statistical Physics 163, 604 (2016).
  * Gillespie [1976] D. T. Gillespie, A general method for numerically simulating the stochastic time evolution of coupled chemical reactions, Journal of Computational Physics 22, 403 (1976).
  * Krapivsky  _et al._ [2010] P. L. Krapivsky, S. Redner, and E. Ben-Naim,  _A kinetic view of statistical physics_ (Cambridge University Press, 2010).
  * Janssen  _et al._ [2004] H.-K. Janssen, F. van Wijland, O. Deloubrière, and U. C. Täuber, Pair contact process with diffusion: Failure of master equation field theory, Physical Review E 70, 056114 (2004).