Exact solution for the force-extension relation of a semiflexible polymer under compression

TL;DR

This paper provides an exact analytical solution for the force-extension behavior of a semiflexible polymer under compression, revealing how thermal fluctuations influence buckling and related thermodynamic properties.

Contribution

It introduces an exact solution using Mathieu functions for the wormlike chain model under compression, detailing the impact of thermal fluctuations on buckling and thermodynamics.

Findings

Thermal fluctuations cause larger elongations than elastic rods at strong compression.

Susceptibility peaks near the Euler buckling force, increasing with persistence length.

Entropy and heat capacity vary non-monotonically with applied load.

Abstract

Exact solutions for the elastic and thermodynamic properties for the wormlike chain model are elaborated in terms of Mathieu functions. The smearing of the classical Euler buckling instability for clamped polymers is analyzed for the force-extension relation. Interestingly, at strong compression forces the thermal fluctuations lead to larger elongations than for the elastic rod. The susceptibility defined as the derivative of the force-extension relation displays a prominent maximum at a force that approaches the critical Euler buckling force as the persistence length is increased. We also evaluate the excess entropy and heat capacity induced by the compresssion and find that they vary non-monotonically with the load. These findings are corroborated by pseudo-Brownian simulations.

Figures

Equations

$${\mathcal{H}}_0$$

$$Z_0({\mathbf{u}}_L,L|{\mathbf{u}}_0,0)$$

$${\mathcal{H}}_{\text{force}}$$

$$\frac{\mathcal{H}}{k_{\text{B}}T}$$

$$Z({\mathbf{u}}_L,L|{\mathbf{u}}_0,0)$$

$$=Z_0⟨\exp \left({\int }_0^L\mathrm{d}sf\mathbf{e}\cdot \mathbf{u}\right)⟩.$$

$$G(T,F)$$

$$\mathrm{d}G$$

$$⟨X⟩$$

$$\chi$$

$$\chi (T,F)$$

$${\left(\frac{\partial S}{\partial F}\right)}_T$$

$$\Delta S(T,F)$$

$$\Delta C_F(T,F)$$

$${\partial }_{s}Z(\mathbf{u},s|{\mathbf{u}}_0,0)$$

$$Z(\mathbf{u},s=0|{\mathbf{u}}_0,0)$$

$${\partial }_{s}Z(\phi ,s|{\phi }_0,0)$$

$$Z(\phi ,s=0|{\phi }_0,0)$$

$$\left[\lambda +f\cos \phi +\frac{1}{{\ell }_{\text{p}}}\frac{{\mathrm{d}}^2}{\mathrm{d}{\phi }^2}\right]z(\phi )$$

$$\left[\frac{{\mathrm{d}}^2}{\mathrm{d}{x}^2}+\left(a-2q\cos (2x)\right)\right]z(x)=0,$$

$${\text{ce}}_{2n}(q,x)$$

$${\text{se}}_{2n+2}(q,x)$$

$$a_{2n}{A}_0^{2n}-q{A}_2^{2n}$$

$$(a_{2n}-4)A_2^{2n}-q(2A_0^{2n}+A_4^{2n})$$

$$(a_{2n}-4m^2)A_{2m}^{2n}-q(A_{2m-2}^{2n}+A_{2m+2}^{2n})$$

$$Z({\phi }_L,L|{\phi }_0,0)=$$

$$\displaystyle\frac{1}{2\pi}\sum_{n=0}^{\infty}\bigl{[}\text{ce}_{2n}(q,\varphi_{0}/2)\text{ce}_{2n}(q,\varphi_{L}/2)e^{-a_{2n}(q)L/4\ell_{\text{p}}}$$

$$\displaystyle+\text{se}_{2n+2}(q,\varphi_{0}/2)\text{se}_{2n+2}(q,\varphi_{L}/2)e^{-b_{2n+2}(q)L/4\ell_{\text{p}}}\bigr{]}.$$

$$Z(L|{\phi }_0,0)={\int }_0^{2\pi }\mathrm{d}{\phi }_{L}Z({\phi }_L,L|{\phi }_0,0),$$

$$=\sum _{n=0}^{\infty }{A}_0^{2n}(q){\text{ce}}_{2n}(q,{\phi }_0/2)\exp [-a_{2n}(q)L/4{\ell }_{\text{p}}],$$

$$Z(L)$$

$$a_0(q)$$

$$a_{2n}(q)$$

$${\text{ce}}_0(q,0)$$

$${\text{ce}}_2(q,0)$$

Full text

I Introduction

The elasticity of biopolymers such as actin, microtubuli, or intermediate filaments is responsible for the mechanical and structural stability of cells, their motility and intracellular transport processes Sackmann (1994); Brangwynne et al. (2006); Bausch and Kroy (2006); Fletcher and Mullins (2010); Lieleg et al. (2010); Nolting et al. (2014). Networks of these cytoskeletal polymers connected by regulatory proteins are exposed to mechanical stresses, thereby exhibiting striking nonlinear elastic behavior MacKintosh et al. (1995); Kroy and Frey (1996); Storm et al. (2005); Claessens et al. (2006); Chaudhuri et al. (2007); Huisman et al. (2008); Carrillo et al. (2013); Razbin et al. (2015); Amuasi et al. (2015); Plagge et al. (2016). An adequate characterization of the mechanical properties of these biological assemblies constitutes a fundamental step toward the future design and synthesis of artificial biomimetic materials Ratner et al. (2004); Ratner and Bryant (2004). Since the macroscopic behavior of these filamentous materials strongly depends on the properties of their constituents, a profound knowledge on the nature of single filaments is necessary to fully understand the elasticity of these networks Hugel and Seitz (2001).

