The conformation of a semiflexible filament in a quenched random potential

TL;DR

This paper investigates the conformational behavior of a semiflexible filament in a quenched random potential, combining theoretical analysis with simulations to understand how cross-linking influences filament shape and energy storage.

Contribution

It introduces a theoretical framework for the conformational statistics of semiflexible filaments in random potentials and compares predictions with Brownian dynamics simulations.

Findings

Curvature distribution depends on potential correlation length and strength.

Theoretical predictions align with simulation results.

Replica field techniques' applicability is discussed.

Abstract

Motivated by the observation of the storage of excess elastic free energy - (prestress) -- in cross linked semiflexible networks, we consider the problem of the conformational statistics of a single semiflexible polymer in a quenched random potential. The random potential, which represents the effect of cross linking to other filaments is assumed to have a finite correlation length $\xi$ and mean strength $V_{0}$. We examine statistical distribution of curvature in filament with thermal persistence length $\ell_{P}$ and length $L_0$ in the limit that $\ell_{P} \gg L_0$. We compare our theoretical predictions to finite element Brownian dynamics simulations. Lastly we comment on the validity of replica field techniques in addressing these questions.

Figures

Equations

$$E[y(x)]={\int }_0^{L_0}dx\left\{\frac{1}{2}\kappa {\ddot{y}}^2+\frac{1}{2}\tau {\dot{y}}^2+V\left(x,y\left(x\right)\right)\right\},$$

$$Z=\int \mathcal{D}y(x)e^{-\beta E[y(x)]}.$$

$${\mathcal{P}}_{\mathrm{V}}(V)=\frac{1}{P_0}\exp \left\{-\frac{1}{2V_0^2}\int d^{2}x\left[(\nabla V{)}^2+{\xi }^{-2}{V}^2\right]\right\}.$$

$${\mathcal{P}}_{\mathrm{F}}(V)\sim \exp \left[-\frac{1}{2V_0^2}\int d^{2}x{\xi }^2({\nabla }^{2}V{)}^2+(\nabla V{)}^2\right].$$

$$P(V_{k_x,k_y})\propto \exp \left\{-\frac{L_x{L}_y}{{\xi }^2}\frac{V_{k_x,k_y}^2}{8V_0^2}\exp \left({\xi }^2(k_x^2+k_y^2)\right)\right\},$$

$$V(x,y)=\frac{\tilde{k}}{2}{\left[y-y_0(x)\right]}^2.$$

$$\frac{E_{\mathrm{flexible}}}{T}={\int }_0^{L_0}dx\left\{y(x)\mathcal{O}y(x)+ky(x)y_0(x)+\frac{k{y}_0^2(x)}{2}\right\},$$

$$\mathcal{O}=\frac{m}{2}{\partial }^2+\frac{k}{2}.$$

$$Z_{\mathrm{flexible}}=\int \mathcal{D}y{e}^{-E_{\mathrm{flexible}}[y(x)]}.$$

$$Z=\frac{1}{\sqrt{\det \mathcal{O}}}\exp \left[-{\int }_0^{L_0}dx\left(-\frac{k^2}{4}{y}_0(x){\mathcal{O}}^{-1}{y}_0(x)+\frac{k{y}_0^2}{2}\right)\right].$$

$$\frac{Z}{Z_{\mathrm{SHO}}}=\exp \left[-{\int }_0^{L_0}dx\left(-\frac{k^2}{4}{y}_0(x){\mathcal{O}}^{-1}{y}_0(x)+\frac{k{y}_0(x{)}^2}{2}\right)\right],$$

$$F=F_{SHO}+T{\int }_0^{L_0}dx\left(-\frac{k^2}{4}{y}_0(x){\mathcal{O}}^{-1}{y}_0(x)+\frac{k{y}_0(x{)}^2}{2}\right).$$

$$\mathcal{P}\left[y_0(x)\right]\sim \exp \left[-\frac{1}{2\xi }{\int }_0^{L_0}dx{\dot{y}}_b^2\right]=e^{-{\int }_0^{L_0}{y}_0\mathcal{G}{y}_{0}dx},$$

$$\mathcal{G}=-\frac{1}{2\xi }{\partial }^2$$

$$\left[F\right]=F_{SHO}+\Delta F,$$

$$\Delta F=-\frac{T{k}^2}{8}\text{mbox}Tr({\mathcal{O}}^{-1}{G}^{-1})+\frac{Tk}{4}\text{mbox}Tr{G}^{-1}.$$

$$\left[F\right]=\frac{L_0\sqrt{V_0\tau }}{4}\left[1+\frac{2T}{\xi \tau }\right].$$

$$L={\int }_0^{L_0}dx\sqrt{1+{\dot{y}}^2}\approx L_0+\frac{1}{2}{\int }_0^{L_0}dx{\dot{y}}^2.$$

$$\Delta L=\left[⟨L⟩\right]-{⟨L⟩|}_{V_0=0}$$

$$\frac{\Delta L}{L_0}=\frac{1}{8}\sqrt{\frac{V_0}{\tau }}\left[1-\frac{2T}{\tau \xi }\right].$$

$$\frac{E}{T}={\int }_0^{L_0}dx\left[\frac{{\ell }_{\mathrm{P}}{\ddot{y}}^2}{2}+\frac{m{\dot{y}}^2}{2}+\frac{k(y(x)-y_0(x){)}^2}{2}\right].$$

$$Z=\int \mathcal{D}y(x)e^{-E/T}.$$

$${\mathcal{O}}_{\kappa }=\frac{{\ell }_{\mathrm{P}}}{2}{\partial }^4-\frac{m}{2}{\partial }^2+\frac{k}{2}$$

$$Z=\frac{1}{\sqrt{\det {\mathcal{O}}_{\kappa }}}{e}^{-{\int }_0^{L_0}dx(-\frac{k^2}{4}{y}_0(x){\mathcal{O}}_{\kappa }^{-1}{y}_0(x)+\frac{k{y}_0(x{)}^2}{2})}.$$

$$F=F_{MHO}+T{\int }_0^{L_0}dx\left(-\frac{k^2}{4}{y}_0(x){\mathcal{O}}_{\kappa }^{-1}{y}_0(x)+\frac{k{y}_0(x{)}^2}{2}\right).$$

$$\left[F\right]=F_{MHO}-\frac{T{k}^2}{8}\text{mbox}Tr({\mathcal{O}}_{\kappa }^{-1}{G}^{-1})+\frac{Tk}{4}\text{mbox}Tr{G}^{-1},$$

$${\mathcal{O}}_{\kappa }=\frac{{\ell }_{\mathrm{P}}}{2}\left[{\partial }^2-{\omega }_1^2\right]\left[{\partial }^2-{\omega }_2^2\right],$$

$${\omega }_{1,2}^2=\frac{m}{2{\ell }_{\mathrm{P}}}\left[1\pm \sqrt{1-\frac{4{\ell }_{\mathrm{P}}k}{m^2}}\right].$$

$$z_n=n\pi /L_0,$$

$$\left[F\right]=F_{MHO}-\frac{\xi k^2}{2{\ell }_{\mathrm{P}}}\sum _{n=0}^{\infty }\frac{1}{z_n^2}\frac{1}{\left(z_n^2+{\omega }_1^2\right)}\frac{1}{\left(z_n^2+{\omega }_2^2\right)}+\frac{\xi k}{2}\sum _{n=0}^{\infty }\frac{1}{z_n^2}.$$

