Semiflexible polymers: Dependence on ensemble and boundary orientations Debasish Chaudhuri Department of Biological Physics, Max Planck Institute for the Physics of Complex Systems, N¨othnitzer Str. 38, 01187 Dresden, Germany ∗ (Dated: February 6, 2008) 7 0 Weshowthatthemechanicalpropertiesofaworm-like-chain(WLC)polymer,ofcontourlengthL 0 andpersistencelengthλsuchthatt=L/λ∼O(1),dependbothontheensembleandtheconstraint 2 on end-orientations. In the Helmholtz ensemble, multiple minima in the free energy near t = 4 n persists for all kinds of orientational boundary conditions. The qualitative features of projected a probability distribution of end to end vector depend crucially on the embedding dimensions. A J mapping of the WLC model, to a quantum particle moving on the surface of an unit sphere, is 8 used to obtain the statistical and mechanical properties of the polymer under various boundary 1 conditions and ensembles. The results show excellent agreement with Monte-Carlo simulations. ] PACSnumbers: 87.15.La,36.20.Ey,87.15.Ya h c e I. INTRODUCTION For instance, microtubule-associatedproteins attachone m or both of their ends to microtubules to arrange them - in microtubule bundles [4]. Again, in gene-regulation at Microtubules and actin polymers constitute the struc- often DNA-binding proteins loop DNA with fixed end t tureofcytoskeletonthatgivesshape,strengthandmotil- s orientations [5, 6, 7]. Thus it becomes important to un- . itytomostofthelivingcells. Theyaresemiflexiblepoly- derstand the statistics and the mechanical properties of t a mers in the sense that their persistence lengths λ are of semiflexible polymers with different possibilities of end m the order of their chain lengths (the statistical contour orientations and ensembles. - lengths) L such that the stiffness parameter t = L/λ is d small and finite. For example, actin, microtubule and During the last decade many single molecule experi- n ments have been performed on semiflexible polymers[3, double stranded DNA (dsDNA) have λ=16.7 µm[1, 2], o 8, 9, 10]. These have been done by using the optical 5.2mm[2]and50nm[3]respectively. Inphysiologicalsit- c tweezers[9],themagnetictweezers[11]andtheAFMs[12]. [ uationLofdsDNAinsideacellmayvaryinbetweenmil- limeters to a meter with average length in human being In the optical tweezer experiments one end of a polymer 2 isattachedtoadielectricbeadwhichis,inturn,trapped 5cm, whereas typical contour lengths of microtubules v ∼ by the light intensity profile of a laser tweezer. In this can be 10µm[4]. The contour length of actin filament 8 ∼ casethedielectricbeadisfreetorotatewithintheoptical canbe aslargeas 100µm[1]. Inthe in vitro experiments, 9 trap. On the other hand, attaching an end of a polymer 2 thecontourlengthsofbio-polymerscanbetailoredchem- toasuper-paramagneticbead,onecanusemagneticfield 8 ically;e.g. intheexperimentdescribedinRef.[1]thecon- gradients to trap the polymer using a magnetic tweezer 0 tour lengths of actin polymer have a distribution up to 6 L = 30µm. For a polyelectrolyte like DNA the persis- setup. Inthis caseonecanrotatethe beadwhile holding 0 itfixedinpositionbychangingthedirectionoftheexter- tence length λ canalso be tuned a little by changing the / nal magnetic field. In the AFM experiments one end of t salt concentration of the medium. The relevant parame- a a polymer is trapped by a functionalized tip of an AFM ter in deciding the mechanical properties is the stiffness m cantilever. Thetwodistinctprocedureswhichcanbefol- parameter t, contour length measured in units of persis- - lowed to measure force-extension are: (a) Both the ends d tence length. While it is obvious that in the thermody- ofthepolymerareheldviathelaserormagnetictweezers n namiclimitoft ,theGibbs(constantforce)andthe →∞ or the AFMs. (b) One end of the polymer is attachedto o Helmholtz (constant extension) ensemble predict identi- c cal properties, the same is not true for real semiflexible asubstratesuchthatthe positionandorientationofthis : end is fixed while the other end is trapped via a laser or v polymerswhicharefarawayfromthislimit. Inbiological magnetic tweezer or an AFM cantilever. i cellsactinfilamentsremaindispersedthroughoutthe cy- X toplasm with higher concentration in the cortex region, While the optical tweezers allow free rotation of di- r just beneath the plasma membrane. microtubules, on electric beads within the trap, thereby, allowing free ori- a the otherhand, haveone end attachedto a microtubule- entations of the polymer end, the magnetic tweezers fix organizingcentre,centrosome,in animalcells. Thus bio- the orientation of the ends and one can study the de- logicallyimportantpolymersmayfloatfreelyormayhave pendence of polymer properties on end-orientations by oneoftheirendsfixed. Eventheendorientationsofpoly- controlled change of the direction of external magnetic mers play a crucial role in many important phenomena. field. In this paper, we call this fixing of orientation of an end of a polymer as grafting. By changing the trap- ping potential fromstiff to soft trapone can gofrom the HelmholtztotheGibbsensemble[13]. Beforeweproceed, ∗Electronicaddress: [email protected] let us first elaborate onhow to fix the ensemble of a me- 2 chanical measurement[13, 14]. In the simplest case we tity that describes statistical property of such polymers can assume that one end of the polymer is trapped in a is the distribution of end-to-end separation. Numerical harmonic well, V(z) = C(z z )2/2 with (0,0,z ) be- simulationsto obtainradialdistributionfunction for dif- 0 0 − ing the position ofthe potential-minimum. The polymer ferent values of t have been reported along with a se- end will undergo continuous thermal motion. One can ries expansion valid in the small t limit[30]. Mean-field use a feedback circuit to shift z to force back the fluc- treatments to incorporate the inextensibilty in an ap- 0 tuating polymer end to its original position. This will proximate way have also been reported[31, 32]. In an ensure a Helmholtz ensemble. This can also be achieved earlier study[16] we investigated the free energy profile by taking C . On the other hand, one can use a of a semiflexible polymer whose ends were free to rotate → ∞ feedback circuit to fix the force C(z z ) by varying in the constant extension ensemble and in the stiffness 0 − − C depending on the position z of the polymer end. This regime of 1 t 10. This work predicted that a clear ≤ ≤ will ensure a Gibbs ensemble. This can also be achieved qualitative signature of semiflexibility would be a non- by taking a vanishingly soft (C 0) trap to infinitely monotonic force extension for stiffnesses around t 4 → ∼ largedistance(z )suchthatwithinthelengthscale in the Helmholtz ensemble. This comes from the multi- 0 →∞ of fluctuation the polymer end feels a constant slope of modality of probability distribution of end to end sepa- the parabolic potential. Surely, in experiments, using a ration. However, this non-monotonicity is absent in the feedback circuit is easier to implement a particular en- Gibbs ensemble[16]. Multiple maxima in the probabil- semble. However,the other procedure is mathematically ity distribution of end to end separation was due to a well defined and one can seek recourseof it to show that competition between entropy, that prefers a maximum the partition function of two ensembles are related by a near zero separation, and energy, that likes an extended Laplacetransform[15]. Thisrelationdoesnotdependon polymer. Aseriesoflaterstudies[15,33,34,35]usedan- the choice of the Hamiltonian for a polymer. An exact alytic techniques to understand the end to end distribu- relation between the two ensembles for worm like chain tionatallstiffnessesincluding the stiffness regimewhere (WLC) model is shown in Sec.II. multimodality was observed. Recently, multimodality is found in transversefluctuations of a grafted polymer us- From the above discussion on possible experiments, it ing simulations [17] and approximate theory[36, 37]. A is clear that there can be three possibilities of boundary Greensfunctiontechniquehasbeendevelopedthattakes conditions in terms of orientations: (a) Free end: Both into account the orientations of the polymer ends [38]. the ends of a polymer can remain free to rotate[15, 16]. The impact of the specific boundary conditions and the (b) One end grafted: Orientationof one end is fixed and comparable length scales of a dsDNA and the beads to the other can take all possible orientations[17]. (c) Both which it is attached in typical force-extension measure- ends grafted: Orientations of both the ends are kept ments have been identified in another recent study [18]. fixed. Thus, in experiments, one can have two possible The WLC model has also been extended to study statis- ensembles and three possible boundary conditions. We tics of end to end separation and loop formation prob- restrictourselvestotheWLCpolymersembeddedintwo ability in dsDNA[39] and to incorporate twist degree of dimensions(2D).Weinvestigatetheprobabilitydistribu- freedom [40, 41, 42, 43]. tion, free energy profile and force extension relation for each of these cases in this paper. We shall see that the The construction of this paper is as follows. In Sec.II properties of a semiflexible polymer depend both on the we present a theoretical technique for exact calculation choiceoftheensembleandtheboundarycondition. Note of the WLC model via a mapping to a quantum parti- that, there can be other possibilities of boundary condi- cle moving on the surface of an unit sphere. This tech- tionse.g. orientationatoneendofapolymercanbe free nique incorporates all the possible end orientations and to rotate ona half-sphere[18]. However,in this paper we predictsresultsinboththe Helmholtz andthe Gibbsen- focus on the three possible boundary orientations listed sembles. In Sec.III we discuss the different discretized above. versions of the WLC model and the Monte-Carlo (MC) simulation procedures followed in this work. In Sec.IV TheWLCmodelisasimplecoarsegrainedwaytocap- wepresentallthe resultsofprobabilitydistributionsand turebendingrigidityofanunstretchablepolymer[19,20] force-extensions etc. obtained from theory and simula- embedded in a thermal environment. Recent single tions. Then, in Sec.V, we summarize our results and molecule experiments in biological physics [3, 8, 9, 10] conclude with some discussions. renewedinterestin this oldmodel ofpolymer physics. It wassuccessfullyemployed [21,22]tomodeldataofforce- extensionexperiments[8]ondsDNA.Mechanicalproper- ties of giant muscle protein titin [23, 24], polysaccharide II. THEORY dextrane [12, 24] and single molecule of xanthane [25] werealsoexplainedusingtheWLCmodel. Duetothein- extensibility constraint, the WLC model is hardto tract In the WLC model a polymer is takenas a continuous analytically except for in the two limits of flexible chain curve denoted by a d-dimensional vector ~r(s) where s is (t )andrigidrod(t 0),aboutwhichperturbative a distance measured over the contour of the curve from →∞ → calculationshavebeendone[26,27,28,29]. Akeyquan- one of its ends. This curve has a bending rigidity and 3 thus the Hamiltonian is given by This immediately gives, = 1/Z˜(~0). We show that N Z˜(~0) is a constant which depends on the constraints on κ L ∂tˆ(s) 2 end orientations. Eq.4 gives, β = ds , (1) H 2 ∂s wmherereistˆ(isn)ex=te∂n~rs(ibs)le/∂is.ei.sZt0tˆh2et=an1(cid:18)g,enβtvise(cid:19)ctthoerainnvdetrhseeptoemly-- P˜(f) = N Ztˆitˆf D[tˆ(s)]e„−(d−41)λR0Lds“∂t∂ˆ(ss)”2+λf R0Ltˆxds« perature. Persistence length is a measure of the dis- = tˆf [tˆ(τ′)]e»−R0t(d−41)“∂∂tˆ(ττ′′)”2−ftˆxffdτ′– tance up to which the consecutive tangent vectors on N Ztˆi D the contour do not bend appreciably and is defined by (5) tˆ(s).tˆ(0) = exp( s/λ). The bending rigidity κ is re- hlated to piersistenc−e length λ via κ=(d 1)λ/2. Thelaststepisobtainedbyreplacingτ′ =s/λandusing − Inthissectionwepresentatheoreticalmethodtosolve the identities κ = (d 1)λ/2 and t = L/λ. Note that, − the WLC model to any desired accuracy[15, 44] for both P˜(f), is the partition function, apart from a multiplica- theHelmholtzandtheGibbsensemblesandallthethree tiveconstant,intheGibbsensemblewheretbehaveslike possible boundary orientations over the entire range of an inverse temperature such that the Gibbs free energy stiffness parameter t. We first present the method for a canbewrittenasG(f)= 1/t lnP˜(f). Nowconsidering free polymer[15]. Then we extend it to calculate proper- τ′ as imaginary time and−replacing τ = iτ′ one gets, − ties of grafted [ one/both end(s) ] polymers. The partition function of a WLC polymer in the P˜(f)= tˆf [tˆ(τ)]e[i 0−itLdτ] ; (6) Helmholtz ensemble is Z(~r) = cexp(−βH) where c N Ztˆi D R denotes a sum over all possible configurations of the polymer that are consistent withPthe inextensibility con- withthe identificationof = (d−1) ∂tˆ(τ) 2+ftˆ asthe straint. The probability distribution of the end to end L 4 ∂τ x vector becomes, P(~r) = Z(~r)/ Ld~rZ(~r) = Z(~r). If Lagrangian,P˜(f)[=Z˜(f)/Z˜(0)]in(cid:16)theab(cid:17)oveexpression thetangentvectorsofthetwoendsofapolymNerareheld is the path integral representation for the propagator of fixed at tˆ and tˆ , the probabiliRty distribution of end to a quantum particle, on the surface of a d- dimensional i f endvectorinconstantextensionensemblecanbewritten sphere, that takes a state tˆi to tˆf . In Schrodinger | i | i in path integral notation as picturethiscanbewrittenastheinnerproductofastate tˆ andanotherstate tˆ evolvedbyimaginarytime it, i f | i | i − tˆf L P(~r)= [tˆ(s)]exp( β ) δd ~r tˆds (2) Z˜(f)= tˆ exp( iHˆ( it))tˆ = tˆ exp( tHˆ)tˆ , (7) i f i f N Ztˆi D − H × −Z0 ! h | − − | i h | − | i where Hˆ is the Hamiltonian operator corresponding to where [tˆ(s)] denotes integration over all possible paths the Lagrangian . in tangDent vector space from the tangent at one end tˆi Once P˜(f) =LZ˜(f)/Z˜(0) is calculated, performing an to the tangent at the other end tˆf. In d-dimensions ~r = inverse Laplace transform one can obtain the projected (r1,r2,...,rd). Recently a pathintegralGreens function probability distribution Px(x). Eq.4 can be written as, formulation has been developed [15] to evaluate the end to end distribution for a free polymer in 3D. We closely 1 P˜(f)= dv exp(tf v )p (v ) (8) follow that method and generalize it to obtain results x x x x Z−1 for various orientation constraints on polymer ends. In particularwefocusonpolymerslivingina2Dembedding where vx = x/L and px(vx) = LPx(x) is a scaling re- space. lation. Note that the Helmholtz free energy is given by The integrated (projected) probability distribution is x(vx) = (1/t)lnpx(vx). Thus Eq.8 gives the relation F − given by, betweentheHelmholtzandtheGibbsensembleforfinite chain (finite t), P (x)= d~rP(~r)δ(r x). (3) x 1 1 − Z exp[ tG(f)]= dvxexp(tf vx)exp[ t x(vx)]. We define the generatingfunctionofP (x)viaaLaplace − Z−1 − F x transform, Inthermodynamiclimitoft ,asteepestdescentap- →∞ proximation of the above integral relation gives G(f) = L (v ) fv , the well known Legendre transform rela- P˜(f)= dxexp(fx/λ)Px(x) (4) Ftioxn.xId−entifyxing iu= tf one can define Fourier trans- Z−L − 1 form relations, p˜ (u)= p (v )exp( iuv )dv and where f is the force in units of k T/λ i.e. f =Fλ/k T x −1 x x − x x B B appliedalongthex-axis. Again,thepartitionfunctionin 1 R∞ the Gibbs ensemble, Z˜(f~) = Ld~rexp(f~.~r/λ)Z(~r)[15]. px(vx)= 2π −∞dup˜x(u)exp(iuvx) (9) Z R 4 have a single maximum at v = 0 for all values of stiff- x 0.6 nessparametert. Nosuchsimplerelationexistsbetween p(v )andp (v )in2D.ThetwodimensionalWLCpoly- x x x 0.5 mer having its ends free to rotate may show more than one maximum in p (v ) and therefore non-monotonicity x x ) 0.4 x inforce-extension. Indeedourcalculationandsimulation v p(x0.3 (seeSec.III)doesshowmultiplemaximainprojecteddis- tribution p (v ) (Fig.1). This is a curious difference be- x x 0.2 tween semiflexible polymers in 2D and 3D. Because of this and the fact that experiments in 2D are possible[1], 0.1 in this work we focus on the 2D WLC polymers. We have already given a general form of Z˜(f) (Eq.7) 0 -1 -0.5 0 0.5 1 which depends on the dimensionality d of the embed- v x ding space. For d = 2, one can assume tˆ= (cosθ,sinθ), leading to = 1/4 θ˙2 + fcosθ . This automatically maintains tLhe ine{xtensibility cons}traint tˆ2 = 1. The an- FIG. 1: (Color online) For a semiflexible polymerin 2D hav- gular momentum p = ∂L = θ˙/2 and thus the corre- ing its ends free to rotate px(vx) ( = py(vy) ) is plotted θ ∂θ˙ at stiffness parameter t = 2. The points are collected from sponding Hamiltonian H = θ˙pθ −L = p2θ −fcosθ. In Monte-Carlo simulation in freely rotating chain model (see planar polar coordinates, replacing p i ∂ one ob- Sec.III). The line is calculated from theory (see Sec.II). The tains the corresponding quantum Hamθil→ton−ian∂θoperator, theory shows excellent agreement with simulation. It clearly Hˆ = ∂2 fcosθ. In this representation of tangent shows bimodality via two maxima in integrated probability −∂θ2 − vectors, distribution at thetwo near complete extensions. Z˜(f) = θ exp( tHˆ)θ i f h | − | i such that P˜(f) = p˜x(u = ift) and the inverse Fourier = φ∗n(θi)φn′(θf)hn|exp(−tHˆ)|n′i , (12) transformcanbewrittenasaninverseLaplacetransform, n,n′ X px(vx)=t21πiZ−ii∞∞dfP˜(f)exp(−tfvx) . (10) tdwhihreeecrtteoitoφannl(aθHs)ain=mEihlntqo.|4θn,ii.aHˆnIf=HeˆxHˆtde0re+nnaoHˆlteIfos=rcae−rii∂s∂gθ2ai2dp−prlfoiectdoosraθl(.oHˆnTghux=s- The simplest way to obtain p (v ), numerically, is to 0 replace f = iu/t in the exprexssioxn for P˜(f) to obtain −∂∂θ22) in presence of a constant external field (HˆI = − fcosθ). The eigenvalues of Hˆ are E = n2 and the p˜x(u) and evaluate the inverse Fourier transform (Eq.9). − 0 n complete set of orthonormal eigenfunctions are given by Up to this pointeverythinghas beentreatedin d- em- φ (θ)=exp(inθ)/√2π wheren=0, 1, 2,..., . In bedding dimensions. Experiments on single polymer can n be performed in three dimensions (3D) as well as in two this basis nHˆI n′ = (f/2)(δn′,n+1±+δ±n′,n−1).±T∞here- h | | i − dimensions (2D). In 3D, polymers are left inside a solu- fore, nHˆ n′ =n2δn′,n (f/2)(δn′,n+1+δn′,n−1). Ifthe h | | i − tion whereas one can float the polymer on a liquid film external force were applied in y- direction Hˆ = fsinθ I to measure its properties in 2D [2]. However, polymers and nHˆ n′ = n2δn′,n (f/2i)(δn′,n+1 δ−n′,n−1). embedded in 2D are more interesting because of the fol- n exhp(| t|Hˆ)in′ canbecalc−ulatedbyexponen−tiatingthe lowing reason. In a free polymer whose end orientations mh a|trix −nHˆ n|′ .iThusonecanfindZ˜(f)andhenceP˜(f) are free to rotate, the system is spherically symmetric h | | i and p (v ). and thus the probability distribution of end to end vec- x x Note that the above formalism can be easily extended tor P(~r) = P(r) where r = ~r . For this system it was | | to find the end to end vector probability distribution shown that in the Helmholtz ensemble in 3D[15], p(v ,v ). A Laplace transform of P(~r) is P˜(f~) = x y 1 dpx Ld~rexp(~r.f~/λ)P(~r). In a similar manner as above one p(v )= (11) x −2πvx dvx can show that P˜(f~) = Z˜(f~)/Z˜(~0) with Z˜(f~) given by R where p(v = r/L) = LdP(r) is the probability distribu- Eq.12 with Hˆ = −∂∂θ22 −fxcosθ−fysinθ. Thus, using tion of end to end distance scaled by the contour length aninverseLaplacetransformonecanfindP(~r)andhence L. In presence of the spherical symmetry of a free WLC p(vx,vy). polymer, this distribution gives the Helmholtz free en- ergy (v) = (1/t)lnp(v)[16]. P(r) is related to the radialFdistribu−tion function S(r) via S(r)=C rd−1P(r) A. Free polymer d where C is the area of a d- dimensional unit sphere. d Since p(v) is a probability distribution, p(v ) 0 and For a polymer which has both its ends free to ro- x ≥ therefore dp /dv 0 for v > 0 thus ruling out mul- tate, integrating Eq.12 over all possible initial and fi- x x x tiple peaks in p (v≤) [13] and showing that p (v ) will nal tangent vectors in rigid rotor basis one gets, Z˜(f)= x x x x 5 2π 0 exp( tHˆ)0 , Z˜(0)=2π and hence III. SIMULATION h | − | i P˜(f)= 0 exp( tHˆ)0 . (13) In this section, we introduce two discretized models h | − | i that we use to simulate semiflexible polymers. Both of ThismeansthatP˜(f)isgivenbythe(0,0)-thelementof these are derived from the WLC model which has been the matrix n exp( tHˆ)n′ . Thus, if the external force used for our theoretical treatment in Sec.II. After in- h | − | i f is appliedinx- direction, rememberingp˜ (u)=P˜(f = troducingthe discretizedmodels weshowhow toimpose x iu/t) one can calculate the inverse Fourier transform thevariousboundaryconditionsonendorientations. We −(Eq.9) to obtain p (v ). In this case, due to spherical perform Monte-Carlo (MC) simulations of these models x x symmetry of a polymer whose ends are free to rotate, to obtain probability distributions in the Helmholtz en- p (v )=p (v ). semble. x x y y One discretized version of the Fokker-Plank equation corresponding to the WLC model is the freely rotating chain (FRC) model[26, 27]. In the FRC model, one con- B. One end grafted siders a polymer as a random walk of N steps each of length b = L/N with one step memory, such that, suc- This symmetrybreaksdownimmediately ifone endof cessive steps are constrained to be at an fixed angle θ the polymer is fixed to a specific direction,namely along with λ=2b/θ2. The continuum WLC model is obtained the x-axis i.e. θ = 0. Then in Eq.(12) integrating over i inthe limitθ, b 0,N keepingλandLfinite. To all possible θ and leaving θ = 0 one obtains Z˜(f) = → →∞ f i simulate a polymerwith ends free to rotate a largenum- n exp( tHˆ)0 in the rigid- rotor basis. Note, for ber of configurations are generated with first step taken nh | − | i this case Z˜(0)=1 and therefore inanyrandomdirection. Whereasifonechoosesthefirst P step to be in some specific direction, this will simulate a P˜(f)= n exp( tHˆ)0 . (14) polymer with one end grafted in that direction. h | − | i A straight forward discretization of the Hamiltonian n X in Eq.1 in 3D (2D) is an 1d Heisenberg (classical XY) model: C. Both ends grafted β = κ N (tˆi−tˆi−1)2 = N ( J tˆi.tˆi−1) (16) Two ends of a polymer can be grafted in infinitely H 2 b − i=1 i=1 X X different ways. Let us fix the orientation of one end along x- direction (θ = 0) and the other end withanearestneighborcouplingJ =κ/bbetween‘spins’ i along any direction θf. Then Eq.(12) gives 2πZ˜(f) = tˆi. We have ignored a constant term in energy. The ap- n,n′ein′θfhn|exp(−tHˆ)|n′i, 2πZ˜(0) = neinθf−tn2 wpriothprJiabt=ecκonfitninitueu.mInlitmhiistmisordeceolvgerraefdtinfogribs→sim0u,laJte→db∞y and hence P P fixing end spins on the 1D chain. If an end is free then P˜(f)= n,n′ein′θfhn|e−tHˆ|n′i . (15) tbhyetehnedenseprignytaaknedseunptroapnyy.oIrnietnhtiastimonodtehla,tbyarfiexainllgowthede P neinθf−tn2 two end-spins, one can easily simulate a polymer with both its ends grafted in some fixed orientations. We fol- P low the normalMetropolis algorithm[46] to perform MC If the external force is in x- direction, the Laplace simulation in this model. transform of Z˜(f), defined in the way described above, We restrict ourselves to two dimensions. In the FRC gives the projected probability distribution in x- direc- model simulations we have used a chain length of N = tion, px(vx). On the other hand, if the external force is 103 and generated around 108 configurations. This sim- in y- direction, the Laplace transform of Z˜(f) gives the ulation does not require equilibration run. Therefore all projected probability distribution in y- direction py(vy), the 108 configurations were used for data collection. In the distribution of transverse fluctuation while one end theXYmodelwehavesimulatedN =50spinsandequi- of the polymer is grafted in x- direction. libratedover106 MCsteps. Afurther106 configurations All the relations derived so far are exact. Since the were generated to collect data. We have averaged over calculation of an infinite order matrix n exp( tHˆ)n′ 103 initial configurations, each of which were randomly h | − | i isnotfeasible,wecalculateitnumerically[45]bytruncat- chosen from nearly minimum energy configurations that inguptoanorderN ,thatcontrolstheaccuracy,limited conform with the boundary conditions. Increasing N, in d only by computational power. Unless otherwise stated, both the models of simulation, do not change the aver- weuseN =11whichalreadygivesverygoodagreement aged data. As a check on the numerics, we compared d with simulated data (see Fig.1 and Sec.IV). The inverse simulation evaluation of r2 and r4 with their exact h i h i Laplace transforms to obtain end to end probability dis- results[16, 38] to obtain agreement within around 0.5%. tributions from P˜(f)s are also done numerically. NoticethatinsimulatingtheFRCmodeloneperforms 6 40 probability distributions obtained at stiffness parameter t=1 t = 2 with that obtained from MC simulation (Sec.III) t=3.33 using the FRC model (see Fig.1). This shows excellent 30 t=10 agreement between theory and simulation. For a free 2 1/[8(1-<v>) ] polymer p (v ) and p (v ) are same due to the spheri- x x y y <v> calsymmetry. Notethat (v )= (1/t)lnp (v )would 20 x x x F − f givea non- monotonicforce-extension f -v due to the x x h i multimodality in p (v ) (Fig.1) via f = (∂ /∂v ). x x x x h i F 10 The force-extension obtained from the projected proba- bilitydistributionp (v )correspondstotheexperimental x x scenario in which the external potential traps the poly- 0 mer end only in the x- direction and the polymer-end is free in y. In general, if the external potential traps 0 0.5 1 the polymer-end in d dimensions (d d) then a d di- <v> r r ≤ r mensional projection( [d d ] dimensional integration) r − of the probability distribution of end to end vector p(~v) gives the appropriate free energy and decides the force- FIG.2: (Coloronline)Noneoftheforce-extensioncurves,ob- tainedintheGibbsensemble,includingthatatt=3.33show extension relation. On the other hand, if the trapping non- monotonic behavior unlike in the Helmholtz ensemble potentialholds apolymer-endinallthed-dimensions,as [16]. ForcesareexpressedinunitsofkBT/λi.e. f =Fλ/kBT. isusuallydoneinmostforce-extensionexperiments,only the end to end vector distribution p(~v) gives the appro- priate Helmholtz free energy that can predict the force- random walk with fixed angle between consecutive steps extension behavior in the Helmholtz ensemble. This un- and does not require to equilibriate. Thus one uses all derstanding is general and does not depend on the spe- the simulated configurations for data collection. On the cificorientationalboundaryconditionsorthedimension- other hand, in simulating XY-model one has to perform ality d of embedding space. This is important to keep in equilibration runs over a large number of steps and an mindwhile analyzingexperimentaldata. Inexperiments averaging over many initial configurations is required. that use the laser tweezers to trap polymer ends in d- Anotherimportantdifferencebetweenthetwosimulation dimensions, ends remain free to rotate and the relevant methodsisthat,intheXY-modelsimulation,ineachMC Helmholtz free energy is obtained from p(v). Ref.[16] step one has to calculate a time consuming exponential predicted multiple minima in this free energy leading to of change in energy, whereas, no such exponential calcu- non-monotonic force-extension in such experiments. lation is required in simulating the FRC model. Thus The Gibbs ensemble: We have already mentioned that simulating the FRC model is clearly much faster, com- the non- monotonic nature of force-extension, a strong putationally. However,implementingthefixedboundary qualitativesignature ofsemiflexibility, is observable only orientationsatboththeendsofapolymerismucheasier in the Helmholtz ensemble and not in the Gibbs en- in the XY-model. semble [16]. In the Gibbs ensemble, the averaged ex- tension comes out to be v = (∂G/∂f) and the re- sponse ∂ v /∂f = t[ v2 h iv 2]− 0. Similar relation IV. RESULTS h i h i−h i ≥ for response function does not exist in the Helmholtz ensemble. Therefore, the force-extension in the Gibbs Once all these theoretical and simulation tools are ensemble has to be monotonic (Fig.2) in contrast to the available, we apply them to bring out the statistical and Helmholtz ensemble. For a polymer with its ends free mechanicalpropertiesofasemiflexiblepolymer. Wehave to rotate, the force extension relations, that have been three different boundary conditions depending on the calculated from theory, at various t are shown in Fig.2. orientational constraints on the polymer ends and two For small forces the polymer shows linear response. At different ensembles. For each case we look at the vari- large and positive force polymer goes to fully extended ous probability densities, ensemble dependence of force- limit beyond which, the inextensibility constraint pre- extension etc. For the case of a polymer with both ends vents further extention. It is possible to do perturbative graftedwefindthatthepropertiesdependontherelative analysis of P˜(f) = 0 exp( tHˆ)0 in the two extreme orientation of the two ends. limits of small andhla|rge f−orces|toi obtain the asymp- toticforce-extensions[47]. Inthesmallforcelimit,fcosθ may be treated as a perturbation about the rigid ro- A. Free polymer tor hamiltonian Hˆ = ∂2/∂θ2. Thus keeping upto 0 − the second order correction to eigen-values we obtain The Helmholtz ensemble: We employthetheoryasde- E0 = f2/2. Within this perturbative approximation − scribedinSec.IItocalculatep (v )andp (v )forapoly- P˜(f)=exp( tE )andthereforeG(f)= 1/t lnP˜(f)= x x y y 0 mer with both its ends free to rotate. We compare the f2/2. Thu−s the force extension relatio−n in this limit − 7 2 0.8 FRC: p (v ) x x FRC: p (v ) y y 1.5 XY: p (v ) x x 0.6 XY: p (v ) y y Theory : p (v ) x x pα 1 Theory : py(vy) v)y0.4 ( 0.65 y p t=0.5 t=0.75 0.5 t=2 0.625 0.2 t=2.8 LMF t=4 0 0.6 -1 -0.5 0 0.5 1 -1 -0.5 0 0.5 1 v α 0 -1 -0.5 0 0.5 1 v y FIG. 3: (Color online) The simulation data for px(vx) and py(vy)fromtheFRCmodelandtheXYmodelarecompared withtheirtheoreticalestimates. Simulationsandcalculations FIG.4: (Coloronline)Forapolymerwithoneendgrafted in were done at t= 2 for a polymer with one end grafted in x- x- direction, the integrated probability distribution py(vy) is direction. plottedatvariousstiffnessest. Att=4thereisasinglemax- imum at vy =0. Decreasing t we see at t=2.8 emergence of two more peaks at nonzero vy apart from the one at vy = 0 (See inset). At t=2 the central peak vanishes, the trimodal is v = ∂G/∂f = f. On the other hand, for large forhceis one−can expand the term cosθ 1 θ2/2 and distribution becomes bimodal. The circles labeled LMF are ≃ − data taken from Ref.[17] at t = 2 and show excellent agree- write Hˆ = tf +tHˆ where Hˆ = ∂2/∂θ2+(1/2)fθ2 0 0 mentwithourtheory. Att=0.75weseere-emergenceofthe − − is the harmonic oscillator Hamiltonian. In the harmonic centralpeakandtrimodalityinpy(vy)(Seeinset,3sarefrom oscillator basis, the ground state eigenvalue E0 = f/2 ourMCsimulation intheFRCmodel at t=0.75.) Thelines and thus the ground state energy corresponding to Hˆ is are calculated from theory. p f+ f/2. Therefore,in a similar manner as in above, − P˜(f)=exp( t f/2+tf)andG(f)= f/2 f. Thus, 5 fboerlvarpg=efor∂cG−es/,p∂tfhe=fo1rce1-e/x(t2e√n2sifo)nwrehliacthpiocnancob−meeisnvoeurttetdo 43 FF((vvxy)) 24 <<ffx>>((vvx)) tfo(ghveit) itnh−eFrige.l2atfiaolnls, ofn−=to1f/[8=(1v−havti)f2]. A0lllimthiteacnudrvteos F2 ><fα 0 y y 1/h[8(i1 v )2] in the f lihmiit. → 1 -2 −h i →∞ 0 -1 -0.5 0 0.5 1 -4-1 -0.5 0 0.5 1 v v α α B. Grafted polymer: One end FIG. 5: (Color online) The left panel shows the Helmholtz The Helmholtz ensemble: Letuscompareourtheoreti- free energies F(vx) and F(vy) of a polymer at t=2 and one calandsimulationestimateofpx(vx)andpy(vy)att=2 end grafted in x- direction. The right panel shows the corre- (Fig.3) for a semiflexible polymer with one end grafted sponding force-extensions in the Helmholtz ensemble. Both in x- direction. The excellent agreement validates both hfxi- vx and hfyi- vy show non- monotonicity and regions of our theory andthe simulationtechniques. In px(vx), the negative slope. Free energies are expressed in units of kBT peak in near complete extension along positive x is due and forces are expressed in unitsof kBT/λ. to the coupling of the end orientation towards this di- rection with large bending energy (also see Fig.12). We then explore, in detail, the transverse fluctuation p (v ) behavior. Theothertwopeaksemergeasentropytriesto y y of this system for different t (Fig.4). At large t(= 10), foldthepolymerandenergyrestrictstheamountofbend- p (v ) has single maximum at v =0. At such low stiff- ing. Since bending in positive and negativey- directions y y y nesses entropy takes over energy contributions. Number are equally likely, the transverse fluctuation shows two of possible configurations and thus entropy gains if end new maxima near v = 0.5 symmetrically positioned y ± toendseparationremainsclosetozero. Thisgivesriseto aroundv =0. Withfurther increaseinstiffness(t=2), y the single central maximum. The emergence of multiple the central entropic peak vanishes (also see Fig.12) and maxima atnonzero v , the multimodality, atlargerstiff- p (v ) becomes bimodal with two maxima (Fig.4). At y y y ness (t=2.8) is due to the entropy- energy competition. even higher stiffness (t = 0.75) the central peak reap- The central peak is due to the entropy driven Gaussian pears, due to a higher bending energy. At t = 0.5 the 8 0.8 stant and the corresponding average force in x- [y-] di- t=0.5 rection is found from the relation fx = ∂ (vx)/∂vx h i F t=1 ( or f = ∂ (v )/∂v ). Notice that, when v [v ] is y y y x y 0.6 t=2 h i F heldconstant,v [v ]remainsfree. This canbe achieved t=5 y x using a trapping potential constant in v [v ] and trap- y x > ping the polymer end in vx [vy]. In Fig.5, we show the x0.4 v Helmholtz free energies (v ) = (1/t)ln p (v ) and < F x − x x (v ) = (1/t)ln p (v ) and the corresponding force y y y F − extension curves in constant extension ensemble. Note 0.2 that unlike the monotonicity obtained in v -f curve y y h i (Fig.6) in the Gibbs ensemble, the f -v curve in Fig.5 y y h i clearly shows non-monotonicity,a signature of semiflexi- 0 -30 -20 -10 0 10 20 30 bility in the Helmholtz ensemble. f The Gibbs ensemble: From our theory we can also ex- y plorethetransverseresponseofapolymerwhichhasone of its ends grafted and a constant force is applied to the other end in a direction transverse to the grafting direc- 1 tion. Assume that the grafting direction is x and a force t=0.5 f is applied in y- direction to study the transverse re- y t=1 sponse. Alinearresponsetheorywasproposedearlier[36] 0.5 t=5 totackle this question. Ourtheorycanpredictthe effect t=0.5: LMF of externally applied force f of arbitrary magnitude on t=1: LMF y <v>y 0 t=5: LMF tphlieedavinerya-gdeirpeocstiitoinoni.se.hvf~xi=ayˆnfdy,hvwyei.haAvseHthIe=forcfeyissinapθ-. − Because one end of the polymer is grafted in x- direc- -0.5 tion we use nHˆI n′ = (fy/2i)(δn′,n+1 δn′,n−1) to h | | i − − evaluate Z˜(f ), whereas to calculate v = (∂G/∂f ) y x x h i − [ or, v = (∂G/∂f ) ], we introduce a small perturb- y y h i − --130 -20 -10 0 10 20 30 ing force δfx [ or, δfy ] in the Hamiltonian matrix to obtainthe partialderivatives. Thus we obtainthe corre- f y spondingforce-extensionsshowninFig.6. Asthegrafted end is oriented in x- direction, we expect, in absence of any external force, v will be maximum and will keep x FIG. 6: (Color online) Average displacements along x- di- h i on reducing due to the bending of the other end gen- rection hvxi and y- direction hvyi as a function of transverse erated by the external force f imposed in y- direction. force(transversetograftingdirectionx)inconstantforceen- y Thus v is expected to be independent of the sign of semble. Lines denote our theoretical calculation while points h xi f . Similarly, v should follow the direction of exter- denote the MC simulation data taken from Ref.[17]. Forces y y h i (Fy) are expressed in units of kBT/λ,i.e. fy =Fyλ/kBT. nal force and therefore is expected to carry the same sign as f . Fig.6 verifies these expectations and shows y very good agreement between our theory and simulated distribution again becomes single peaked at v = 0 as datatakenfromRef.[17]. Itisinterestingtonotethat,in y bending energy takes over entropy and the polymer be- theHelmholtzensemble,themultimodalityinprobability comes more like a rigid rod. However, even at very high distributionpredictsnon-monotonicityinforce-extension stiffnessliket=0.5thesinglepeakeddistributionp (v ) relation. However, as expected, this non-monotonicity y y isquitebroadunderliningtheinfluenceofentropicfluctu- does not survive in the Gibbs ensemble. ations. Notice thatwe haveplottedMC datatakenfrom Ref.[17] for the XY model simulation at t = 2 (Fig.4). This shows very good agreement with our theory. Infact C. Grafted polymer: Both ends all the simulated data from Ref.[17] at different t show excellent agreement with our theoretical predictions. In The Helmholtz ensemble: Let us first fix the orienta- the inset of Fig.4, we have magnified the multimodality tions of the polymer at both its ends along x- axis and at t = 2.8 and t = 0.75. We have also plotted our FRC comparep (v ) andp (v )obtainedfromourXY model x x y y model simulation data at t = 0.75 and obtained very simulationand our theory (Fig.7). The very goodagree- good agreement. ment validates both our theory and simulation. Then, At this point, it is instructive to look at the force we go on to explore the properties of this system using extension behavior in the Helmholtz ensemble, the en- the theory developed in Sec.II-C. Let us fix the orienta- semble in which p (v ) and p (v ) have been calcu- tion at one end in x- direction (θ = 0) and that in the y y x x i lated above. In it the extension v [v ] is held con- otherend(θ )canbevariedtostudythechangeintrans- x y f 9 4 t=1 p (v ) x x t=2 p (v ) t=4 3 y y 1 t=10 Theory: p (v ) x x Theory: p (v ) ) y y vy pα2 p(y 0.5 1 0 -1 -0.5 0 0.5 1 0 v -1 -0.5 0 0.5 1 y v α 1 FIG. 7: (Color online) The simulation data for px(vx) and θ = 0 py(vy)fromtheXYmodelsimulationsofaWLCpolymerare θ = π/4 θ = π/2 compared with their theoretical estimates. Simulations and 0.8 θ = 3π/4 calculationsweredoneatt=2foraWLCpolymerwithboth θ = 7π/8 θ = π its endsgrafted in x-direction. 0.6 ) y v ( 0.88 y p 0.4 versefluctuationp (v ). Tobeginwith,letusfindp (v ) y y y y for different stiffness parameters t with θ = 0 (Fig.8). 0.84 f 0.2 The height of the central peak shows non-monotonicity – with increase in t from t = 1 the height of the central 0.8 0 0.4 peak first decreases up to t = 2 and then eventually it 0 -1 -0.5 0 0.5 1 increases again. The initial decrease in peak height is v y due to the fact that with increase in t, i.e. with lower- ing in stiffness, the other end of the polymer (relative to the first end) starts to sweep larger distances from the FIG. 8: (Color online) The upper panel shows py(vy) for a x- axis. With further increase in t (t =4), the height of polymerwith both ends grafted along x-direction at various the maximum increases (also see Fig.12). From Fig.12, stiffness parameters t. They always show single maximum. notice that at t=4 multimodality appears in the distri- In lower panel, py(vy) is plotted for various relative angles θ butionofendtoendvector. Thenewentropicmaximum between the orientations of the two ends at t=4. The inset at ~v =~0 contributes towards increasing the peak height magnifies the emergence of bimodality at θ=π/4. in p (v ) at v = 0. Though, in p(v ,v ) multimodality y y y x y is present (see Fig.12) at t = 4, after integration over probability weights along x- direction the projected dis- inate to bring down the peak positions to lesser v with y tributionpy(vy)becomesunimodal. Thusmultimodality respect to that attained at θ = π/2. However, entropy in the probability distribution of end to end vector does alwaysplayacrucialrole,e.g. atθ =π/4,p (v )showsa y y notguaranteemultimodalityinprojectedprobabilitydis- double peak around v =0.37. At θ =π the two ends of y tributions. thepolymerarekeptanti-parallel. Noticethat,asθ =π To see the impact of change in relative angle of graft- and θ = π are physically same, at θ =π, energetically, − ing, we fix one end along x- axis and rotate the orienta- vy = 0.64 are equally likely. Entropy would like the ± tion of the other end and find out the transverse fluctu- two ends to bend to vy = 0. Competition between en- ation p (v ) at t = 4 (Fig.8). At θ θ = 0 the fluctu- ergy and entropy leads to almost a constantdistribution y y f ation is unimodal with the maximu≡m at vy = 0. With up to |vy|∼0.5. The behavior of py(vy) for −π ≤θ ≤0 increaseinθtheorientationoftheotherendrotatesfrom is mirror symmetric about vy = 0 with respect to the positive x- axis towards positive y- axis. Energetically behavior of py(vy) in the region 0 θ π. ≤ ≤ thepolymergainsthemost,ifitbendsalongtheperime- The Gibbs ensemble: We then work in the constant terofacircle. Therefore,energetically,atanyθ,thepeak force ensemble by applying a force f~ = yˆf on an end y ofp (v )wouldliketobeatv =(1 cosθ )/θ. Thusat oriented along any direction θ to x-axis while the other y y y − θ = 0,π/4,π/2,3π/4,7π/8,π the peak of p (v ) should end is oriented along x- direction. We find out the cor- y y beatv =0,0.37,0.64,0.72,0.69,0.64respectively. Fig.8 responding responses, v -f and v -f to this force y x y y y h i h i shows that the peak positions almost follow these values (Fig.9) in the similar manner as has been done in the up to θ =π/2, above which entropic contributions dom- last subsection for the case of a polymer with one end 10 1 0.8 t=1 t=2 0.8 t=3 0.6 t=4 t=5 > 0.6 t=10 vx0.4 θ=0 >x sin θ/θ < θ=π/4 <v θ=π/2 0.4 0.2 θ=3π/4 0.2 0 -10 -5 0 5 10 0 -4 -2 0 2 4 f θ y 0.8 0.8 0.6 0.6 0.4 0.4 0.2 0.2 > > vy 0 vy 0 t=1 <-0.2 θ=0 < -0.2 t=2 θ=π/4 t=3 -0.4 θθ==π3π/2/4 -0.4 tt==45 -0.6 θ=π -0.6 t(=1-1c0os θ)/θ -0.8 -0.8 -10 -5 0 5 10 -3 -2 -1 0 1 2 3 θ f y FIG. 10: (Color online) The upperpanelshows thevariation FIG. 9: (Color online) Average displacements hvxi (upper ofhvxiasafunctionofθandthelowerpanelshowsthevaria- panel) and hvyi (lower panel) as a function of a force fy for tionofhvyiasafunctionofθ. hvxiandhvyiarecalculatedfor a polymer having one end grafted along x- direction and the stiffness parameters t = 1,2,3,4,5,10. The thick solid line, otherinanangleθtothex-direction. Forcesareexpressedin in both the plots, show the expected behavior coming from units of kBT/λ. All the force-extension curves are obtained energetics ignoring the entropy. at t=1. elastic constant ∂f /∂ v near f = 0 (linear response) y y y h i grafted. If θ = 0, the force extensions carry the same islargeratθ =0ascomparedtoatθ =π; i.e. thetrans- qualitative features as for a single end grafted polymer verseresponseofasemiflexiblepolymerwithparallelend at all t (see θ = 0 curves for t = 1 in Fig.9). Therefore, orientationsismorerigidthanwithanti-parallelendori- instead of showing the t dependence of force-extension entations. To see the impact of the change in relative behavior, we show the θ dependence of force extensions angleθ,indetail,wecalculate vx and vy aswevaryθ h i h i at t = 1. The peak in v -f curve shifts to f < 0 as (Fig.10) keeping external force at zero. Bending energy x y y h i θ is increasedup to π/2 above which it againshifts back would like vx =sinθ/θ and vy = (1 cosθ)/θ. Note h i h i − towards fy = 0. With increase in θ, hvxi decreases, as that at θ → 0, energetically, hvxi → 1 and hvyi → θ/2. with these boundary orientations the polymer is forced Again, at θ π bending energy requires vx 0 and →± h i→ to close in x- and open up in y- direction. However, for v 2/π though entropy likes ~v ~0. Thus at y h i → ± h i → θ π entropy likes v 0. For θ < π/2 small small t, the approach of v -θ curve to sinθ/θ is much y x → ± h i → h i negative f leads to unfolding thereby increasing v . better than approach of v -θ to (1 cosθ)/θ (Fig.10). y x y h i h i − Whereas for θ >π/2 the effect of negative force is oppo- It should be noted that the angle θ in this study de- site – it helps the polymer to get folded to reduce v . notes a relative angle of bending between the two end x h i At θ = π, v always remains zero. The responses for orientationsofaWLC polymer. This shouldnotbe con- x h i negative θ are reflection symmetric about f = 0. The fused with the twist angle as in Ref.[41]. In an earlier y foldingbehaviorisalsoapparentfrom v -f curves. Up study[38] the impact of changing θ on the averagedroot y y h i toθ =π/2theresponseshiftstowardspositive v asthe mean squared end to end vector has been obtained. In y h i polymerlikesto openupiny-directionduetothe bend- this section we have shown the impact of changing θ on ing energy cost. However, for large θ entropy wins and projected probability distribution, averaged end to end at θ =π, v -f curve, again, goes through origin. The distance ( v , v ) and force-extension relations. y y x y h i h i h i