Due to their semiflexibility, biopolymers show a peculiar response to external forces Ott et al. (1993); Broedersz and MacKintosh (2014); Razbin et al. (2016), where the static and dynamic properties are dominated by their enthalpic elasticity, but conformational entropy still plays a significant role. Experimental techniques such as optical Ashkin (1997); Mehta et al. (1999) and magnetic tweezers Gosse and Croquette (2002), transmission electron microscopy Kuzumaki and Mitsuda (2006), and (atomic) force spectroscopy Janshoff et al. (2000); Hugel and Seitz (2001) have been used to quantify the mechanical properties in terms of force-extension relations Saleh (2015). These capture the nonlinear effects emerging in the stretching and buckling of semiflexible polymers such as DNA Marko and Siggia (1995); Bouchiat et al. (1999); Bustamante et al. (2000), actin Liu and Pollack (2002), and synthetic carbon nanotubes Kuzumaki and Mitsuda (2006), but also single molecules such as titin Kellermayer et al. (1997), which is an important component in striated muscle tissues, and collagen Sun et al. (2002), present in, e.g., skin and bones.

To characterize theoretically the physics of the nonlinear elastic behavior of semiflexible polymers, the wormlike chain (WLC), also referred to as Kratky-Porod model Kratky and Porod (1949) has been analyzed in terms of the end-to-end distribution function in the weakly-bending approximation  Wilhelm and Frey (1996) as well as in terms of the exact low-order moments Hamprecht et al. (2004); Spakowitz et al. (2004); Spakowitz and Wang (2005); Mehraeen et al. (2008). The WLC model has been shown to be a reliable model for semiflexible polymers, as it reproduces, e.g., the mechanical behavior of DNA Marko and Siggia (1995) and actin Liu and Pollack (2002), while for microtubuli internal shearing of adjacent protofilaments leads to an apparent length dependence of the elastic moduli Taute et al. (2008); Pampaloni et al. (2006). In these studies, approximate solutions of the force-extension relation of stretched polymers have been compared to experimental data in order to determine mechanical properties such as the persistence length.

The approximations have been complemented by exact expressions for stretched polymers in the plane Prasad et al. (2005), while force-extension relations of polymers under a compressive load have so far only been analyzed within the regime of stiff polymers Baczynski et al. (2007); Emanuel et al. (2007); Lee et al. (2007); Bedi and Mao (2015). Therefore, a characterization of the buckling behavior for a broad range of semiflexible polymers is needed to predict elastic properties within experimental observations.

Here, we provide an exact analytical solution for the force-extension relation of a compressed semiflexible polymer within the framework of the WLC model and compare the results for polymers with different rigidity to computer simulations. Boundary conditions imposed within the experimental setup play a crucial role for the response of polymers to external forces, as has been discussed earlier for stretched polymers Ghosh et al. (2007). To quantify this effect, we also compare the force-extension relations with two clamped ends, one free and one clamped end, and two free ends. Thermodynamic properties, such as the in principle measurable excess heat capacity of the semiflexible polymers induced by the compression, are also discussed.

II The wormlike chain model

For the elastic properties of a semiflexible polymer, we rely on the WLC model Kratky and Porod (1949), where the bending energy of an inextensible semiflexible polymer is expressed by its squared curvature,

[TABLE]

Here, $s$ is the arc length, $\mathbf{u}(s)=\mathrm{d}\mathbf{r}(s)/\mathrm{d}s$ denotes the tangent vector of the polymer along its contour $\mathbf{r}(s)$, $\kappa$ the bending stiffness of the polymer, and $L$ its contour length. The corresponding partition sum $Z_{0}(\mathbf{u}_{L},L|\mathbf{u}_{0},0)$ of the polymer with initial orientation $\mathbf{u}_{0}$ and final orientation $\mathbf{u}_{L}$ is obtained as a sum of Boltzmann weights over all possible configurations,

[TABLE]

where the inextensibility constraint $|\mathbf{u}(s)|=1$ has to be fulfilled. In addition, if the polymer is subject to a constant external force $\mathbf{F}$, the stretching energy,

[TABLE]

has to be added to the bending energy. Here, the force acts along a fixed direction $\mathbf{e}$ with strength $F$ such that $\mathbf{F}=F\mathbf{e}$. Collecting terms, the full Hamiltonian of the system reads

[TABLE]

where $f=F/k_{\text{B}}T$ is the reduced force with units of an inverse length. The persistence length $\ell_{\text{p}}=\kappa/k_{\text{B}}T$ for 3D, respectively $\ell_{\text{p}}=2\kappa/k_{\text{B}}T$ for 2D, corresponds to the decay length of the tangent-tangent correlations of the polymer and permits a classification of polymers into three categories: polymers with $\ell_{\text{p}}/L\ll 1$ are referred to as flexible and therefore more coil-like structured, $\ell_{\text{p}}/L\sim 1$ are semiflexible, and $\ell_{\text{p}}/L\gg 1$ are stiff, hence, rodlike Doi and Edwards (1986); Broedersz and MacKintosh (2014). The partition sum $Z(\mathbf{u}_{L},L|\mathbf{u}_{0},0)$ of such a polymer subject to an external force is given by the path integral over all weighted, accessible chain configurations. In particular, except for the prefactor $Z_{0}\equiv Z_{0}(\mathbf{u}_{L},L|\mathbf{u}_{0},0)$, the partition sum can be viewed as the generating function of a semiflexible polymer:

[TABLE]

The corresponding thermodynamic potential is the Gibbs free energy defined by

[TABLE]

with fundamental relation

[TABLE]

where $S=S(T,F)$ denotes the entropy of the polymer and $\langle X\rangle:=\int_{0}^{L}\mathrm{d}s\ \langle\mathbf{e}\cdot\mathbf{u}(s)\rangle$ is the mean end-to-end distance projected onto the direction of the force.

Hence, the mean projected end-to-end distance is obtained by

[TABLE]

In particular, the response of the polymer to the applied force is encoded in the isothermal susceptibility,

[TABLE]

which characterizes the strength of the response with respect to the applied force. Furthermore it serves as an indicator for the buckling force of semiflexible polymers in analogy to the critical Euler buckling force for stiff rods Landau and Lifshitz (1986).

We can also obtain the linear susceptibility $\chi$ in terms of the fluctuation-response theorem Doi and Edwards (1986) for the semiflexible polymer,

[TABLE]

In addition to these elastic properties, we obtain the change in the entropy with respect to the force at constant temperature by using the Maxwell relation:

[TABLE]

Thus, the excess entropy $\Delta S(T,F)=S(T,F)-S(T,F=0)$ can then be evaluated by

[TABLE]