Full text

I Introduction

Semiflexible polymer networks are well known to trap prestress in their formation. Cross linking molecules typically lock-in curved and thus elastically stressed states of the filaments. As a consequence of this being an out of equilibrium process, the cross linkers may, however, trap more thermal energy – $k_{\rm B}T$ – per bending mode. The result is that cross-linking during network formation typically traps excess free energy, which then slowly bleeds out of the system. That relaxation process appears to lead to large, nonequilibrium stress fluctuations, and is associated with the glassy power-law rheology of the network at very low frequencies. Such low-frequency power-law rheology has been observed in both simulation Müller et al. (2014) and experiment Broedersz et al. (2010). Living cells are similarly observed to have soft, glassy, power-law rheology Sollich et al. (1997), albeit with a distinctly different power law exponent Lenormand et al. (2004); Bursac et al. (2005). The large nonequilibrium stress fluctuations have currently been observed solely in simulation and we suggest that they should make an intriguing target for future experiments.

The underlying dynamics of transiently cross linked semiflexible networks is likely to be fundamental to the mechanics of both the active biological and passive in vitro systems. Simulations suggest that both the stress fluctuations and this characteristic power-law rheology of transiently cross linked networks are associated with the reorganization of progressively larger sections of the network occurring on progressively longer time scales. Understanding these dynamics presents a theoretical challenge.

In an effort to better understand the excess free energy trapped in such networks, we consider the problem of a single semiflexible filament at temperature $T$ interacting with a quenched random potential. We hope that, in this single filament model, the quenched pinning potential mimics the effect of random cross linking sites that mechanically couple the filament in question to the surrounding network. In particular, we will examine the role of spatial correlations in the pinning potential, which, at least loosely speaking, introduces an effective mesh size of the surrounding network. We analyze the disordered-averaged conformational fluctuations of the filament and the elastic energy stored in the system as a function of persistence length of the filament and the correlation length of the potential. We compare these predictions to the results of large-scale finite element Brownian dynamics simulations of such a semiflexible filament in a random potential.

The study of the statistical mechanics of a single stiff filament in a random potential recalls a number of related systems in which a low-dimensional elastic object interacts with a quenched pinning potential. Examples include disorder-pinned domain walls between symmetry-equivalent ground states Kuntz and Sethna (2000), vortex lines in superconductors Marchetti and Nelson (1990); Hwa et al. (1993); Ertaş and Kardar (1994), and the three-phase contact line associated with the spreading of fluids as they wet a disordered substrate Narayan and Fisher (1993); Schäffer and Wong (2000). The key distinction between these systems and the one of current interest is the presence of a bending term in the elastic Hamiltonian of the filament. In cases where this term dominates the statistical weights of various filament conformations, i.e., when filament tension is sufficiently small or when examining bends on short enough lengths, we expect to obtain results distinct from those obtained for these previously studied systems.

A filament interacting with a quenched random potential can be characterized by three lengths: the filament length $L_{0}$, the thermal persistence length $\ell_{\rm P}$, and the correlation length $\xi$ of the potential. We will focus on the case of stiff filaments (although we introduce alternatives) in which the persistence length is typically longer than that of the filament itself: $L_{0}<\ell_{\rm P}$. In this case, there are still two distinct limits. One might imagine that either $\ell_{\rm P}\ll\xi$, in which case the filament should be flexible enough to follow the twists and turns of the local minima of the random potential. Alternatively one may consider the case where $\ell_{\rm P}\gg\xi$ and the filament is so stiff that the elastic energy cost for following the valleys of the potential becomes prohibitive.

Similarly, the problem is endowed with two energy scales: the thermal energy $T$ (we work in units where Boltzmann’s constant is set equal to unity) and the typical energy scale of the pinning potential $V_{0}$. The potential has dimensions of energy per length. The inverse length scale $\nu=V_{0}/T$ must control the states of the filament, allowing one to examine both “strongly pinned” on scales where $L_{0}\nu\gg 1$ versus “weakly pinned” $L_{0}\nu\ll 1$ cases. The former is more interesting for the system under consideration.

There are two quantities that will provide insight into the ensemble of filament configurations. These are the effective persistence length of a filament in the pinning potential. This length differs from the usual tangent vector correlation length of the filament due to its interaction with the pinning potential. It also provides a clear target for experimental studies of labeled filaments in networks. The effect of this extra filament bending imposed by the pinning potential (or the network in experiment) is that the filaments store excess elastic energy. We propose that this excess free energy is a prediction for the prestress in networks. The single filament model thus makes two predictions for nonequilibrium networks.

The remainder of this article is organized as follows. In section II we describe the filament Hamiltonian. We then provide in section III an analytical calculation of the averaged free energy of the filament in the strong-pinning limit of a random potential, described above. We start with the case of the completely flexible filament, i.e., the filament without bending energy, and then move forward to the more general case of the semiflexible filament. We then turn in section V to numerical simulations of the system using finite element Brownian dynamics of a geometrically exact Simo-Reissner beam model. We then summarize our results and their implications for pre-stress in networks in section VI, where we conclude with a proposal for new experiments. Finally, we provide an appendix addressing the applicability of the replica method to the problem explaining that it provides unphysical results for many of the measurable quantities.

II The model

We consider a filament in two dimensions that is anchored at both ends. We introduce a coordinate system in which these anchoring sites are at positions $(0,0)$ and $(L_{0},0)$ respectively. The anchoring sites are assumed to be able to generate arbitrary constraint forces necessary to hold the filament at this points, but to provide no constraint torques. The directed filament is further assumed to be free of overhangs, allowing its configuration to be described by a function $y(x)$. A representation for a numerical simulation (to be described in section V) of a few filaments (green lines) interacting with the potential (heat map) is shown in Fig. 1. We will allow the arclength of the filament between these two anchoring points to vary. In short, the anchoring points are reservoirs of extra length. A simple mechanical model of this situation can be thought of as follows. The filament fluctuates on a table whose height topography in a uniform gravitation field gives the pinning potential. The anchoring points may be thought of as holes in this table through which more filament may enter or exit the table’s surface. Weights may also be added below the table to enforce a fixed tension on the filament. We should notice that for the simulation we use a slightly different model, with one end fixed and the other free to move in one direction (see section V). We assume that for the small conformation the difference between these two models is negligible.

In the small bending limit ($\frac{dy}{dx}\ll 1$), which should be valid for filaments much shorter than their persistence length, we can write the filament’s energy functional as

[TABLE]

where $\kappa$ is the filament’s bending modulus defining a thermal persistence length $\ell_{\rm P}=\kappa/T$. $\tau$ is the tension imposed on the filament. $V(x,y)$ is the quenched random potential (with dimensions of energy per length), described in more detail below. We have introduced the notation $\dot{y}=\frac{dy}{dx}$. Then the classical partition function for such a filament at temperature $T=1/\beta$ is given by the path integral over all trajectories of the filament weighted by a Boltzmann factor obtained from Eq. 1:

[TABLE]

We will later consider averages of the free energy, obtained from Eq. 2 in the usual way: $F=-T\ln Z$, over an ensemble of random potentials.

Turning to the pinning potential, we can consider two rather simple forms for its probability distribution that allow for a finite spatial correlation length, but do not break rotational invariance. The first is inspired by the massive scalar field Lagrangian

[TABLE]

