Globular Structures of a Helix-Coil Copolymer: Self-Consistent Treatment. C. Nowak, V.G. Rostiashvili, and T.A. Vilgis Max-Planck-Institut fu¨r Polymerforschung, Ackermannweg 10 , 55128 Mainz, Germany Aself-consistentfieldtheorywasdevelopedinthegrand-canonicalensembleformulationtostudy transitions in a helix-coil multiblock globule. Helical and coil parts are treated as stiff rods and 7 self-avoidingwalksofvariablelengthscorrespondingly. Theresultingfield-theorytakes,inaddition 0 totheconventionalZimm-Braggparameters,alsothree-dimensionalinteractiontermsintoaccount. 0 Theappropriatedifferentialequationswhichdeterminetheself-consistent fieldsweresolvednumer- 2 ically with finite element method. Three different phase states are found: open chain, amorphous n globule and nematic liquid-crystalline (LC) globule. The LC-globule formation is driven by the a interplay between the hydrophobic helical segments attraction and the anisotropic globule surface J energyofanentropicnature. Thefullphasediagramofthehelix-coilcopolymerwascalculatedand 9 thoroughlydiscussed. Thesuggested theoryshowsaclearinterplaybetweensecondaryandtertiary 2 structuresin globular homopolypeptides. ] PACS numbers: 61.30.Vx Polymer liquid crystals; 87.14.Ee Proteins; 87.15.-v Biopolymers: structure and t f physicalproperties o s . t I. INTRODUCTION a m - Significant work has been done to investigate the helix-coil transition theoretically1 and computationally2,3. One d of the first and well-known approaches is the Zimm-Bragg theory4 which is designed for single-stranded polypeptide n chains. Itconsidersaone-dimensionalIsingtype modelinwhicheachsegmentcanbe inoneofthe twostates: helical o state (stiff) or coil state (flexible). The helix is stabilized by hydrogen bonds which generates an energy gain of ǫ c − [ for each segment in the helical state. This energy gain is partly balanced by an entropy loss ∆S. The free energy difference between the helical state and the coil state is therefore given by ∆f = ǫ+T∆S fo−r each segment. 2 In an α-helix hydrogen bonds can be only formed between the ith and the (i+−3)th peptide group. The formation v of a hydrogen bond between the first and the third peptide group fixes the conformation of three groups. The next 1 bond betweenthe secondand the fourthgroupfurnishes the same energygainbut fixes only one new groupandthus 6 6 leads to a much smaller entropy loss. The formation of an α-helix ( as a simple element of the secondary structure) 2 is therefore a cooperative process and the formation of a junction between helix and coil is energetically unfavorable. 1 This can be modelled by an energy penalty µ for each junction between two successive groups of helical and flexible J 6 segments. Similarargumentsholdforallkindofhelices,suchthatthehelixformationisalwaysacooperativeprocess. 0 It is convenient to define the following fugacities / t a s e−β∆f, σ e−2βµJ, (1.1) m ≡ ≡ - where s gives the statistical weight of a helical segment compared to a coil segment. The cooperativity parameter d σ gives the statistical weight of a junction point. σ = 1 corresponds to µ = 0 and therefore to a non-cooperative n J o system. σ 0 corresponds to µJ and therefore to a totally cooperative system, i.e. either the entire chain c forms one b→ig helix or no helix is for→me∞d at all. In most helix forming biopolymers σ is roughly 10−3 10−4 (see e.g. : ref.1). − v i Using a transfer matrix method1 the one dimensional model can be solved exactly. The s - dependance of the X fractionofstiffhelicalsegmentΘ =N /N hasatypicalsigmoid-type form. With increasingcooperativity(i.e. as R R r the parameter σ is reduced) the transition becomes sharper and sharper. This rather sharp crossovertransition (due a to the cooperativity effect) is also observed experimentally, for instance in polybenzylglutamate4. Severalextensionshavebeenmadetothisone-dimensionalIsingtypemodel. Ithasbeenshown5 thatthetransition becomeslesscooperative,ifthehydrogen-bondingabilityofthesolventistakenintoaccount. Thehelix-coiltransition in grafted chains was studied6 as well as the effect of an external applied force on the transition7. Ofspecialinterestistheapplicationofone-dimensionalmodelstoproteins,seeforinstance8,9. Astudyofthehelix- coil transition including long-range electrostatic interactions10 can, to some extent, explain the amount and location of helical segments in globular proteins. However, to understand how α-helices (or generally secondary structure elements) are formed in the folding process of a protein and how this influences the compaction and formation of tertiary structure (and vice versa), it is necessary to combine the one-dimensional physics of the helix-coil transition withtheinteractionsofsegmentswhichapproacheachotherinthethree-dimensionalspace. Thisenablesadescription Typeset by REVTEX 2 of the interplay between the secondary structure and the mesoscopic structure formation of the entire chain (known as tertiary structure) . In proteins the stiff helical parts are often hydrophobic since the hydrogen bonds stabilizing the helix composition (as opposed to coil parts) are mainly saturated (due to bond formation between the ith and the (i+3)th peptide groups). This hydrophobicity causes an additional attractive interaction between stiff parts and drives the protein into a compact globular phase. Statistical analysis of the data from 41 globular proteins in native and partially folded conformational states11 showed a strong correlation between the amount of secondary structure elements and compactness of the proteins. This indicates that the formation of secondary structure (for instance α-helices) and the hydrophobic collapse into a compact globule occur simultaneously. This problem has been partially discussed withincomputersimulationsofglobularproteins12,13,14,15. Amongotherthings,ithasbeenshown13,14,15 (asopposed to the earlier findings by Dill et al.12) that the compactness itself, driven by the hydrophobicity, is insufficient to generateanyappreciablesecondarystructure. Itwasnecessarytointroducealocalconformationalpropensitytoward α - helix formation. The interplay of compaction and secondary structure leads to the formation of the specific three-dimensional tertiary structure. Computational and experimental studies of this mechanism can for instance be found in ref.16,17. Tostudytheinterplayofhelix-coil(orstiff-flexible)transitionandcollapsetransitionofthepolymerintoacompact globule,wehavedevelopedanapproachwhichcombinesvariablecompositionwiththree-dimensionalexcludedvolume interactions using self-consistent field theory. The paper is organized as follows. In Section II we have covered the generalself-consistentfieldtheoryofahelix-coilcopolymericglobule. Thefinalequationsarewrittendowninaform which is convenient for the numerical solution. Section III is devoted to the analysis of this solution. The formation of different globule structures (e.g. amorphous and liquid - crystalline (LC) globules) are studied in details. Among other things we argue that the LC - globule formation is mainly driven by the globule surface tension anisotropy. Finally, in Section IV we sketch the main results and compare them with the appropriate findings of some other authors. II. THEORY In this section the derivation of the field theory is outlined. For a more detailed description see18. Regarding all three dimensionalinteractionsandentropiccontributionsthe helicesaremodelled asstiff rods. Hence forone specific microscopiccompositionofhelicalandflexiblepartsthesystemlookslikearod-coilmultiblockcopolymerasshownin Fig.(1). The multiblock copolymermaybe composedofK rod-coilblocks. The conformationofthenthrod-coilblock is given by the vector-function r (s) describing the contour of the coil, by the vector R which gives the position of n n the junction point between rod and coil and the unit vector u describes the orientation of the rod, see Fig.(1). The n u r (f N ) n n n n R r (0) n n r (s) n u n - 1 R n - 1 FIG.1: Rod-coilmultiblockcopolymerparameterization. Twosuccessiverod-coilblockswithordinalnumbersn−1andnare shown. Theflexiblechaincoordinates usethevectorsetr(s). Theorientation ofthenth rodisdenotedasu andthejunction n point between rod and coil is given byR . n length of the nth rod-coilblock in units of the segment length b is given by N . The fraction of the flexible segments n (coil) is givenby f and the fraction of stiff segments (rod) thus by 1 f . The microscopic flexible segment density n n − 3 ρˆ (r) and stiff segment density ρˆ (r) can now be defined as follows. C R K fnNn ρˆ (r)= dsδ(r r (s)), (2.1) C n − n=0 Z X 0 K (1−fn)Nn ρˆ (r)= dsδ(r R u sb). (2.2) R n n − − n=0 Z X 0 In addition an orientation density Sˆij(r) is introduced which is sensitive to the collective orientation of the system K (1−fn)Nn 1 Sˆij(r)= dsδ(r R u sb) uiuj δij . (2.3) n n − − − 3 n=0 Z (cid:20) (cid:21) X 0 The interaction Hamiltonian of the model can thus be written as v βH [ρˆ ,ρˆ ] = χ d3rρˆ (r)ρˆ (r)+ d3r[ρˆ (r)+ρˆ (r)]2 int C R R R C R 2 Z Z w + d3r[ρˆ (r)+ρˆ (r)]3+g d3rTr Sˆij(r)Sˆij(r) , (2.4) C R 3! Z Z h i where v and w control the strength of the global two- and three-body interactions. The parameter χ controls a selective two-body interaction between the stiff segments which is caused by the hydrophobic nature of the helical parts of the chain (cf. Introduction). The last term is a alignment interaction between the rods of the standard Maier-Saupe form19. The representation of this type has a wide use in the polymeric liquid-crystal context20. The canonical partition function of the entire system can be written as the functional integral over the collective coil, rod and orientation densities, ρ (r), ρ (r) and Sij(r) respectively, i.e. C R K Z( N ,K) = ρ (r) ρ (r) Sij(r) Dr (s)d3R d2u n C R n n n { } D D D Z Z Z Z n=1 Y δ(ρ (r) ρˆ (r))δ(ρ (r) ρˆ (r))δ(Sij(r) Sˆij(r)) C C R R × − − − 3 fnNn ∂r 2 δ(u 1)δ(r (f N ) R )exp ds n × | n|− n n n − n −2b2 ∂s Z0 (cid:18) (cid:19) v exp χ d3rρ (r)ρ (r) d3r[ρ (r)+ρ (r)]2 R R C R × − − 2 (cid:26) Z Z w d3r[ρ (r)+ρ (r)]3 g d3rTr Sij(r)Sij(r) . (2.5) C R − 3! − Z Z (cid:27) (cid:2) (cid:3) Bymakinguseoftheintegralrepresentationfortheδ-functioninEq.(2.5)(whichresultsinanappearanceofexternal fields h (r), h (r) and hij(r)) the partition function of the entire system reads C R S Z( N ,K) = ρ (r) ρ (r) DSij(r) h (r) h (r) Dhij(r) { n} D C D R D C D R S Z Z Z Z Z Z v exp χ d3rρ (r)ρ (r) d3r[ρ (r)+ρ (r)]2 R R C R × − − 2 (cid:26) Z Z w d3r[ρ (r)+ρ (r)]3 g d3rTr Sij(r)Sij(r) C R − 3! − Z Z (cid:2) (cid:3) + i d3rρ (r)h (r)+i d3rρ (r)h (r)+i d3rSij(r)hij(r) C C R R S Z Z Z (cid:27) Z(0)( N ,K;[h ],[h ], hij ) , (2.6) × { n} C R S h i 4 where the partition function Z(0) of the non - interacting system in the external fields h , h , hij. is given by C R S Z(0)( N ,K;[h ],[h ], hij ) { n} C R S K h i = r (s)d3R d2u δ(u 1)δ(r (f N ) R ) n n n n n n n n D | |− − Z n=1 Y 3 fnNn ∂r 2 fnNn exp ds n i dsh (r (s)) × (−2b2 Z0 (cid:18) ∂s (cid:19) − Z0 C n (1−fn)Nn i dsh (R +u s) R n n − Z0 (1−fn)Nn 1 i dshij(R +u s) uiuj δij . (2.7) − Z0 S n n (cid:18) n n− 3 (cid:19)) The composition of the system is assumed to be equilibrated with respect to the total number of stiff (N ) and R flexible (N ) segments as well as to the number of junction points (N ) between stiff and flexible parts . Therefore C J the description canbe reduced, i.e. N ,K N ,N ,N , and it is more convenientto switch to a grandcanonical n R C J partition function Z of the interactin{g p}olym→er written in terms of the grandcanonical partitionfunction Z(0) of the non-interacting polymer in the external fields h , h , hij, i.e. C R S Z(µ,ǫ,σ) = ρ (r) ρ (r) DSij(r) h (r) h (r) Dhij(r) D C D R D C D R S Z Z Z Z Z Z exp βH [ρ ,ρ ]+i d3rρ (r)h (r)+i d3rρ (r)h (r)+i d3rSij(r)hij(r) × − int C R C C R R S (cid:26) Z Z Z (cid:27) Z(0)(µ,ǫ,σ;[h ],[h ], hij ). (2.8) × C R S h i In Eq.(2.8) µ denotes the chemical potential conjugated to the whole number of segments and ǫ is the energy gain − of a helical segments compared to coil segment. The meaning of the cooperativity parameter σ was discussed in the introduction. The grand canonical partition function Z(0) of the non-interacting system can be derived by using the polymeric correlation function Ξ(0), i.e. Z(0)(µ,ǫ,σ;[h ],[h ], hij ) = d1d1′ Ξ(0)(1,1′;µ,ǫ,σ;[h ],[h ], hij ). The poly- C R S C R S mericcorrelationfunctionΞ(0)(1,1′;µ,ǫ,σ;[h ],[h ],hhiji)givRestheconditionalunnormalizedprohbabiilityoffinding C R S ′ the first segment of the copolymer at the coordinate 1hproivided that the last segment is at 1. The symbol 1 stands either for r or for (r ,u ) depending on whether the first segment is a flexible or a stiff one. The same holds for 1 1 1 ′ the coordinate 1 of the last segment. The polymeric correlation function is therefore the partition function of the ′ multiblock helix-coil copolymer with two ends fixed at 1 and 1. By using the equations of motion for the rod and coil Green functions - G and G - the polymeric correlation coil rod function Ξ(0) can be constructed. The inverse Greens function operators are given by18 G−1 = δ(r r′) βµ b2 2+ih (r) (2.9) coil − − 6 ∇r C (cid:18) (cid:19) Gb−1 = δ(r r′) β(µ ǫ) b2(u·∇R)2 +ih (r)+i↔h (r):↔P , (2.10) rod − − − β(µ ǫ) R S (cid:18) − (cid:19) ↔ where the tensor P is givben by Pij uiuj δij/3. In Eqs.(2.9) - (2.10) we have set the elementary unit lengths in ≡ − coilandrodpartsequaltoeachother(namelytob)forsimplicity. Thisaffectssomequantitativevaluesbutdefinitely does not alter qualitative predictions of this paper. Knowing the inverse Green function operators the grand canonical polymeric correlation function can now be represented as the following Gaussian 2-dimensional path integral 1 Ξ(0)(1,1′;µ,ǫ,σ;[h ],[h ], hij )= ψ ϕψ(1)ϕ(1′) C R S Θ D D 1h i ψZ(3) T G−1 σ1/2 ψ(4) ×exp(−2Z d3d4(cid:18)ϕ(3) (cid:19) (cid:18)−σr1o/d2 −G−co1il (cid:19)(cid:18) ϕ(4) (cid:19)), b (2.11) b 5 where 1 ψ(3) T G−1 σ1/2 ψ(4) Θ=Z DψDϕexp(−2Z d3d4(cid:18)ϕ(3) (cid:19) (cid:18)−σr1o/d2 −G−co1il (cid:19)(cid:18) ϕ(4) (cid:19)). b (2.12) b ′ The choice of ψ(1)ϕ(1) under the path integral yields one rod-coil unit as a basic building block, see Fig.(2) and Eq.(2.14). This choice assigns the coordinates (1) and (1′) to (r,u) and (r′) respectively. Similarly, ψ(3) and ψ(4) are shorthand notations for ψ(r ,u ) and ψ(r ,u ). 3 3 4 4 The inversion of the 2 2-matrix in Eq.(2.11) reads × G−1 σ1/2 −1 1 G−1 σ1/2 σr1o/d2 −G−1 = G−1 G−1 σ σ1c/o2il G−1 . (2.13) (cid:18)−b coil (cid:19) rod∗ coil− (cid:18) b rod (cid:19) The calculation of the path integral in Ebq.(2.11) yieldbs the fbollowing result b σ1/2 Ξ(0)(1,1′;µ,ǫ,σ;[h ],[h ], hij ) = C R S G−1 G−1 σ h i rod∗ coil− = σ1/2G G 1+σG G b robd∗ coil∗ rod∗ coil h + σ2Grbod Gcboil Grobd Gcboil+..b. . ∗ ∗ ∗ i (2.14) b b b b The asterisk in Eq.(2.14) is a shorthand notation for a convolution of Green function operators. The geometric progressionwith a convolution as a binary relation has a clear pictorial representation, see Fig.(2). The first term of Ξ ( = 1, 1’ ; [h ], [h ], [h ] ) µ, σ, ε C R S + + . . . (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) FIG. 2: Pictorial representation of the geometric progression in Eq.(2.14). One rod-coil unit constitutes the basic building block. this series has the analytical expression σ1/2G G , whereas the ratio is equal to σG G , see also Fig.(3). rod coil rod coil The first term represents one rod-coil unit with∗one junction between rod and coil, hence∗the factor σ1/2. If an additional building block is added, two morebjunctiobns are created, hence the factor σ inbthe rabtio. This series gives a correct representation of the grand canonical polymeric correlationfunction. ^ ^ = σ (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) G G (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) * (cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:1) rod coil ^ ^ G G rod coil σ FIG. 3: Pictorial representation of the series ratio in Eq.(2.14): each bar (rod), zigzag line (coil) and fat dot correspond to Gbrod,Gbcoil , and σ1/2 respectively. 6 The denominator in Eq.(2.11) can be avoided by introducing de Gennes’ n 0 trick21. Consider the two n- → component vector fields ψ ,ϕ , where α=1,2,...n. Then Eq.(2.11) can be formally rewritten as α α { } n Ξ(0)(1,1′;µ,ǫ,σ;[h ],[h ], hij )= lim ψ ϕ ψ (1)ϕ (1′) C R S n→0 D αD α 1 1 h i αY=1Z 1 n ψ (3) T G−1 σ1/2 ψ (4) ×exp −2 d3d4 ϕαα(3) σr1o/d2 −G−1 ϕαα(4) . n Z αX=1(cid:18) (cid:19) (cid:18)−b coil (cid:19)(cid:18) (cid:19)o (2.15) b Integrations over the external fields h , h , hij, over the densities ρ (r),ρ (r), Sij(r) and over all endpoints 1,1′ C R S C R { } yield the final field theoretic representation of the grand canonical partition function n Z(µ,ǫ,σ) = lim ψ ϕ d3rd2uψ (r,u) d3r′ϕ (r′) α α 1 1 n→0 D D α=1Z (cid:20)Z (cid:21) (cid:20)Z (cid:21) Y 1 n b2(u )2 exp d3rd2uψ (r,u) β(µ ǫ) ·∇r ψ (r,u) α α × (−2 − − β(µ ǫ) αX=1Z h − i 1 n b2 d3rϕ (r) βµ 2 ϕ (r) − 2 α − 6 ∇r α αX=1Z h i n 2 χ d3r d2uψ2(r,u) − 4 " α # Z α=1Z X n n 2 v − 8Z d3r"α=1Z d2uψα2(r,u)+α′=1ϕ2α′(r)# X X n n 3 w − 48Z d3r"α=1Z d2uψα2(r,u)+α′=1ϕ2α′(r)# X X n + σ1/2 d3rd2uψ (r,u)ϕ (r) α α α=1Z X n n g − 4α=1α′=1Z d3rd2ud2u′P2(u·u′)ψα2(r,u)ψα2′(r,u′)), (2.16) X X where the second Legendre polynomial is given by 1 P (θ)= 3cos2θ 1 . (2.17) 2 2 − (cid:0) (cid:1) To evaluate the partition function Z(µ,ǫ,σ) the self-consistent field approximation is used22. In the self-consistent field approximationfluctuations are neglected and the functional integral over the fields in Eq.(2.16) is integratedby steepest descent. The saddle point solutions for ϕ (r) andψ (r,u) arechosensuchthat the effective grandpotential α α keeps the full symmetry of the Hamiltonian in replica space, i.e. it is invariantunder rotations in replica space. This is the case for (see ref.21) ψ (r,u) = n ψ(r,u) α α ϕ (r) = n ϕ(r), (2.18) α α where n is a unit vector such that n n2 =1. α=1 α In order to make the problem more tractable we expand ψ(r,u) in spherical harmonics, ψ(r,u) = ψ (r)Y (u), (see e.g. ref.2P3). Since the solution for ψ(r,u) must respect uniaxial and cylindrical symme- l,m lm lm try, the expansion reduces to Legendre polynomials, i.e. m 0. This expansion is truncated to the lowest nontrivial oPrder24 ≡ 1 1/2 5 1/2 ψ(r,u) ψ (r)+ ψ (r)P (u n ), (2.19) 0 2 2 z ≈ 4π 4π · (cid:18) (cid:19) (cid:18) (cid:19) 7 where n is the main direction along which the rods in the core of the globule are aligned, if the system forms an anisotropicglobulewithalignedrods. Iftheyarenotalignedandthesystemformsanamorphousglobule,onlythefirst term in the expansion in Eq.(2.19) differs from zero. On the other hand, the second term in Eq.(2.19) is responsible forthenematicLC-order. Themaindirectionofalignmentncanbechosenarbitrarilywithoutlossofgenerality,since a change in alignment direction only corresponds to a rotation of the complete globule in the laboratory coordinate frame. Thus we choose n directed along the z-axis of the (x,y,z) laboratory frame, i.e. n = n . This expansion z allows us to perform the u-integrations. The resulting effective saddle point grand potential Ω(µ,ǫ,σ) can now be calculated. It is given in terms of the saddle point fields ψ (r), ψ (r) and ϕ(r) by 0 2 β(µ ǫ) βΩ(µ,ǫ,σ) = − d3r ψ2(r)+ψ2(r) 2 0 2 Z b2 (cid:2) (cid:3) d3r 35ψ (r) 2ψ (r) − 210β(µ ǫ) 0 ∇r 0 − Z n + 14√5ψ (r) 2∂2 ∂2 ∂2 ψ (r) 0 z − x− y 2 h i + ψ (r) 25∂2+25∂2+55∂2 ψ (r) 2 x y z 2 1 h b2 i o + d3rϕ(r) βµ 2 ϕ(r) 2 − 6 ∇r Z h i + χ d3r ψ2(r)+ψ2(r) 2 4 0 2 Z + v d3r (cid:2)ϕ2(r)+ψ2(r)(cid:3)+ψ2(r) 2 8 0 2 Z + w d3r(cid:2) ϕ2(r)+ψ2(r)+ψ2(r(cid:3)) 3 48 0 2 Z (cid:2) (cid:3) 2√πσ d3rϕ(r)ψ (r) 0 − Z g 2 + d3r ψ (r) 7ψ (r)+√5ψ (r) . (2.20) 2 0 2 245 Z n h io The coil and rod densities as well as the orientation density are given by the following relations 1 ρ (r) = ϕ2(r) C 2 1 1 ρ (r) = d2uψ2(r,u) ψ2(r)+ψ2(r) R 2 ≃ 2 0 2 Z 1 (cid:2) (cid:3) S(r) Szz(r) = d2uP (cosθ)ψ2(r,u) 2 ≡ 3 Z 2 ψ (r) ψ (r)+√5ψ (r) . (2.21) 2 0 2 ≃ √5 h i Functional minimization of Eq.(2.20) with respect to ψ (r), ψ (r) and ϕ(r) yields a set of three coupled partial 0 2 differential equations, which are solved numerically with the finite element tool kit Gascoigne25. The results are extensively discussed in the next section. III. RESULTS In this section the numerical results for the full set of equations describing a rod-coil multiblock copolymer with a variable composition of stiff and flexible segments are presented. The segment length is set to b = 1. In addition, all energies such as ǫ, µ and also the saddle point free energy F are given in units of k T. In this section this will B not be indicated in order to avoid complicated notation. Eq.(2.20) describes the copolymer in a grand-canonical representation. In the grand-canonical ensemble the total number of segments of the polymer N is not fixed but its mean value is determined by equilibrium conditions. A real polymer has a fixed length. In order to ensure this fixed length N, the chemical potential µ is - for each set of physical parameters (v,w,χ,g,ǫ,σ) - tuned such that the equilibrium value of N is equal to the desired one. The 8 totalnumberofsegmentsN =N +N iscalculatedbynumericalintegrationoverthe rodandcoildensitiesgiven coil rod by Eq.(2.21) . For a given set of parameters, N(µ) can be computed24 and a typical example of this curve is shown in Fig.(4). For µ 0 the total number of segments N diverges. This corresponds to the N µ−1 behavior which is → ∼ 2 1.8 5 1.6 0 1 / N1.4 1.2 1 0.8 0.095 0.1 0.105 0.11 0.115 0.12 µ FIG. 4: N as a function of µ for w = 1, v = −0.2 and σ = 10−4. The dotted curve corresponds to ǫ= 0.08 and χ = 0. The continuouscurvecorresponds to ǫ=0.1 and χ=0.0138. well-knownfor a Θ-solventchain26. The divergence of N at a specific value µ on the righthand side of the minimum corresponds to a fully collapsed infinite globule and has been discussed first by Kholodenko and Freed21. These two branches of N(µ) meet each other in the minimum of N(µ) which can be associated with the coil-globule transition point, as will be shown in the next section. Since N is always fixed by tuning µ, it is possible to distinguish from a plotlike the one showninFig.(4)whether the systemis left ofthe coil-globuletransitionpoint (i.e. inthe openchain regime) or right of the transition point (i.e. in the globular regime). However, for fixed N, it is necessary to choose one of the two branches. Since this work focuses on the study of globular structures, the numerical calculations are alwaysrestrictedto the right branchincluding the minimum. The self-consistentfield theory is only expected to give good results for this branch since fluctuations are neglected. The three-body interaction parameter w is chosen to be w = 1 throughout the entire paper and the two-body interaction constant v is always negative to ensure that the system stays in the globular regime (up to the transition point). As a first step the pure homopolymer globule case21 has been studied in order to test the numerical routines. For this end one should simply set ψ = ψ = 0 in Eq.(2.20). Here we are not going to discuss these results in detail . 0 2 It is pertinent only to note that the resulting numerical solution fully supports the following well known theoretical results. The critical value of the two-body interaction constant, v, scales as v N−1/2, whereas the maximal c globule density in the critical point behaves as ρ N−1/2. | | ∼ crit ∼ Ashortremarkontheterminologythatwillbeusedbelowisnecessaryatthispoint. Thetermsphaseandtransition will be used frequently although the system is a polymer of finite length. All transitions are therefore crossovers of finite width with continuous order parameters. It is nevertheless common now as applied to ”soft matter” to use the theterm”phase”todistinguishthedifferentstructuralstatesofapolymerandtorefertothecrossoverbetweenthese states as a ”transition”. A. Coil-globule transition The rod-coilmultiblock copolymershowsa coil-globuletransitionsimilar to the one ofa homopolymer. To demon- stratethis, the interactionswhicharespecific for the stiffsegmentsaresetto zero,thatis χ=g =0. Inaddition, the energygainper stiff segmentis setto zero(ǫ=0)andit is assumedthat there is nocooperativityin the formationof stiff segments (σ =1). The two-body interaction constant v is varied. The transition point between coil and globule 9 is defined to be the minimum of the N(µ)-curve as it was explained above. The length of the polymer is fixed at N =550. Fig.(5)shows the profile ofthe totaldensityρ(r) ofthe copolymeras a functionof radialdistance fromthe center. The total density ρ(r) at each point is given by 1 1 1 ρ(r)=ρ (r)+ρ (r)= ϕ2(r)+ ψ2(r)+ ψ2(r). (3.1) C R 2 2 0 2 2 As can be seen from Fig.(5), the density profile becomes broader with decreasing v . At v = 0.5 the copolymer is | | − 1 0.8 = −0 5 v . 0.6 = −0 3 v . ρ 0.4 = −0 077 v . 0.2 5 10 15 20 25 r FIG. 5: This plot shows the radial density profile of the entire copolymer in radial direction for different values of v. The dashed curvefor v=−0.077 corresponds to thecoil-globule transition point. deepin the globularstate with abig plateauofalmostconstantdensity anda rathersmallsurfacelayerofdecreasing density. At the transition point v = 0.077 the plateau basically vanished and the surface layer becomes very broad. − To further illustrate the structural change, Fig.(6) shows a colour-coded plot of the local density in ̺-z space, where ̺ denotes the radialdirectionandz the axialdirectionincylindricalcoordinates. The center ofthe globuleis located at the bottom left corner of eachpicture. The pictures on the left show the localdensity of flexible segments and the pictures on the right correspondto the local density of stiff segments. Red indicates high density and dark blue zero density. In the upper two plots the copolymer is at the transition point (v = 0.077). In the lower two plots it is − deep in the globular state (v = 0.5). − Amuchclearerindicationthatv = 0.077correspondsindeedtoatransitionpointcanbeseenfromFig.(7),where − the globule radius is plotted versus v. The globule radius is defined as the point R in radial direction at which the density ρ(R) has decreased to ρ(R) = 10−3ρ , where ρ is the maximum density at the center of the globule. The 0 0 radius R shows a rapid increase when v = 0.077 is approached. Note, that the copolymer is finite (here: N=550) − and therefore all transitions are crossoversas discussed above. B. Helix-coil transition In this section the fraction of stiff segments is investigated as a function of the energy gain per stiff segment ǫ. This is similar to the helix-coil transition described by the Zimm-Bragg model4. A major difference is that the model used here is a three-dimensional model of the polymer including interactions, whereas the Zimm-Bragg model and its extensions are one-dimensional models with no three-dimensional interactions and no explicit entropy term . On the other hand, the Zimm-Bragg model can be solved exactly, whilst the self-consistent field treatment of the three-dimensional model is a mean-field approach which neglects fluctuations. Two different regimes will be discussed in the following: a low cooperativity regime with σ in the range 0.05 1 and a high cooperativity regime with σ in the range 7 10−3 10−4. Remember, that σ =1 means no cooperativ−ity · − and σ =0 means total cooperativity. Throughoutthis sectionχ and g are setto zero. There are therefore no specific 10 FIG.6: Thedensityofflexibleandstiff(helical)segmentsareshownontheleftandontherightpanelsrespectively. v=−0.077 corresponds to the upperplots and v=−0.5 to thelower plots. 40 35 30 25 R 20 15 10 0.1 0.2 0.3 0.4 0.5 v − FIG. 7: The globule radius R is plotted as a function of v. interactions between the stiff segments. The only interactions are attractive two-body interactions and repulsive three-body interactions between all segments. Fig.(8) shows how the fraction of stiff segments Θ depends on ǫ for R differentvaluesofσ intherange0.05 1. Evenforsmallcooperativity,the slopeclearlydepends onσ andgetslarger − with increasing cooperativity (decreasing σ).