Similarly, the experimentally accessible excess heat capacity $\Delta C_{F}(T,F)=C_{F}(T,F)-C_{F}(T,F=0)$ is obtained by

[TABLE]

II.1 Analytic solution

To determine the elastic and thermodynamic properties, the partition sum has to be computed by solving for the path integral [Eq. (5)]. Therefore, we discretize the path integral in terms of the arc-length Kleinert (2009) and find that the partition sum obeys an equation of the Schrödinger type, describing the evolution of the partition sum along the contour of the polymer,

[TABLE]

where $\Delta_{\mathbf{u}}$ denotes the angular part of the Laplacian, and with initial condition

[TABLE]

such that the delta function $\delta(\cdot,\cdot)$ enforces both orientations to coincide.

Here we restrict the discussion to a polymer confined to a plane. Consequently, the inextensibility constraint of the orientation can be parametrized in terms of a polar angle, $\mathbf{u}=(\cos\varphi,\sin\varphi)^{T}$, were the angle $\varphi$ is measured with respect to the direction of the applied force. Thus, the evolution of the partition sum along the arc length reads

[TABLE]

subject to the initial condition

[TABLE]

The Fokker-Planck equation [Eq. (17)] is reminiscent of the Schrödinger equation of a quantum pendulum Aldrovandi and Ferreira (1980) and can be solved analytically by separation of variables in terms of appropriate angular eigenfunctions, $\exp(-\lambda s)z(\varphi)$. Inserting this ansatz into Eq. (17), we obtain

[TABLE]

Here we focus on a polymer exposed to a compressive load, in particular, we rewrite the reduced force as $f=-|f|$. A change of variables $x=\varphi/2$ and rearranging terms, leads to the equation

[TABLE]

which is known as the Mathieu equation DLMF ; Olver et al. (2010) with deformation parameter $q=2|f|\ell_{\text{p}}$ and eigenvalue $a=4\ell_{\text{p}}\lambda$. Thus, the general solution is expressed as a linear combination of $\pi$-periodic even and odd Mathieu functions $\text{ce}_{2n}(q,x)$ and $\text{se}_{2n+2}(q,x)$ with corresponding eigenvalues $a_{2n}\equiv a_{2n}(q)=4\ell_{\text{p}}\lambda_{n}$ and $b_{2n+2}\equiv b_{2n+2}(q)=4\ell_{\text{p}}\lambda_{n}$ DLMF ; Olver et al. (2010), respectively. Note that in the case of a polymer under tension, where the sign of the force is positive, $f=|f|$, a different change of variables $x=\pi/2-\varphi/2$ has to be used in order to arrive at the Mathieu equation Prasad et al. (2005).

The Mathieu functions are essentially deformed sines and cosines:

[TABLE]

The coefficients $A_{2m}^{2n}(q)$ and $B_{2m+2}^{2n+2}(q)$ are determined by recurrence relations, which result from inserting the Fourier series into Eq. (20):

[TABLE]

and similar relations hold for the coefficients $B_{2m+2}^{2n+2}(q)$ DLMF ; Olver et al. (2010). The Mathieu functions constitute a complete, orthogonal and normalized set of eigenfunctions: $\int_{0}^{2\pi}\mathrm{d}x\ \text{ce}_{2n}(q,x)\text{ce}_{2m}(q,x)=\delta_{nm}\pi$ and similarly for $\text{se}_{2n+2}(q,x)$. Hence, the full solution of Eq. (17) in terms of the eigenfunctions reads

[TABLE]

Since the eigenvalues are ordered with increasing $n$, $a_{0}(q)<b_{2}(q)<a_{2}(q)<b_{4}(q)<...$, the series expansion converges and can be evaluated numerically, which is discussed in more detail in Appendix A.

We impose boundary conditions, that reflect different experimental setups accounting for clamped or free orientations at the ends of the polymer. The presented partition sum [Eq. (26)] represents the case of a clamped polymer with given initial and final orientation. This can be complemented by two or more scenarios. First, we look at a half-clamped polymer, where the initial orientation is clamped and the final one is free. The corresponding partition sum $Z(L|\varphi,0)$ is obtained via integration over all final angles,

[TABLE]

where $A_{0}^{2n}(q)$ results from the integral over the even Mathieu functions and the odd Mathieu functions do not contribute anymore.

Further, integrating also over the initial angle the partition sum corresponding to a free polymer reduces to

[TABLE]

Interestingly, the same result occurs for the semiflexible polymer under tension Prasad et al. (2005), which reflects that the polymer is free to align with the external force. To evaluate the partition sums in Eqs. (26) and (27) numerically, we rely on an implementation of the Mathieu functions in a computer algebra system Wolfram Research (2016). However, to determine the Fourier coefficients $A_{0}^{2n}(q)$ in Eqs. (27) and (28) it is more efficient to solve numerically the eigenvalue problem of the recurrence relations [Eq. (23)-(25)].

II.2 Linear response

The linear response $\chi(T,F=0)$ of the semiflexible polymer to the compression force can also be obtained exactly by expanding the full solution in Eqs. (26), (27), and (28) to second order in the force. For a clamped polymer we use the expansions of the Mathieu functions for $\varphi_{0}=\varphi_{L}=0$ and the corresponding eigenvalues up to second order in the deformation parameter $q$ from Refs. DLMF ; Olver et al. (2010):

[TABLE]

Note that for $\varphi_{0}=\varphi_{L}=0$ the odd Mathieu functions do not contribute to the partition sum. Inserting these expansions into the full solution [Eq. (26)] one can derive both the projected mean end-to-end distance as well as the force-free susceptibility [Eq. (11)]. The explicit expressions are rather lengthy for arbitrary parameters $L$ and $\ell_{\text{p}}$ and will not be shown here, but we will discuss them for a clamped polymer in Sec. III, Fig. 3.

For the half-clamped polymer ($\varphi_{0}=0$) only the first three terms, $n=0,1,2,$ contribute to the order considered to the partition sum in Eq. (27). Here, we use the expansions up to second order in $q$ of the Fourier coefficients $A_{0}^{0}(q)=(1-q^{2}/16)/\sqrt{2}$, $A_{0}^{2}(q)=q/4$ and $A_{0}^{4}(q)=q^{2}/96$ obtained from Refs. DLMF ; Olver et al. (2010), and inserted these together with the expansions from Eqs.(29)-(33) into the full solution [Eq. (27)]. Thus, the partition sum reads