Here $\xi$ sets the correlation length and $V_{0}$ the energy scale of the potential. It is straightforward to see that this model generates an ensemble of random potentials in which the amplitude of each Fourier mode is selected as an independent Gaussian random variable from a distribution with zero mean and a width that depends on the magnitude of the wavenumber $k=|{\bf k}|$.

We may also consider a related problem in which the correlation length is assigned to the force rather than to the potential. Since the force is the (negative) derivative of the potential, one can obtain the necessary Gaussian probability distribution for force by introducing another derivative in Eq. 3. We obtain

[TABLE]

The first version of the potential in which the Fourier amplitudes of the scalar potential are Gaussian distributed generates quite large pinning forces at short correlation lengths since the slope of the potential is $V_{0}/\xi$. We speculate that the second version of the potential in which the spectrum of pinning forces on the filament is Gaussian distributed is a better approximation of the physical problem since one may then manipulate the correlation length of the potential (which is our approximation to the mesh size of a filament network) without changing the scale of forces to which the filament is subjected. We return to this point in our discussion. Hereafter we refer to the first type of random potential as the energy controlled distribution while the second will be called the force controlled distribution.

These distributions can also be expressed in terms of a probability distribution for the Fourier components of the potentials. Working in terms of those Fourier modes, we can also introduce another distribution with an exponential suppression of the higher Fourier modes:

[TABLE]

where the rectangular system has an area $A=L_{x}L_{y}$. This particular form of the random potential is not convenient for analytic calculations, but generates more numerically stable simulations. See Fig. 2 for examples of random potentials selected from these distributions. The pinning forces generated from these potentials are shown in Fig. 3. The energy controlled potential produces a complex force landscape with very short ranged correlations. We do not reproduce that vector field here.

III The Valley Approximation

There is a considerable simplification to be found if we assume that $L\nu\gg 1$. In this strong pinning limit, the ensemble of filament configurations is dominated by states where the filament is confined to the valleys of the pinning potential. We indeed observe this in simulations. There we also see instances in which filament crosses from one valley to another over a saddle point of the potential. We will address escapes from one valley to the next over saddle points later. Excluding such saddle points for now, we assume that the potential is roughly constant along the bottom of the valley and that the curvature of the potential in the direction orthogonal to the path along the valley floor is also constant. Thus, the local form of the pinning potential is given by

[TABLE]

We have introduced a curvature (spring constant) scale $\tilde{k}=\frac{V_{0}}{\xi^{2}}$. The path of the valley minimum $y_{0}(x)$ remains a random curve. To analyze the effect of the quenched distribution of such paths, we may either calculate physical quantities of interest for an arbitrary curve $y_{0}(x)$ and then average, or use the replica trick to handle the average over the potential simultaneously with the thermal averaging. The replica trick, however, provides unphysical results as we explain in appendix C.

III.1 Flexible filaments with the energy controlled distribution

We now explore the valley approximation first for a flexible polymer by setting $\kappa=0$. We require a finite tension so the polymer’s path can be considered to be nearly straight. In this limit the energy of the polymer in the random potential may be written as

[TABLE]

where we have scaled the parameters by temperature $T$: $m=\beta\tau$ and $k=\beta\tilde{k}$. We have also introduced the differential operator

[TABLE]

The partition function is given by the integral

[TABLE]

Calculating this partition sum is equivalent to performing the Euclidean path integral (with $x$ being the time-like coordinate) for a quantum particle with a mass $m$ in a harmonic potential with spring constant $k$. In the analogous quantum problem the effect of the pinning potential is to introduce a time-dependent force $-ky_{0}(x)$. Since the integral is Gaussian, we can obtain a closed form solution for the partition function for a particular path of the valley floor:

[TABLE]

If we choose $y_{0}(x)=0$ everywhere, Eq. 10 reduces to the partition function of an unforced oscillator. From this we obtain the prefactor, leading us to write

[TABLE]

where $Z_{\rm SHO}$ is the well-known result for the simple harmonic oscillator Kleinert (1995). It is then straightforward to compute the free energy $F$ of the flexible chain in a particular realization of the random potential, the one whose valley follows the path $y_{0}(x)$. This free energy is given by

[TABLE]

We now average the free energy in Eq. 12 over a distribution of the paths of the valley floor that is consistent with our previous random potential distribution given by Eq. 3. To obtain this we weight the paths of the valley floor $y_{0}(x)$ by

[TABLE]

where differential operator

[TABLE]

incorporates the correlation length of the pinning potential $\xi$. This weight is analogous to the Euclidean path integral of a free particle. The distribution of quenched potentials determined by Eqs. 6 and 13 is not identical to that given by Eq. 3. But the statistical weight associated with valley floors of these potentials has the same spatial correlations as those valley floors determined by the original potential distribution in Eq. 3.

The averaging over the distribution of these valley floors we obtain a correction to the simple harmonic oscillator (SHO) free energy

[TABLE]

where

[TABLE]

Here and throughout this article, we use the squared brackets $\left[\cdot\right]$ to indicate averages with respect to the quenched random potential. The angled brackets $\langle\cdot\rangle$ represent thermal averages. We compute these traces by diagonalizing the two relevant operators – see appendix A. We find that for long filaments (see appendix A) the disorder-averaged free energy is

[TABLE]

The correlation length of the quenched potential and the tension set a natural energy scale that controls the free energy correction. When that potential is sufficiently heterogeneous so that its valleys are quite tortuous on the scale of a Pincus blob Pincus (1976) $\xi<T/\tau$, the pinning potential has a significant effect on the free energy. The result above – Eq. 17 – cannot be extended to arbitrarily small tensions $\tau\rightarrow 0$ since our assumption that $x(y)$ is a well-defined function breaks down. In effect our use of the Monge gauge fails to adequately represent the polymer’s shape.

Finally, we compute the contribution to the length, as compared with the length of the filament in the absence of the pinning potential. We recall that the contour length of the filament

[TABLE]

where the Taylor expansion of the integrand is justified by the fact that the filament is nearly straight when under sufficient tension. The quantity of interest is thermal average $\langle L\rangle$. As the integrand is nonnegative, this average $\langle L\rangle$ is necessarily longer than the separation of the end points $L_{0}$. In fact, in the absence of a bending modulus (as assumed here), $\langle L\rangle$ is divergent; there is an infinite amount of contour length trapped in the high wavenumber modes of deformation. It is thus useful to define the change in the contour length of the polymer due to the quenched pinning potential. We introduce

[TABLE]

as the difference between the thermal average of the filament’s contour length when averaged again over ensemble of pinning potentials at fixed $V_{0}$ and $\xi$ and the same quantity with no pinning potential, i.e., the problem obtained by setting $V_{0}=0$. This difference remains finite. We find that

[TABLE]

Once again we see that the corrections due to the pinning potential enter through the dimensionless ratio of the Pincus blob length to the correlation length of the potential. The negative sign in Eq. 20 may appear to be counterintuitive. One might imagine that a shorter correlation length would, in fact, create more transverse undulations in the filament as it tries to follow the more sinuous potential minimum. A shorter correlation length $\xi$ would then be expected to increase $\Delta L$. It does not. The tensed flexible polymer has large undulations on length scales below that of the Pincus blob. The effect of a decreasing potential correlation length at fixed $V_{0}$ is to increase the curvature of the potential, making the harmonic constraint forces on the filament stronger. These larger forces straighten out the filament on scales below the Pincus blob length, decreasing $\Delta L$ by straightening it out on these small scales.