[TABLE]

which results in a projected mean end-to-end distance

[TABLE]

and linear susceptibility

[TABLE]

For the free polymer only the first two terms contribute and the partition sum in Eq. (28) is evaluated similar to the partition sum of a half-clamped polymer [Eq. (34)]. It reduces to

[TABLE]

Here the projected mean end-to-end distance vanishes for zero force after averaging over the initial and final angles,

[TABLE]

and the linear susceptibility reads

[TABLE]

II.3 Pseudodynamics

To validate our analytical results, we also simulate the semiflexible polymer under compression. To generate a set of representative conformations of the polymer in equilibrium, we rely on Brownian dynamics simulations that yield the canonical ensemble as stationary state. Since we are not interested in the dynamic properties here, features such as hydrodynamic interactions, anisotropic friction, and the overall translational motion are ignored. Thus the purpose of the pseudodynamics is merely to reproduce the equilibrium properties.

We first discretize equidistantly the contour $L$ of the polymer in terms of the positions of the beads $\{\mathbf{R}_{i}\}_{i=0}^{N}$ and corresponding tangent vectors $\{\mathbf{u}_{i}\}_{i=0}^{N-1}$, where $\mathbf{u}_{i}=(\mathbf{R}_{i+1}-\mathbf{R}_{i})N/L$ is of unit length, $|\mathbf{u}_{i}|=1$.

The pseudodynamics for the orientation of the $i$th bead is prescribed by the Langevin equation in the It$\bar{\text{o}}$ sense,

[TABLE]

for $i=1,\dotsc,N-2$ and appropriately modified Langevin equations at the boundary, $i=0,N-1$. Here $\mathbf{u}_{i}^{\perp}(t)$ is the unit orientation rotated clockwise by an angle of $\pi/2$, i.e. $\mathbf{u}_{i}^{\perp}(t)\cdot\mathbf{u}_{i}(t)=0$ and $\det[\left(\mathbf{u}_{i}(t),\mathbf{u}_{i}^{\perp}(t)\right)]=1$. Furthermore $\omega_{i}(t)$ is a Gaussian white noise process with zero mean and variance $\langle\omega_{i}(t)\omega_{j}(t^{\prime})\rangle=\delta_{ij}\delta(t-t^{\prime})$ for $i,j=0,\dotsc,N-1$. The scaled rotational diffusion of the orientations $\hat{D}_{\text{rot}}$ sets only the time scale in our case. We show in Appendix B that the Langevin equations [Eq. (40)] in fact reproduce the canonical ensemble of the semiflexible polymer problem.

The boundary conditions for a clamped polymer are set by fixed initial and final orientation chosen into the direction antiparallel to the compression force, in particular, $\mathbf{u}_{0}(t)=\mathbf{u}_{N-1}(t)=-\mathbf{e}$ for all times $t$. A half-clamped polymer fulfills the same boundary condition for the initial orientation, $\mathbf{u}_{0}(t)=-\mathbf{e}$, whereas the final orientation can freely rotate.

Simulation snapshots [Fig. 1] reveal that a semiflexible polymer shows qualitatively different buckling behavior with respect to the boundary conditions. Since clamped polymers cannot rotate their ends, they exhibit an S-shaped configuration under a compressive load. In contrast, if the second end of the polymer is free to rotate, it can align with the direction of the applied force and a hook-shape configuration is observed. Considering two free ends the polymer will just rotate under the force and thus its behavior is expected to be identical to that of a polymer under tension.

To obtain reliable statistics in equilibrated configurations we have conducted $120$ realizations of a polymer with $N=300$ segments and a time step of $10^{-5}/\hat{D}_{\text{rot}}$ over a time horizon of $10^{3}/\hat{D}_{\text{rot}}$. Measurements are taken after the polymer has reached equilibrium. For a free polymer the simulations have been performed over a longer time horizon of $10^{4}/\hat{D}_{\text{rot}}$, however, with fewer realizations, since the free rotation of the polymer requires long transients to reach equilibrium at low forces.

In addition we have performed Monte Carlo simulations, similar to Ref. Baczynski et al. (2007), to crosscheck the results. They confirm our analytical predictions for clamped polymers, but they start to fail in terms of convergence and statistics when refining the discretization to capture the semiflexibility of the polymers and also the simulation of a free polymer has been demanding.

II.4 Classical Euler buckling

A stiff rod under compression at small forces does not yield at all but starts to buckle at a critical force that is known as the Euler buckling instability Landau and Lifshitz (1986). The critical force is determined by the rigidity $\kappa$ and the contour length $L$ by

[TABLE]

where $\gamma$ accounts for the different boundary conditions. In the case of a clamped rod, $\gamma=1$, with one end fixed and one translationally free; however, if one end of the rod is clamped and one orientationally and translationally free, $\gamma=2$ Landau and Lifshitz (1986).

For small temperatures the Fokker-Planck equation [Eq. (17)] of the partition sum can be approximated for a stiff rod by an Eikonal approximation Kleinert (2009). Here we set the partition sum $Z=\exp(-G/k_{\text{B}}T)$, neglect terms of higher order in the inverse temperature, and consequently obtain the approximate equation for the Gibbs free energy,

[TABLE]

Using the method of characteristics leads to the equation of motion for the angle,

[TABLE]

which is as anticipated reminiscent of the equation of motion of a classical pendulum Aldrovandi and Ferreira (1980). The same equation results as Euler-Lagrange equation minimizing the total energy [Eq. (4)]. Starting from this equation, the force-extension relation for $|F|>F_{c}$ of a clamped and half-clamped rod is well known Landau and Lifshitz (1986); Baczynski et al. (2007),

[TABLE]

These results will be compared to the elastic properties of semiflexible polymers confined to two dimensions.

III Elastic properties

In this section we provide a discussion of the analytic solution and simulation results for the force-extension relation and the associated susceptibility.

*Clamped polymer: *The mean projected end-to-end distance $\langle X\rangle$ of a clamped polymer is obtained by numerical differentiation of the Gibbs free energy [Eq. (9)]. We find a smooth monotone crossover from an almost stretched to a buckled configuration with increasing compression force, in contrast to the classical Euler buckling instability, where the first yielding takes place at the critical force $F_{c}$ [Fig. 2 (a)]. Due to thermal fluctuations already for vanishing forces, semiflexible polymers display undulations, hence, their projected end-to-end distance is shorter than their contour length, $\langle X\rangle<L$. As anticipated, in the regime of small forces higher flexibility leads to shorter mean end-to-end distances, which approach zero for completely flexible polymers, $\ell_{\text{p}}/L\ll 1$. Qualitatively, the linear response behavior appears to be correct up to forces comparable to the classical Euler buckling force (see Figs. 2  (a) and (b)). Interestingly, for increasing forces, $|F|>F_{c}$, the end-to-end distance of semiflexible polymers intersects with that of a stiff rod, indicating a stiffening of the polymer due to thermal fluctuations. A similar trend occurs when comparing different persistence lengths, hence, initially more flexible polymers yield more strongly while for large forces they resist harder than stiffer polymers.

Furthermore we compare the susceptibility $\chi$ of semiflexible polymers in response to the compression force to the susceptibility of a stiff rod [Fig. 2 (d)]. The latter vanishes for forces smaller than the critical Euler buckling force, reflecting that a classical rod does not yield. Directly at the transition the classical susceptibility assumes the value $\chi(F_{c}+0)=2L/F_{c}$ and decreases monotonically for larger forces. In contrast the susceptibility of semiflexible polymers initially increases starting from $\chi(T,F=0)$ as predicted from the linear response (not shown), then displays a maximum at a force $|F_{\text{max}}|$, which exceeds the Euler buckling force, $|F_{\text{max}}|>F_{c}$. Increasing the stiffness, this force approaches the critical buckling force $|F_{\text{max}}|\downarrow F_{c}$, which suggests to use $|F_{\text{max}}|$ as a measure to charactzerize the smooth transition for fluctuating semiflexible polymers. For simplicity, we refer to $|F_{\text{max}}|$ as the buckling force also in the case of non-zero temperature. Interestingly, for stronger forces the projected end-to-end distance becomes even negative (see Fig. 3), which reflects the S-shaped configuration in Fig. 1.

Half-clamped polymer: Here, the force-extension relation behaves similarly to the one for clamped polymers (see Fig. 2 (b)). Nevertheless half-clamped polymers yield more strongly even after taking into account that the classical buckling force is reduced by a factor of $\gamma^{2}=4$ due to the different boundary conditions. In fact, the data demonstrate that the semiflexibility becomes more important for half-clamped polymers.

Qualitatively, the monotonous transition from the straightened to the buckled configuration of these polymers remains unchanged, intersections with the end-to-end distance of a stiff rod persist also here. The corresponding susceptibilities [Fig. 2 (e)] are less pronounced in comparison to clamped polymers of the same stiffness [Fig. 2 (d)]. Also, here we characterize the transition by a buckling force $|F_{\text{max}}|$ taken as the maximum of the susceptibility. For half-clamped ends the approach $|F_{\text{max}}|\downarrow F_{c}$ is slower than for clamped ends.

Free polymer: In contrast to clamped polymers, the force-extension relation of a free polymer exhibits qualitatively different behavior (see Fig. 2 (c)). Direct inspection of the configurations (not shown) reveals that the polymer essentially rotates and aligns with the direction of the force. Thus, rather than compressing, the force in fact stretches the polymer in the reverse direction, in particular, the projected end-to-end distances are all negative. The analytic solution for the case of pulling has been provided before Prasad et al. (2005), however, only an approximate solution for rather flexible polymers of length $L>4\ell_{\text{p}}$ has been discussed in more detail. Thermal undulations prevent the polymer to be in a completely straight configuration. Stiffer polymers react to smaller forces than more flexible polymers, and extend nearly to their full contour length. For free ends the force-extension relations approach monotonically the case of a fully aligned rod as the stiffness increases, in contrast to the other two boundary conditions. Similarly, the susceptibility [Fig. 2 (f)], is always maximal at vanishing forces, where it assumes the value $\chi(T,F=0)$ as predicted by linear response theory, and is larger for increasing stiffness. Correspondingly, no characteristic buckling force is extracted for this case, since after all the polymer stretches rather than buckles in this case.

Approximate Solution: An approximate solution of the force-extension relation in the regime of rather stiff polymers has been presented by Baczynski et al. Baczynski et al. (2007) by integrating out short wavelength fluctuations. Their results predict qualitatively similar behavior for clamped polymers in terms of the smooth buckling transition. Yet, already in the regime for small forces $|F|\lesssim F_{c}$ deviations to our analytic theory become apparent (see Fig. 3). The differences decrease as stiffer polymers are considered. Their approximate solution for stiff polymers exhibits far weaker buckling than the predicted analytic solution and in addition the force-extension relation remains positive for all forces. This finding is not in agreement with the analytic solution for semiflexible polymers as well as with the behavior of a classical stiff rod, which both predict a negative projected end-to-end distance for large forces, $|F|\gg F_{c}$. Hence, the effect, that the polymer assumes and S-shaped configuration in response to the applied force is not captured within the approximate theory.

For both the analytic and the approximate solution, the force-extension curves of a stiff rod and the semiflexible polymers intersect. The corresponding intersection force of the force-extension has been computed in Ref. Baczynski et al. (2007) and confirms that this intersection takes place in the vicinity of the Euler-buckling instability. In contrast to the critical buckling force $|F_{\text{max}}|$ defined by the maximal susceptibility, the critical buckling force $|F_{\text{ap}}^{c}|$ of Ref. Baczynski et al. (2007) has been defined as the onset of the buckling of the polymer, shown in Fig. 3. Similar to our analysis, it exceeds the critical Euler buckling force, but always remains smaller than the maximal buckling force $|F_{\text{ap}}^{c}|<|F_{\text{max}}|$. Furthermore, the susceptibility extracted from the approximate theory deviates significantly from our analytic solution, yet, the force at the maximal susceptibility remains similar (see Fig. 2(d)).