The first term in Eq. 20 increases the $\Delta L$ of the filament with increasing potential (or decreasing $\tau$). This reflects the expected effect of the potential. In a stronger potential the filament is forced to follow more precisely the tortuous valley of the potential minimum and thereby use more arclength. We return to this idea in the case of semiflexible filaments. From this analysis, we see that one may wish to consider an ensemble of random potentials for which the typical scale of the pinning forces remains fixed even as the correlation length is changed.

III.2 Semiflexible filaments: energy controlled distribution

We now include the bending modulus in the filament Hamiltonian. Using the same valley-based approximation for the random potential, the energy of the semiflexible filament takes the form

[TABLE]

We introduce $\ell_{\rm P}=\beta\kappa$, the normal thermal persistence length of the filament. The partition sum is then given by

[TABLE]

In Eq. 22 we restrict the paths by imposing boundary conditions such that the filament begins and ends at $x=0,L_{0}$ respectively. Moreover, it starts and ends at zero tangent angle with respect to the mean direction (along the $x$ axis): $\dot{y}(0)=\dot{y}(L_{0})=0$.

Following our previous procedure, we may formally integrate over all paths $y(x)$ by introducing the inverse of the differential operator

[TABLE]

and write the partition function as

[TABLE]

The prefactor in the above equation is the partition function of a modified harmonic oscillator (MHO), previously discussed in Ref. Kleinert (1995). Leaving the details of that aside for the moment, we write free energy of the system as

[TABLE]

We must now average this free energy over the ensemble of paths taken by the local potential minimum. We average over $y_{0}(x)$ in an ensemble where we weight each such path by Eq. 13. Doing so, we obtain

[TABLE]

where $G$ is defined by Eq. 14. Once again, we are required to compute the traces of the relevant operators, defined by Eqs. 14 and 23. The analogous exercise for the flexible filament was relegated to appendix A. We expand on that discussion for the case of semiflexible filaments here.

First we note that one can factor ${\mathcal{O}}_{\kappa}$ into two commuting operators

[TABLE]

where we have introduced the (potentially complex) frequencies:

[TABLE]

Expanding in a complete eigenbasis of the two operators those product make up ${\mathcal{O}}_{\kappa}$, we obtain a form of the disorder averaged free energy in terms of an infinite sum over the (discrete) eigenvalues of ${\mathcal{O}}_{\kappa}$ indexed by

[TABLE]

where $n=0,1,\dots$. The resulting free energy is

[TABLE]

The MHO free energy is given in the appendix B. If we take the limit of a vanishing bending modulus ($\ell_{\rm P}\rightarrow 0$), we find that $\omega_{1}\rightarrow\infty$ and $\omega_{2}\rightarrow\sqrt{k/m}$. This returns us to the previous calculated free energy (up to a constant) for the flexible filament in the disordered potential – see appendix A. As the bending modulus increases from zero, the two frequencies become complex when $\ell_{\rm P}\geq m^{2}/4k$ – see Eq. 28 . The free energy, however, remains real since $\omega_{1,2}^{2}$ are complex conjugates.

Finally, we note that the product of fractions in the above summation can be broken up into a set of three independently convergent sums. This allows one to write the disorder-averaged free energy in terms of a sum of cotangents of $\omega_{1,2}^{2}L_{0}$. We consider that solution in the limit of long and stiff filaments $\omega_{1,2}L_{0}\gg 1$, obtaining a simple algebraic expression

[TABLE]

Using Eq. 28, we reintroduce the original model parameters. That result is most succinctly expressed in terms of the dimensionless parameter

[TABLE]

In terms of $\phi$ we find the disorder-averaged free energy to be

[TABLE]

Physically, we see that this dimensionless quantity $\phi$ is the ratio of the Pincus blob size to the correlation length (as was observed in the flexible polymer case) multiplied by a correction factor that incorporates the strength of the random potential. Specifically, we see that this correction factor is given by the square-root of the ratio of two length scales in the problem: the persistence length $\ell_{P}$ of the filament and the arclength of the filament $\nu^{-1}=T/V_{0}$ required for the potential energy of pinning to equal thermal energy. We may interpret the effect of a finite persistence length as extending the size of the underlying Pincus blobs. This result extends the disorder-averaged free energy of a flexible polymer to the semiflexible regime, and may be compared to the previous result found in Eq. 17. The effect of the finite bending modulus enters in both the difference between $F_{\rm SHO}$ and $F_{\rm MHO}$ in the first term on the right hand side of Eq. 33 and in second term on the right hand side of the above equation, where $b$ enters solely through the dimensionless parameter $\phi$ defined in Eq. 32.

We can now calculate the effect of the pinning potential on the length of the chain between the two pinning points. Following the same approach and using the same definitions as used for the flexible filament, we now find

[TABLE]

The pinning potential always increases the arclength, but its contribution becomes smaller as the filament’s bending modulus is increased, i.e., for larger $\phi$. Using the result for the MHO and reintroducing the original model parameters, we write

[TABLE]

The first term in the above expression is proportional to the ratio of the Pincus blob size $\zeta=T/\tau$ to the bending length $\ell_{b}=\sqrt{\kappa/\tau}$, which sets the cross over length between a regime where tension dominates the statistical ensemble of filament configurations at longer lengths and the bending modulus at shorter lengths. The second term in the above expression is proportional to the ratio of the same Pincus blob size and the length scale set by the pinning potential $\nu^{-1}=T/V_{0}$. In the limit of very stiff filaments so that $\ell_{P}\nu\gg 1$ and $\phi\gg 1$ we see that the first term provides a negative change in length, as observed for the flexible polymer, but the second term increases the filament arclength with a contribution $\sim\sqrt{\zeta\nu}\phi^{-1/2}$. The interpretation is similar to that given regarding the flexible polymer. The first term produces filament straightening as discussed there. The second term provides for the increase of the arclength due to the tortuosity of of the potential minima.

III.3 Semiflexible filaments: force controlled distribution

We now consider the semiflexible filament to be interacting with the smoother, force controlled distribution. The calculation proceeds in a manner analogous to the previous two sections. The key difference is that the differential operator appearing in the statistical weight of valley paths ${\cal G}$ is replaced by a new one accounting for the fact that the typical curves along bottom of the valleys of the potential now have their own persistence length. Thus the differential operator ${\cal G}$ is replaced by

[TABLE]

Following the methods outlined above, we arrive at an expression for the free energy of the filament in the quenched potential, written in terms of a sum over the eigenvalues of the differential operators ${\cal O}_{\kappa}$ and ${\cal G}_{\rm F}$. We find that the disorder-averaged free energy is given by

[TABLE]

where the first term is the free energy of the modified harmonic oscillator, as discussed in appendix B. The second term is the correction due to complex geometry of the valleys of the potential. It is given by

[TABLE]

where $\Omega=1/\xi$. The frequencies $\omega_{1,2}$ are defined in Eq. 28 and $z_{n}$ is defined in Eq. 29. This result simplifies considerably when we examine the limit of very short or very long filaments. If the former case, where $\omega_{1,2}L_{0}\ll 1$, only zeroth term survives and we get

[TABLE]

In the latter case, where $\omega_{1,2}L_{0}\gg 1$, we may change the summation to an integration over $z=\frac{\pi n}{L_{0}}$, obtaining

[TABLE]