*Pulled polymer: *We have also evaluated numerically for the first time the force-extension relation of a pulled polymer. The case of free ends under pulling is up to a sign identical to compression, as discussed above (see Fig. 2(c)). In Fig. 4 we evaluate the analytic solution provided by Ref. Prasad et al. (2005) and compare it to the weakly-bending approximations Smith et al. (1992); Marko and Siggia (1995) for a clamped polymer,

[TABLE]

and a half-clamped polymer,

[TABLE]

see Appendix C for a derivation. These approximations accurately reproduce the force-extension relation of pulled stiff polymers. However, the weakly-bending assumption is not fulfilled for more flexible polymers at small forces, which display deviations from the exact solution (see Fig. 4). The force-extension relation of a (half-) clamped polymer exhibits qualitatively similar behavior to that of a free polymer, and aligns along the direction of the applied force.

IV Thermodynamic properties

In addition to the elastic behavior, we discuss the thermodynamic properties in terms of the excess entropy

[TABLE]

and the excess heat capacity $\Delta C_{F}=T(\partial\Delta S/\partial T)_{F}$ in terms of our exact solution.

For clamped polymers, we observe an increase of the excess entropy $\Delta S$ at small forces and a decrease at forces exceeding the Euler buckling force $|F|\gtrsim F_{c}$ (see Fig. 5 (a)). Thus, by applying a small force the number of accessible configurations of the polymer is increased, which reflects that by the S-shaped configurations the fluctuations in the direction perpendicular to the force become more important. Yet, for strong forces the excess entropy decreases again and eventually even becomes negative. In fact, the simulations show that for these strong forces the clamping suppresses the undulations of the contour.

The temperature dependence of the excess entropy is encoded in the heat capacity $C_{F}(T,F)$. We find that the excess heat capacity is negative for forces $|F|\lesssim F_{c}$, whereas it becomes positive and displays a maximum for larger forces; see Fig. 5 (a) (inset). Thus, the initial increase of fluctuations upon compression is less important for more flexible polymers, while the suppression for large forces is more relevant for stiffer polymers.

Qualitatively similar results are obtained for half-clamped polymers both for the excess entropy and the excess heat capacity (not shown).

In contrast, free polymers behave qualitatively differently in terms of their thermodynamic properties. In particular, the excess entropy $\Delta S$ always decreases monotonically with increasing force (see Fig. 5 (b)). Here the polymer reverses direction and is effectively under tension, such that the force suppresses the thermal undulations and the polymer straightens out. The corresponding excess heat capacity $\Delta C_{F}$ increases and displays a prominent maximum before approaching a constant value for very large forces. As corroborated by the simulations this reflects that for more flexible polymers the thermal fluctuations are less suppressed than for stiffer ones (see Fig. 5 (b) (inset)).

V Summary and Conclusion

We have derived an analytic expression for the partition sum of a semiflexible polymer under compression within the framework of the WLC model for arbitrary stiffnesses. The partition sum is represented as an expansion in terms of Mathieu functions, which are the eigenfunctions of the associated Fokker-Planck equation. Since the eigenvalues form an ascending series to infinity, the infinite sum of decaying exponentials can be evaluated numerically (see Appendix A). Elastic properties such as the force-extension relation, the susceptibility, and thermodynamic properties in terms of the excess entropy and excess heat capacity of the polymers with different persistence lengths, have been obtained as derivatives of the exact Gibbs free energy. These properties have been discussed for different boundary conditions reflecting various experimental setups.

Our results for the force-extension relation of (half-) clamped polymers predict a smooth crossover from an equilibrated, almost straight rod to a buckled configuration. Interestingly, the force-extension curves intersect in the vicinity of the classical Euler buckling instability. Therefore, stiffer polymers resist more strongly for small forces, yet at forces larger than the buckling force they yield more strongly than flexible ones. Furthermore, the susceptibilities display a maximum in the vicinity of the classical Euler buckling force and by the fluctuation-response theorem this implies that the fluctuations of the projected mean end-to-end distance become most important here. Our analysis reveals that for large forces the weakly-bending approximation breaks down, as is apparent already for a classical elastic rod. Therefore the weakly-bending approximation remains valid only for large persistence lengths $\ell_{\text{p}}\gg L$ but is also restricted to small compression forces $|F|\lesssim F_{c}$.

In contrast, a free polymer cannot resist a compression force, rather it rotates and aligns with the applied force. Therefore, our solution for this case reproduces the well-studied setup of a semiflexible polymer under tension.

Apart from the weakly-bending approximation Baczynski et al. (2007), also different approximation schemes have been elaborated earlier to understand the behavior of rather stiff polymers close to the buckling transition. These studies have incorporated corrections by thermal fluctuation to the theory of the classical Euler buckling instability Lee et al. (2007); Emanuel et al. (2007), which capture quantitatively the force-extension relation for forces much smaller and much larger than the critical Euler buckling force $F_{c}$ and compare it to simulations of semiflexible polymers modeled by a bead-spring chain Emanuel et al. (2007). Including the lowest-order quartic mode in the fluctuations predict qualitatively the same behavior as our analytical theory for planar inextensible polymers close to the critical Euler buckling instability Lee et al. (2007). In addition, the buckling behavior of extensible polymers is also predicted to delay the buckling transition within the regime of small thermal fluctuations Bedi and Mao (2015). Alternatively to existing literature, one can also start from the probability distribution of a semiflexible polymer, which has been elaborated in the weakly-bending regime Wilhelm and Frey (1996) and within a mean-field approach for filament inextensibility Blundell and Terentjev (2009a), to compute approximately the force-extension relation.

However, these approximations Wilhelm and Frey (1996); Baczynski et al. (2007); Lee et al. (2007); Emanuel et al. (2007); Bedi and Mao (2015) remain valid for stiff polymers or away from the buckling transition only, as thermal fluctuations increase drastically for more flexible polymers at the maximum of the susceptibility. In contrast, our analytical theory provides the force-extension relation for the full range of semiflexible polymers.