The summation for the general case can also be performed, producing a quite lengthy expression that do not reproduce here. Performing the integral in Eq. 40, we find

[TABLE]

We now wish to compute the average arclength of the filaments to observe the effect of the pinning potential upon their ensemble of conformations. This calculation involves taking the derivative of the free energy with respect to the tension, which is conjugate to the length. First, we write the free energy difference in terms of two auxiliary functions $f$ and $g$ whose argument is $r=m/\ell_{\rm P}$:

[TABLE]

where

[TABLE]

and $s=\sqrt{\frac{k}{\ell_{\rm P}}}$. Then taking the appropriate derivative, we compute the difference in excess arclength between the filament in the confining potential and the same filament without it

[TABLE]

The first term, which is the excess arclength of the semiflexible filament in a uniform harmonic potential, is calculated in appendix B.

We plot the $\Delta L$ versus applied tension in Fig. 4. There we see the decrease in excess arclength of the filament with increasing tension $\tau$. There are two regimes characterized by a different power laws $\Delta L\sim\tau^{a}+\mbox{const}$ (there is a finite constant in the low tension regime) in the low and high tension regimes, referred to as weak and strong in the figure caption. In the weak tension regime, the filament is better able to track the valley of the potential minimum. As the applied tension is increased, the ensemble of filament configurations becomes restricted to straighter ones that cannot follow these valleys with high fidelity.

In comparing the theoretical calculation to the numerics, we freely adjusted the curvature and correlation length of the valley to obtain the fit.

IV Prestress and excess free energy of the pinned filament

From the results for the excess free energy, we have a prediction for the increase of the energy density of a network due to the effect of pinning. This enhancement of the energy density should be interpreted as the observed pre-stress found in biopolymer networks. Pre-stress can be considered to be stored in at least two separate manners. First, there should be excess bending and stretching energy of the filaments as they are pinned by the network (represented by the quenched potential in our analysis) into configurations that store more than $T/2$ per bending mode – the amount of energy storage consistent with the equipartition theorem. Secondly, there may be excess energy stored in the strain associated with the cross links themselves. We cannot directly measure this quantity in our model. To obtain the excess energy stored in the filament due to the pinning potential, we simply compute the disorder average of the thermal expectation value of the squared amplitude of each Fourier mode of the filament.

IV.1 Strong tension

First we work under the assumption of strong tension in which the filament is nearly straight. We start by expanding the undulations of the filament

[TABLE]

In terms of these Fourier modes the energy in a particular realization of the force-controlled pinning ensemble is

[TABLE]

where $f_{n}$ are the Fourier modes of the pinning force so that

[TABLE]

In this Fourier expansion of the pinning force, we assume that the potential is defined in a box of size $[L_{x}\times L_{y}]$ with periodic boundary conditions and we set $L_{x}=L_{0}$. Now, we may use the equipartition theorem to demand that each Fourier mode stores $T/2$ energy. This leads to

[TABLE]

We now average the above result using the force-controlled distribution for the potential from Eq. 4 in Fourier representation

[TABLE]

where Fourier modes are defined in a standard way

[TABLE]

We now express $f_{n}$ from Eq. 48 in terms of the gradients of $V_{k_{x},k_{y}}$

[TABLE]

with $k_{x}=\frac{\pi n}{L_{0}}=z_{n}$ Using Eq. 52 and Eq. 50 we obtain

[TABLE]

Taking $L_{y}$ to be large we replace the summation with integration to obtain

[TABLE]

Since only the second term in Eq. 49 depends on the pinning potential the average over the quenched disorder yields an expression for the mean excess bending energy of the filament

[TABLE]

in which the second term contains all the information about the filament’s interaction with the pinning potential. That second term is

[TABLE]

The first term is merely the standard result from the equipartition theorem for a semiflexible filament MacKintosh et al. (1995). It follows similarly that contribution from the potential to the mean energies stored in filament tension and in the pinning potential are given by

[TABLE]

In Fig. 5 we plot the bending energy stored in the filament as a function of dimensionless wavenumber $q=\xi z_{n}$ for a variety of pinning potential strengths (at fixed external tension) (upper panel A) and a variety of tensions at a fixed value of the strength of the pinning potential (lower panel B). From dimensional analysis we note that there is a single scale that sets the strength of the pinning potential

[TABLE]

When this quantity is large $\tilde{E}\gg 1$ we expect the pinning potential to control the statistical ensemble of the filament configurations. Conversely, we expect high tension to straighten out the filament so that it cannot follow the local potential minima. Higher tension leads to both straighter typical filament configurations and configurations for which the effect of pinning becomes harder to distinguish against a background of thermal undulations. This transition between pinning dominated states of the filament and thermally dominated ones is wavenumber dependent. At sufficiently high wavenumbers $q=\xi z_{n}>q^{\star}$ the modes of the filament are generically freed from the pinning potential.

In Fig. 6 we observe this transition from strongly pinned modes, where the pinning contribution to the bending energy $E^{\rm pin}$ is greater than the thermal component $E^{\rm T}$ for $q<q^{\star}$ to effectively unpinned ones.

Since the work done by tension to extend the filament is $\tau\Delta L$, we can use Eqs. 56,57 to compute the excess length

[TABLE]

The sum is rapidly convergent so the first terms will dominate. In the limit of a short correlation length $\xi/L_{0}\ll 1$, $L_{0}\sqrt{\frac{\tau}{\kappa}}\gg 1$ we get

[TABLE]

The summation results in

[TABLE]

If instead we consider the case of the large bending we should take the limit $L_{0}\sqrt{\frac{\tau}{\kappa}}\ll 1$ and

[TABLE]

For the case of large correlation length, $\xi/L_{0}\gg 1$, large tension $L_{0}\sqrt{\frac{\tau}{\kappa}}\gg 1$ we get

[TABLE]

which results in

[TABLE]

where $\zeta(x)$ is Riemann $\zeta$-function. Since the correlation length is much larger than that the filament, we expect that the filament feels an essentially uniform force field, much like that of a hanging rope in a gravitational field. Indeed, the result for the excess length in this classical problem is $\Delta L=\frac{g^{2}L_{0}^{3}}{6\tau^{2}}$, which demonstrates the same power law dependence on tension $\tau$ and the separation of the end points $L_{0}$.

IV.2 Weak tension

For case of weak tension where the filament can better follow the potential minima, we are free to use the valley approximation. In this case the filament’s energy is given by Eq. 21. Doing the same calculation for the variance of each Fourier mode of the filament in a valley whose bottom curve is described by the Fourier modes of $y_{0}(x)$, $y_{0}^{n}$, we find

[TABLE]

Once again, the first term is independent of the disorder in the valley, but does depend on the curvature of the potential. This result corresponds to the case of a semiflexible filament in a straight parabolic potential Smith (2001). The second term corrects this result for the tortuosity of the valley. To compute this correction we note that

[TABLE]

From this result and from energy function – see Eq. 21, we immediately find that

[TABLE]

Again the first term represents the bending energy associated with the semiflexible filament in a straight parabolic potential. The potential decreases the bending energy, as is physically reasonable. The potential suppresses the normal thermal undulations of the filament. The second term, however, represents an increase in the mean bending energy associated with the curvature of the potential valleys. For a fixed local mean curvature of the pinning potential, $V_{0}/\xi^{2}=\mbox{const}$, the dependence of this tortuosity correction is nonmonotonic in wavenumber. We examine this in Fig. 7. There we see that the wavenumber dependent bending energy goes through a local maximum at the scale of the potential’s correlation length,i.e., where $z_{n}\xi=1$.

V Simulations

The established computational framework for the Brownian dynamics of semiflexible filaments proposed in Cyron and Wall (2012) and used e. g. in Müller et al. (2014) has been extended to account for forces resulting from the random potential field.

V.1 The finite element Brownian dynamics simulation framework

A single filament is modeled by nonlinear, geometrically exact, 3D Simo-Reissner beam theory and discretized in space using finite elements Jelenic and Crisfield (1999). In terms of the structural rigidity of the filament, we thus account for axial, torsional, bending, and shear deformation. To model the Brownian motion, we include viscous drag as well as thermal forces, each distributed along the entire filament length. More precisely, viscous forces and moments are computed assuming a quiescent background fluid and individual damping coefficients for translations parallel and perpendicular to the filament axis as well as rotation around the filament axis, respectively. Thermal forces are determined from the stochastic Wiener process in accordance with the fluctuation-dissipation theorem. Finally, an Implicit-Euler scheme is used to discretize in time and a Newton-Raphson algorithm solves the resulting nonlinear system of equations. Further details on this simulation framework including all formulae can be found in Cyron and Wall (2012).

V.2 Incorporation of the background potential field

As described in section I, the potential field $V$ acting on the filament has dimensions of energy per length. Its contribution to the virtual work required for the weak, variational formulation of the problem can be stated as

[TABLE]

where $\mathbf{r}\in\mathbb{R}^{3}$ is the centerline position and $s\in[0,L_{0}]$ denotes the arclength coordinate in the stress-free reference configuration of the filament. Subsequent discretization of the admissible centerline variations $\delta\mathbf{r}$ according to the finite element method yields the contributions to the discrete element force vector. We apply the trapezoidal rule on each finite element to numerically evaluate the integral along the filament. Regardless of the fact, that we only consider planar problems throughout this article, the entire simulation framework as well as Eq. 69 is capable of modeling arbitrary filament configurations in 3D. Note that the potential exerts forces on the filament, however, as it models a surrounding network, it is independent of the filament motion. This is commonly denoted as one-way coupling.

To mention the most important algorithmic details: In a pre-processing step, we use a random number generator, apply a discrete Fourier transformation, and finally a finite difference scheme to arrive at the force field $-\nabla V$ on a sufficiently fine grid in the entire simulation domain. In each iteration, we then interpolate these grid values to compute $-\nabla V$ at the current position of each node and evaluate Eq. 69 element-wise.

V.3 Simulation setup and results

The simulation setup consists of a single filament of length $L_{0}=20\mathrm{\mu m}$ and persistence length $L_{p}\approx 14\mathrm{\mu m}$. Its initial, stress-free reference configuration is chosen straight and parallel to the global $x$-axis, as shown in Fig. 8. By means of Dirichlet boundary conditions, the filament is constrained to the $xy$-plane and simply supported, i. e., free to rotate at both ends, however only free to move in $x$-direction at one endpoint. Its circular cross-section is specified by the area $A=1.9\times 10^{-5}\mathrm{\mu m}^{2}$, area moment of inertia $I=4.3{\times}10^{-12}\mathrm{\mu m}^{4}$ and polar moment of inertia $I_{p}=8.6{\times}10^{-12}\mathrm{\mu m}^{4}$. The material is defined by the Young’s modulus $E=1.3\times 10^{10}\mathrm{pN/\mu m^{2}}$ and the Poisson ratio $\nu=0.3$. Temperature is set to $T{\,=\,}293\mathrm{K}$ and the dynamic viscosity of the quiescent background fluid to $\eta{\,=\,}10^{-3}\,\mathrm{Pa\,s}$. The filament is discretized in space using $400$ linear beam finite elements and the time step size is chosen as $\Delta t=0.01\mathrm{s}$.

Two variants of the potential field have been considered in simulations. First, the potential with exponential suppression of high Fourier modes defined by Eq. 5 and second, the force controlled distribution defined by Eq. 4.

V.3.1 Results for the potential with exponential suppression of high wavenumber modes

For the results presented already along with the theoretical prediction in Fig. 4, we applied the potential from Eq. 5, using $3096$ Fourier modes, a correlation length of $\xi\approx 1.6\mathrm{\mu m}$ and $V_{0}=0.175\mathrm{pN}$. In addition, a single point force $\tau=10^{-2}\ldots 50\mathrm{pN}$ was applied in global $x$-direction to the (right) endpoint of the filament that is free to move in this direction. Each simulation was run for $5\times 10^{4}$ time steps. To speed up simulations, we made use of parallelization and simulated systems of five filaments at a vertical spacing much larger than $\xi$ and without any interactions between the filaments. Each data point in Fig. 4 results from the statistical ensemble of $70$ to $100$ such systems with five filaments each, depending on the deviation in results that was higher for the small tension values. Finally, the excess arclength $\Delta L$ is obtained from simulation data as the negative displacement of the (right) filament endpoint in $x$-direction.

V.3.2 Results for the force controlled distribution

While the simulations using the potential with the exponential suppression of the high wavenumber modes are more robust, the theory assumes a force controlled distribution. Therefore, unlike for the previous variant of the potential distribution, there are no adjustable parameters necessary to directly compare these numerical results to the theory. We consider two parameter sets: one representing the case of large tension and another one representing the case of large bending and zero tension.

The first parameter set is given as $V_{0}=1/256^{2}$ pN, $\tau=0.006$ pN, $\xi=1\mu$m, $L_{0}=5\mu$m. The simulation gives us $\langle\Delta L\rangle=2.0\pm 0.6\times 10^{-5}\mu$m, while Eq. 61 predicts $\Delta L\approx 1.3\times 10^{-5}\mu$m. The second parameter set is $V_{0}=1/256^{2}$ pN, $\tau=0$, $\xi=1\mu$m, $L_{0}=5\mu$m, $\kappa\approx 0.0125\,\mbox{pN}\,\mu\mbox{m}^{2}$. Here, the result from our simulations is $\langle\Delta L\rangle=2.2\pm 0.4\times 10^{-5}\mu$m, while Eq. 63 predicts $\Delta L\approx 1.2\times 10^{-5}\mu$m. We performed sixteen simulation runs for each parameter set.

VI Discussion

We have examined the statistical mechanics of a single semiflexible filament in a quenched pinning potential as a model for studying how the network environment changes the typical stored elastic energy of filaments and leads to prestress. Based on these calculations we propose that there are two experimental quantities for which we may make predictions even with our single filament model. The first is that we expect the pinning environment of the network to impose a different (and nonequilibrium) statistical weight to filament configurations. One way to parameterize this difference between the ensemble of filament configurations in a network and of a filament in isolation is that the effective persistence length of the network filament will no longer be $\ell_{\rm P}=\kappa/T$.

Using our results for the disordered-averaged Fourier modes of the filament’s undulations, we may directly compute the tangent tangent correlations. We find that

[TABLE]

where the nonthermal persistence length $\tilde{\ell}_{\rm P}$ is given by

[TABLE]

where $\omega_{1,2}$ are the eigenvalues introduced in Eq. 28. We examine the dependence of the effective, disorder-influenced persistence as a function of tension and potential correlation length in Fig. 9.