Also for pulled semiflexible polymers a lot of work has been done in elaborating the force-extension relation Meng and Terentjev (2017), ranging from weakly-bending approximations Smith et al. (1992); Marko and Siggia (1995) to various approaches to account for the inextensibility of a wormlike chain MacKintosh et al. (1995); Ha and Thirumalai (1997); Palmer and Boyce (2008); Blundell and Terentjev (2009b). Here, we have evaluated for the first time the force-extension relation of pulled polymers using the analytical solution provided by Ref. Prasad et al. (2005), where up to now only approximations of the fore-extension relations of rather flexible polymers have been discussed. As we have shown in Appendix A, for flexible polymers the numerical evaluation of the partition sum reduces to a few terms, and can therefore be approximated due to the properties of the eigenvalues by the first term only, as has been observed earlier Prasad et al. (2005). However, for increasing stiffness of the polymer more terms contribute to the partition sum, the approximation fails and the full solution is required to accurately reproduce these elastic properties. In the limiting case of stiff polymers, we find excellent agreement with the weakly-bending approximation for all forces (see Fig. 4).

We anticipate, that our results serve as a reference for experimentally measured force-extension relations, and allow determining the persistence length for a broad range of semiflexible polymers. We have provided tables (see Supplemental Material [SeeSupplementalMaterialat\url{http://link.aps.org/supplemental/10.1103/PhysRevE.95.052501}fortablesoftheforce-extensionrelations]supplement) containing the numerical evaluation of the force-extension relation for semiflexible polymers subject to compression as well as pulling forces, in order to make the developed theory accessible and applicable to experimental observations.

Mathematically, the partition sum can be viewed as the Laplace transform of the distribution of the projected end-to-end distance. Therefore, the distribution can be obtained in principle by an inverse Laplace transform; however, this is numerically an ill-posed problem. Consequently, one should rather consider the characteristic function, i.e., the Fourier transform of the distribution, which is obtained formally by analytic continuation in the force to complex values. The Mathieu functions and the eigenvalues in this case can still be obtained numerically by solving a matrix eigenvalue problem. This approach has been applied recently for the mathematical analog of a self-propelled particle Kurzthaler et al. (2016).

Our solution method is not restricted to the plane but can be extended to three dimensions. There the Mathieu functions of the eigenvalue problem are replaced by generalized spheroidal wave functions Leitmann et al. (2016); Kurzthaler et al. (2016). Similarly, our approach can be extended in principle to account for a spontaneous curvature such that the classical reference system consists of a circular arc.

The equilibrium single-polymer behavior under tension or compression plays a crucial role as input for elastic properties of networks, and the regime of strong compression close to buckling may be useful to characterize the mechanical stability of such networks MacKintosh et al. (1995); Kroy and Frey (1996); Storm et al. (2005); Claessens et al. (2006); Chaudhuri et al. (2007); Huisman et al. (2008); Carrillo et al. (2013); Razbin et al. (2015); Amuasi et al. (2015); Plagge et al. (2016). In particular, networks are often modeled as entangled solutions of single semiflexible polymers with entanglement points, where they cross or loop each other, branching points or cross-links MacKintosh et al. (1995); Kroy and Frey (1996); Storm et al. (2005); Huisman et al. (2008); Carrillo et al. (2013); Razbin et al. (2015); Amuasi et al. (2015); Plagge et al. (2016). Already in equilibrium, single polymers experience a stretching or compression force induced by the surrounding network, and consequently the mean distance between, for example, two cross-links differs from the contour length of the polymers. Here, our theoretical predictions for the elastic properties of semiflexible polymers subject to different boundary conditions, which depend on the properties of the entanglement points or cross-links, allow adequate modeling and predict the equilibrium properties of polymer networks. Furthermore, exposing the network to mechanical stresses strongly depends on the force-extension relation of these individual filaments and crucially determines their elastic response. Up to now only approximate force-extension relations have been applied in the analysis of networks MacKintosh et al. (1995); Kroy and Frey (1996); Storm et al. (2005); Huisman et al. (2008); Carrillo et al. (2013); Amuasi et al. (2015), whereas our characterization of (half-)clamped and free semiflexible polymers allows analyzing the elasticity of inhomogeneous networks composed of filaments with arbitrary stiffnesses. In particular, the clamped boundary conditions might serve as a starting point for modeling polymers connected by crosslinks, where their orientation is fixed due to the specific bonds. Differently, weaker entanglements can be accounted for by applying half-clamped boundary conditions, where one end of the polymer is allowed to rotate freely. Hence, the theoretical predictions can be incorporated by using the tables provided in the Supplemental Material sup in further modeling studies of networks, which would not alter the computational cost of simulations significantly.

Similarly, the equilibrium behavior in free space constitutes the reference for the relevant case of semiflexible polymers immersed in a dense crowded medium Schöbl et al. (2014); Keshavarz et al. (2016); Leitmann et al. (2016). Furthermore, the equilibrium properties should serve as the basis also for dynamic studies Hallatschek et al. (2007a, b); Obermayer et al. (2009), in particular, the pseudodynamics we have derived here constitutes a convenient starting point to investigate the relaxation of undulation modes in the regime where the weakly-bending approximation is no longer valid.

Acknowledgements.

We thank Sebastian Leitmann and Klaus Kroy for many discussions, and Sebastian Leitmann for a critical reading of the manuscript. This work has been supported by Deutsche Forschungsgemeinschaft (DFG) via Contract No. FR1418/5-1.

Appendix A Numerical evaluation

Here, we discuss the numerical evaluation of the partition sums in Eqs. (26), (27), and (28), which are sums of relaxing exponentials with respect to the length $L$. Yet, the coefficients, i.e. the Mathieu functions ($\text{ce}_{2n}(q,x)$ and $\text{se}_{2n+2}(q,x)$), as well as the exponents, i.e. the eigenvalues of the Mathieu functions ($a_{2n}(q)$ and $b_{2n+2}(q)$), are non-monotonic functions with respect to the force $F$; see Fig. 6.

In particular, we find, that the zeroth eigenvalue, $a_{0}(q)$, remains negative for all forces, whereas the second, $a_{2}(q)$, and fourth eigenvalues, $a_{4}(q)$, change sign within the parameter range considered. The other eigenvalues remain positive and increase with increasing mode $n$ for all relevant forces (see Fig. 6 (b)). Consequently, the higher modes of the partition sum become exponentially suppressed, and therefore induce a natural cut-off of the infinite series [Eqs. (26), (27), and (28)].

Moreover, with increasing mode $n$ the Mathieu functions become more important at intermediate and larger forces (see Fig. 6 (a)), and therefore many modes are expected to contribute to the partition sum. Interestingly, we observe that the number of modes necessary to achieve a desired accuracy of the elastic properties highly depends on the force.

As an example, we consider the mean end-to-end distance $\langle X\rangle$, and determine the number of modes necessary to obtain the accuracy $|\langle X\rangle_{n+1}-\langle X\rangle_{n}|/L<10^{-8}$, where $\langle X\rangle_{n}$ denotes the mean end-to-end distance with $n$ included modes. For a clamped and a half-clamped polymer, the number of modes increases and assumes a maximum in the vicinity of the critical Euler buckling force, and decreases again for larger forces; see Figs. 7(a) and (b). Similarly, for a free polymer the number of modes also increases, yet at smaller forces, and decreases for large forces [Fig. 7(c)].

In addition, the number of modes rises drastically for stiffer polymers, which makes numerical evaluation more tedious. Furthermore, the number of modes for half-clamped boundary conditions is largest, whereas about half the modes are required for a free polymer.

Appendix B Derivation of the pseudo dynamics

In general, equations of motion in terms of stochastic differential equations can be derived from the Fokker-Planck equation of the stochastic process Gardiner (2009). At this stage, only the equilibrium distribution of the polymer subject to an external compression force is known and an equation of motion for the probability density needs to be formulated. Therefore, we construct a Fokker-Planck equation starting with the requirement that the known equilibrium distribution of the WLC model coincides with the stationary distribution of the stochastic process. Therefore, we discretize the polymer equidistantly in terms of the positions of the beads $\{\mathbf{R}_{i}\}_{i=0}^{N}$ and corresponding tangent vectors $\{\mathbf{u}_{i}\}_{i=0}^{N-1}$, where $\mathbf{u}_{i}=(\mathbf{R}_{i+1}-\mathbf{R}_{i})N/L$ with unit length $|\mathbf{u}_{i}|=1$. The discretized Hamiltonian takes the form

[TABLE]

where $\hat{\ell}_{\text{p}}=\ell_{\text{p}}N/L$ and $\hat{f}=fL/N$ denote the scaled persistence length and force ($\hat{f}>0$ for tension and $\hat{f}<0$ for compression), respectively.

To fulfill the inextensibility constraint, we parametrize the orientation $\mathbf{u}_{i}=(\cos\varphi_{i},\sin\varphi_{i})^{T}$ in terms of the polar coordinates $\{\varphi_{i}\}\equiv\{\varphi_{i}\}_{i=0}^{N-1}$ measured relative to the fixed direction $\mathbf{e}$ of the force. Therefore, the probability density in equilibrium is given by

[TABLE]

such that $\int\left[\prod_{i=0}^{N-1}\mathrm{d}\varphi_{i}\right]\mathbb{P}_{\text{eq}}(\{\varphi_{i}\})=1$.

The discretized Hamiltonian in polar coordinates reads

[TABLE]

Next we derive the Fokker-Planck equation describing the time evolution of the conditional probability density $\mathbb{P}\equiv\mathbb{P}(\{\varphi_{i}\},t|\{\varphi_{i}^{0}\})$ to find a discretized polymer with orientations $\{\varphi_{i}\}$ at time $t$, given that it was oriented with $\{\varphi_{i}^{0}\}$ at time $t=0$. It is expressed by

[TABLE]

where the velocity of the probability current is obtained by the friction law, $U^{i}(\{\varphi_{i}\})=\sum_{k}K^{ik}(\{\varphi_{i}\})F_{k}(\{\varphi_{i}\})$. We choose for simplicity $\mathbf{K}=N\xi_{r}^{-1}\mathbb{I}$, which denotes the rotational mobility tensor and $\mathbf{F}$ the forces,

[TABLE]

Here, the first term corresponds to the mechanical forces while the second accounts for the thermal Brownian forces. Thus one verifies that the equilibrium distribution $\mathbb{P}_{\text{eq}}$ corresponds to the stationary distribution, $\partial_{t}\mathbb{P}_{\text{eq}}=0$. Collecting results we find for the time evolution of the conditional probability density,

[TABLE]

with scaled rotational diffusion coefficient $\hat{D}_{\text{rot}}=Nk_{\text{B}}T/\xi_{\text{r}}$ and initial condition

[TABLE]

Starting from the Fokker-Planck equation [Eq. 53] standard methods Gardiner (2009) are used to obtain the Langevin equations for the angles, governing the pseudodynamics of a semiflexible, discretized polymer,

[TABLE]

where $\mathrm{d}\omega_{i}(t)$ is the increment of a Gaussian white noise process $\omega_{i}(t)$ with zero mean $\langle\omega_{i}(t)\rangle=0$ and delta-correlated variance $\langle\omega_{i}(t)\omega_{j}(s)\rangle=\delta_{ij}\delta(t-s)$ for $i,j=0,\ldots,N-1$. These equations can be transformed by It$\bar{\text{o}}$’s lemma Gardiner (2009) to the equations for the tangent vectors in Eq. (40).

Appendix C Force-extension relation for a pulled polymer in the weakly-bending regime

In the weakly-bending regime of a pulled polymer, $f=|f|$, we can approximate Eq. (17) by

[TABLE]

and use the Gaussian ansatz

[TABLE]

Then the inverse variance $M(s)$ and the normalization $\Gamma(s)$ have to fulfill the equations of motion

[TABLE]

Using the initial condition $Z(\varphi,s=0|\varphi_{0}=0,0)=\delta(\varphi)$, the solutions of these differential equations read

[TABLE]

To obtain the partition sum for a half-clamped polymer with contour length $L$, we average over the final orientation,

[TABLE]

which is correct if the angular fluctuations are small, $M(s)\gg 1$. This is fulfilled for stiff polymers $L/\ell_{\text{p}}\ll 1$ or strong pulling $\ell_{\text{p}}|f|\gg 1$. Finally, we obtain the force-extension relations for (half-) clamped polymers in the weakly-bending regime by taking the derivative of the Gibbs free energy [Eq. (7)] with respect to the force (see Sec. II).