For the case of sufficiently strong tension or weak pinning potential $\tau\gg V_{0}$, we find that the tangent correlations do not decay exponentially with separation along the filament. Instead they decay as a polylogarythmic function of the separation: $\mbox{Li}_{2}(e^{i\pi|x_{1}-x_{2}|/L})$. In that case, no persistence length can be defined. See Fig. 10 for the behavior of the correlation function

Turning to the case of prestress, we offer a prediction for the mean excess free energy density of the network due to trapped elastic deformations of its constituent filaments. To arrive at this prediction within our one filament model, we assume the excess free energy density may be computed by summing the excess bending and tension energy $\Delta E^{\rm bend}+\Delta E^{\rm ten}$ of a filament due its interaction with the pinning potential and then dividing that quantity by the area occupied by that filament. In our two dimensional calculation this is simply $L_{0}\times\xi$. To make a prediction for a three dimensional network we assume that there are two independent polarization states of the filament’s undulations (which is reasonable for small bending angles) resulting in a prediction

[TABLE]

We expect that this quantity should set the scale for the anomalous nonequilibrium stress fluctuations observed in transiently cross linked networks of semiflexible filaments.

There are limits to the single filament description of the collective phenomenon of network structure and excess free energy. Following the example of mean field models of magnets, one might imagine finding a self consistent description of network structure in which the statistics of the pinning potential are determined by the calculated properties of a filament in that potential in a sort of Weiss molecular field description. We leave such self-consistent calculations to future work.

We also note that our description of a single filament is inadequate for studying the crumpling of filaments in a even stronger pinning environments. States of larger deformation, as might be expected in flexible polymers without sufficient tension, cannot be described by our framework; their description requires more complex and inherently nonlinear elasticity. Fortunately, nature provides numerous examples of semiflexible protein filaments for which our analysis should be sufficient. Transiently cross linked networks of such stiff filaments are an arena for the study of the role of quenched disorder on their ensemble of shapes and elastic energy storage.

Acknowledgements.

AJL and VMS acknowledge partial support from NSF-DMR-1709785. AJL and VMS are thankful to Robijn Bruinsma and Jonathan Kernes for fruitful discussions. WAW and MJG acknowledge partial support from BaCaTeC.

Appendix A Calculating the traces

In this appendix we discuss the calculation of the trace of the operator $\mathcal{O}^{-1}\mathcal{G}^{-1}$ obtained in Eq. 16. The operators $\mathcal{O}$ and $\mathcal{G}$ were introduced in Eqs. 8 and  14 respectively. To calculate the traces, we need the complete spectrum of the two differential operators in question. In other words, we must solve

[TABLE]

subject to the fixed-end boundary conditions $\Psi(0)=\Psi(L_{0})$. The solution is a constant and the set of standing waves satisfying the boundary conditions: $\Psi(y)=A\sin(z_{n}y)$, where, as discussed in the main text, $z_{n}=n\pi/L_{0}$ for nonnegative integer $n$. Defining $\omega=\sqrt{\frac{k}{m}}$ we may write the trace as

[TABLE]

Since both operators can be simultaneously diagonalized, it is also possible to write

[TABLE]

Using the above results and Eq. 16, one finds directly that

[TABLE]

The second term in the product occurring in the summand removes the singularity associated with the $n=0$ term.

[TABLE]

The summation can be done in closed form leading to the disorder-averaged free energy:

[TABLE]

Appendix B The Modified Harmonic Oscillator

In this appendix we review the partition function of the modified harmonic oscillator. The modified harmonic oscillator Hamiltonian has an addition $p^{4}$ where $p\rightarrow-i\partial$ is the momentum operator. This leads to a time-independent Schrödinger operator of the form $\mathcal{O}_{\kappa}$ introduced to discuss the semiflexible filament in the text. A more complete exposition of this problem can be found in Refs. Kleinert (1995) and Feynman and Hibbs (1965). We begin by factoring the differential operator $\mathcal{O}_{\kappa}$, defined in Eq. 23 into a product of two commuting second order differential operators as shown in Eq. 27:

[TABLE]

where

[TABLE]

It immediately follows that the partition sum for the modified harmonic oscillator is given by a product of the partition functions of two separate harmonic oscillators $Z_{1}$ and $Z_{2}$,

[TABLE]

with frequencies $\omega_{1}$ and $\omega_{2}$ respectively. The free energy immediately follows. From Eq. 81, it is the sum of the free energies of the two harmonic oscillators introduced above in Eq. 79,

[TABLE]

The remainder of this appendix uses the above result.

We first compute the free energy of the filament in a confining potential (finite $k$). Since the mass term in the equation of motion represents tension $m=\tau/T$ – see Eq. 8– and since tension is conjugate to arclength in the Hamiltonian, the derivative of the free energy with respect to mass gives us mean arclength that we seek.

[TABLE]

There is one complication. The free energy is divergent in the limit that $k\rightarrow 0$. We require the arclength of the unconfined filament in order to compute $\Delta L$, as defined in Eq. 19. To address this problem, we compute the mean arclength at finite $k$, take the necessary derivative, and then take the limit $k\rightarrow 0$, which then provides a finite result.

We first recall the (possibly complex) frequencies $\omega_{1,2}$ and use them to compute derivatives of the free energy in Eq. 82.

[TABLE]

We know that $\omega_{1}^{2}+\omega_{2}^{2}=\frac{m}{\ell_{\rm P}}$ and that $\omega_{1}\omega_{2}=\sqrt{\frac{k}{\ell_{\rm P}}}$. When there is no confining potential $k=0$, one frequency vanishes. This implies that the spectrum of the corresponding operator includes a zero eigenvalue, which, as discussed above, will be problematic for the analysis. In preparation for taking the $k\rightarrow 0$ limit we keep only the lowest terms in $k$ here. Expanding to lowest order in $k$ we find that

[TABLE]

which means

[TABLE]

Solving for $\omega_{1,2}$ we find that

[TABLE]

These results allow us to take the appropriate $m$ derivatives:

[TABLE]

By assuming that the filament is long we may take $\omega_{1}L_{0}\gg 1$. We also assume that $k$ may be sufficiently small so that $\omega_{2}L_{0}\ll 1$ (we will take it to zero shortly). The free energy of the filament is then

[TABLE]

Taking the derivative as shown in Eq. 83, we obtain the mean length

[TABLE]

For large $L_{0}$ the last term can be ignored, leaving us with the $k=0$ result.

We now turn to the case where $k$ remains finite. Then for long filaments we have $\omega_{1}L\gg 1$ and $\omega_{2}L\gg 1$, allowing us to write $F_{MHO}$, defined in Eq. 82 as

[TABLE]

Then, taking the $m$ derivative as above, we obtain the mean arclength of the filament with a confining potential $k\neq 0$:

[TABLE]

Subtracting the equivalent quantity for the $k=0$ case, shown in Eq. 92, we obtain an expression for the change in excess mean arclength due to the presence of the confining potential

[TABLE]

Returning the original, physical parameters of the semiflexible filament model, this expression becomes

[TABLE]

Appendix C Replica Trick

The core idea of the replica trick is the following. We want to find $\left[\ln Z\right]$ where the average is over the ensemble of quenched potentials. We may rewrite this in the following form

[TABLE]

In many cases the straightforward calculation of the left hand side is intractable, while it is possible to compute on the right hand side for positive integer $R$ and then take the limit. The reader is referred to Ref. Mezard et al. (1987) for further details of this approach. We first demonstrate the utility of the replica method in a simplified model. We assume a Gaussian Hamiltonian so that the partition function has the form

[TABLE]

where $A$,$B$ are operators acting on the variable $y$. The variable $y$ is affected by a quenched potential, $y_{0}$ whose statistics are fixed by the Gaussian distribution

[TABLE]

controlled by another operator $C$. This problem is designed to be sufficiently simple that it can be solved without the replica trick. It is also directly related to valley approximation, as the reader may confirm. A straightforward calculation of the partition function, followed by taking the logarithm and then averaging the resulting free energy over the quenched potential distribution leads to

[TABLE]

If we now repeat the calculation using the replica trick, we first replicate the Hamiltonian to form $Z^{R}$:

[TABLE]

where the replicated variables are indexed by $r$: $y_{r},\,r=1,\ldots,R$. We now average the replicated partition sum over the quenched potential (notice we are swapping the order of the two averages) to obtain

[TABLE]

Introducing $\vec{y}=\left(y_{0},y_{1}...y_{R}\right)$ we can write Eq. 102 in the form

[TABLE]

where the matrix $\mathcal{A}$ can be written in the block form

[TABLE]

The $\tilde{A}$ block is the $R\times R$ matrix coupling the replicated variables. The $\tilde{B}$ block is an $R\times 1$ matrix coupling the $R$ replicated variables to the quenched disorder $y_{0}$. The remaining calculation involves performing the Gaussian integrals over $\vec{y}$. We find

[TABLE]

The determinant of the block matrix can be written as

[TABLE]

and we observe that, in our model, the variables $y_{r}$, $r=1,\ldots,R$ are noninteracting so that the $\det{\tilde{A}}=\det A^{R}$. As a consequence,

[TABLE]

By taking the $R\rightarrow 0$, we reproduce the exact result, Eq. 100.

Now, we consider using the replica trick for the problem at hand. Specifically, we examine a flexible polymer at temperature $T$ interacting with a quenched delta-correlated random potential. The polymer has tension $\tau$ and length $L$. We will examine the average potential energy of the filament in the pinning potential. For later computational convenience, we introduce a coupling constant $\alpha$ controlling the interaction of the polymer with pinning potential.

The partition function is

[TABLE]

We can calculate the average potential energy as

[TABLE]

After this thermal average, we take the average of the potential over the ensemble of pinning potentials

[TABLE]

We assume that averaging and differentiation commute. We expect that the mean potential energy should be negative. Indeed, the random potential will pick negative and positive values with the same probability. The polymer, however, will prefer the negative regions to the positive.

We now perform the necessary averages using the replica trick with, as before, the replica index being $r=1..R$. The replicated partition sum is now

[TABLE]

The average over the pinning potential, once again, couples the replicas. We first write that average over the distribution of $V(x,y)$,

[TABLE]

Since ${\cal P}(V)$ is a Gaussian with fixed variance $\sigma$ we can directly perform the Gaussian integral over $V$. We also notice that the delta functions $\delta(y-y_{r}(x))$ allow us to do the integral over $y$ in the exponent. Doing both steps we obtain

[TABLE]

Eq. 113 now appears to be the partition function of $R$ particles interacting through a delta-potential. It can be written in the form

[TABLE]

where the multi-particle Hamiltonian is

[TABLE]

Notice that in Eq. 114 that the length of the polymer $L_{0}$ plays the role of inverse temperature $1/T$ for the fictitious particles. Of course, $L_{0}\neq\beta$, which is the inverse temperature of the polymer. We also note that $x$ plays the role of time in this dynamical system of interacting particles.

Since the sum contains the term $i=j$, there is a positive infinite constant $C$ in the Hamiltonian. Separating these terms explicitly we rewrite the Hamiltonian as

[TABLE]

Reflecting on the fact that the path integral in Eq. 114 is analogous to the quantum transition amplitude for the system of fictitious particles  Polyakov (1987) and that the long polymer limit $L_{0}\rightarrow\infty$ corresponds to the zero temperature limit for that system, we may focus on the ground state wavefunction Kardar (1987, 2007). To find the ground state wavefunction, we use the Bethe ansatz:

[TABLE]

The constant $C_{0}$ is chosen to satisfy the normalization condition. The constant $\kappa$ is set so that the discontinuity of the first derivative cancels the delta-function, i.e., $\frac{2}{\beta\tau}\kappa=2\sigma\beta^{2}\alpha^{2}$ .

We transform the sum by choosing a particular ordering of the $R$ particles on the $y$ line: $y_{1}<y_{2}<\ldots<y_{R}$. In this ordering the sum in Eq. 117 is particularly simple:

[TABLE]

From this we obtain the ground state energy

[TABLE]

Then the partition function is

[TABLE]

and the free energy is

[TABLE]

Taking the limit we get

[TABLE]

Using Eq. 110, we obtain the potential energy of the thermalized polymer in the random pinning potential, averaged over an ensemble of such potentials:

[TABLE]

This potential energy is positive definite, which is clearly unphysical. The result shows that the replica trick fails to find the averaged potential energy of the polymer in the pinning potential. The problem may lie in either the introduction of the infinite constant $C$, or in the replica trick itself.

We conjecture that the failure of the replica trick in this case is related to the failure of the analytical continuation inherent in the method. To expand on this idea we consider an example related to directed polymers (without a bending energy) in a complex random potential with zero correlation length ($\xi\rightarrow 0$), as explored by Zhang Zhang (1989). In that work [Eq. (11) of that article] the author asserts that

[TABLE]

where $P(F)$ is the probability distribution of free energies of the directed polymer due to the statistical ensemble of quenched potentials in its environment and $n$ is the replica index. The left hand side of this expression is $\left[Z^{n}\right]$. Now $P(F)$ should decay exponentially or faster for negative $F$, otherwise the integral will be divergent. The result on the right hand side, arrived at through a saddle point evaluation of the integral, is exact for all positive integer $n$ in the limit of large $t$. The goal (as in our case too) is to examine the $n\rightarrow 0$ limit.

To examine this point let us perform the integration in Eq. 124 by first Taylor expanding the exponential and then integrating term by term.

[TABLE]

From Ref. Kolokolov and Korshunov (2009) we know that the large positive $F$ behavior of $P(F)$ is $\sim\exp(-F^{\frac{5}{2}})$, ensuring the convergence of these integral and justifying the swapping of the order of summation and integration. The sum on the right hand side can be reorganized in terms of powers of the replica index $n$ so that each term in the summand takes the form $d_{k}n^{k}$. In this form that sum can be combined with integral on the left to write an infinite sum that is now equal to zero.

[TABLE]

where

[TABLE]

If we were allowed to assert that each term in the sum Eq. 126 separately vanished, then we would be able to calculate all moments of the free energy directly from the known coefficients $d_{k}$. This process for $k=1$ is the essence of the replica trick as used to compute the mean free energy. Unfortunately, this conclusion is not necessarily valid. A simple counter example can be obtain from the Taylor expansion of the function $F(z)=f(z)\sin(\pi z)$, where $f(z)$ is an analytic function. Let $a_{k}$ be Taylor expansion coefficients for $F(z)$. Then it is easy to see that if we choose $c_{k}=\frac{a_{k}}{k!}$ we will get Eq. 126, but not all $c_{k}$ are zero. The same issue appears in Ref. Kardar (1987).