Identifying transitions in finite systems by means of partition function zeros and microcanonical inflection-point analysis: A comparison for elastic flexible polymers
Julio C. S. Rocha, Stefan Schnabel, David P. Landau, Michael Bachmann
aa r X i v : . [ phy s i c s . d a t a - a n ] A ug Identifying transitions in finite systems by means of partition function zeros andmicrocanonical inflection-point analysis: A comparison for elastic flexible polymers
Julio C. S. Rocha, ∗ Stefan Schnabel, † David P. Landau, ‡ and Michael Bachmann
1, 3, 4, § Center for Simulational Physics, The University of Georgia, Athens, Georgia 30602, USA Institut f¨ur Theoretische Physik and Centre for Theoretical Sciences (NTZ),Universit¨at Leipzig, Postfach 100920, D-04009 Leipzig, Germany Instituto de F´ısica, Universidade Federal de Mato Grosso, Cuiab´a (MT), Brazil Departamento de F´ısica, Universidade Federal de Minas Gerais, Belo Horizonte (MG), Brazil (Dated: June 19, 2018)For the estimation of transition points of finite elastic, flexible polymers with chain lengths from13 to 309 monomers, we compare systematically transition temperatures obtained by the Fisherpartition function zeros approach with recent results from microcanonical inflection-point analysis.These methods rely on accurate numerical estimates of the density of states, which have beenobtained by advanced multicanonical Monte Carlo sampling techniques. Both the Fisher zerosmethod and microcanonical inflection-point analysis yield very similar results and enable the uniqueidentification of transition points in finite systems, which is typically impossible in the conventionalcanonical analysis of thermodynamic quantities.
PACS numbers: 05.10.-a,05.70.Fh,82.35.Lr
I. INTRODUCTION
Phase transitions are among the most fascinating phe-nomena in nature and huge efforts have been made tounderstand the features that characterize these coopera-tive processes for many different systems in a general andsystematic way. Strictly speaking, thermodynamic phasetransitions occur only in the thermodynamic limit, i.e.,for infinitely large systems. However, recent growing in-terest has also involved finite systems. Prominent repre-sentatives for such systems are finite polymer chains and,in particular, proteins. Because of surprisingly manifestcommon properties of transitions in finite and infinitesystems, the question arose to what extent the relation-ship between “pseudo-transitions” in finite systems andtheir infinite-system counterparts can be stressed. It iswell known that the precise determination of the locationof transitions in finite systems is typically ambiguous anddifferent fluctuating quantities suggest different points inparameter space as transition points. In the thermody-namic limit, scale freedom would let this space collapseto a single unique transition point. However, most con-temporary problems in soft condensed matter and tech-nology are apparently of small size, for which the ther-modynamic limit is not applicable at all. For this reason,it is necessary to verify if the methods of statistical anal-ysis that have been developed for infinitely large systemsand have proven to be so extremely successful in thesecases can be employed for, or adapted to, finite systems ∗ E-mail: [email protected] † E-mail: [email protected] ‡ § as well.Another important aspect is the fact that computersimulations open a completely new view on statisticalphysics, as only the most recently developed computa-tional methods and algorithms enable the accurate studyof fundamental statistical quantities that could hardly beapproached by theoretical methods in the course of theestablishment of the theory of complex phenomena andphase transitions in the past decades. One such quantityis the density of states g ( E ), i.e., the number of systemconfigurations within a given energy interval. Its loga-rithm can be associated with the entropy of the systemin energy space, S ( E ) = k B ln g ( E ), and the first deriva-tive with respect to energy yields the inverse temperature β ( E ) = dS ( E ) /dE . It has been shown recently that thecareful analysis of inflection points of this quantity re-veals all transitions in the system uniquely and withoutany ambiguity [1]. Since in this approach the tempera-ture is considered to be a derived quantity and a functionof energy, this method is a representative of microcanon-ical statistical analysis.In this paper, we will also make use of the densityof states, but we are going to interpret its features in acanonical way by considering the partition function Z ( T )of the system as a function of the (canonical) tempera-ture T . The thermodynamic potential associated withthe canonical ensemble (we consider fixed system size N and volume V ) is the free energy F ( T ) = − k B T ln Z ( T ).Thermodynamic phase transitions are located in temper-ature space, where a derivative of F of a certain orderexhibits a singularity [2–6]. Examples are the canonicalentropy S ( T ) = − ( dF ( T ) /dT ) N,V and response quanti-ties such as the heat capacity C V = T ( dS ( T ) /dT ) N,V = − T ( d F ( T ) /dT ) N,V . Yang and Lee were the first to re-late catastrophic singularities to partition function zerosin the grand canonical ensemble by introducing complexfugacities [7]. Fisher evolved this idea for the canonicalpartition function by introducing a complex temperatureplane [8].There is extensive literature on applications of suchmethods to various physical systems such as spin models(see, e.g., Refs. [9–11]), proteins [12, 13], and to poly-mers [14, 15]. Most applications of the partition functionzero analysis method are considered to be alternative ap-proaches to scaling properties near phase transitions inlarge systems. However, this method is also promisingfor the identification and characterization of analogs ofphase transitions in finite systems, in particular in finitelinear polymer chains that are known to exhibit a va-riety of structural transitions which sensitively dependon the chain length [1, 16–18]. The understanding ofthese structure formation processes is relevant from bothfundamental scientific and applied technological perspec-tives of molecular building-block systems.Typically these processes are accompanied by nucle-ation transitions, where crystalline shapes form from aliquid or vapor phase. Crystalline or glass-like structuresof single polymer chains can serve as the basic elementsof larger assemblies on nanoscopic scales; and beyondthat, the crystallization behavior exhibits strong simi-larities to the cluster formation of colloidal (or atomic)particles [17]. The nucleation is governed by finite-sizeand surface effects, where functionalization is based onthe individual structural properties of small moleculesforming large-scale composites [17]. These effects canbe analyzed by means of microcanonical thermodynam-ics [19], in which case transition properties can be de-rived directly and systematically from the caloric entropycurve [1]. This approach has been successfully appliedto a variety of structural transitions in macromolecularsystems such as folding [1, 20–22], aggregation [23], andadsorption processes of polymers and proteins [24, 25].One particular problem that has gained increased inter-est recently is the influence of the interaction range onthe stability of structural phases [21, 26]. This has beenaddressed by means of systematic microcanonical analy-ses in discrete and continuous polymer models.In principle, once the density of states g ( E ) is given,the partition function can easily be calculated and its ze-ros identified. However, examples of systems for which g ( E ) can be calculated exactly, or quite accurately bytheoretical methods, are very rare. It requires sophis-ticated numerical methods such as generalized-ensembleMonte Carlo sampling that allow for accurate estimatesof g ( E ). Among the most popular methods are mul-ticanonical sampling [27, 28] and the Wang-Landaumethod [29]. These methods are capable of scanning theentire phase space effectively in a single simulation.Compared to recent studies on partition function zeroanalyses of polymers such as Ref. [14], we here employa more realistic coarse-grained model for elastic, flexiblepolymers with continuous, distance-dependent monomer-monomer interactions based on van der Waals forces. Re-cently developed sophisticated simulation methodologiesspecific to this model [30] enable a very precise estimation of fundamental statistical quantities such as the densityof states. This is essential for the careful identificationof low-entropy phases that include liquid-solid and solid-solid transitions. For finite systems, these transitions arestrongly affected by finite-size effects, which are of partic-ular interest in this comparative study of advanced sta-tistical analysis methods. One major question is whetherthe partition function zeros method, which is effectivelya canonical approach, is capable of revealing the sameintricate details of these effects as the microcanonicalinflection-point analysis [1]. For this purpose, we system-atically analyze the canonical partition function zeros forall chain lengths ranging from 13 to 309 monomers in thismodel and identify and classify all structural transitions.Since the finite-size effects in the solid phases are surfaceeffects specific to the explicit chain length, transitionsin between them do not exhibit obvious scaling proper-ties [1, 17, 31]. Therefore, scaling considerations are notin the focus of this study.This paper is organized as follows: In Sec. II, we reviewthe partition function zeros approach and describe thenumerical methods used for the estimation of the densityof states and for the identification of the Fisher zeros.This section also includes a brief discussion of the mi-crocanonical inflection-point analysis. The results of ourstudy are presented in Sec. III, where we first discuss thedifferent scenarios in the liquid-solid and solid-solid tran-sition regimes thoroughly by investigating the zero mapsfor four representative examples that differ in the pro-cesses of Mackay and anti-Mackay overlayer formation.We then generalize and summarize the results obtainedby the zeros method for all polymers with chain lengthsup to 309 monomers and compare with former resultsobtained by microcanonical inflection-point analysis [1].The paper is concluded by a summary in Sec. IV. II. METHODS AND MODELA. Partition function zeros and thermodynamics
We consider a polymer system in thermal equilibriumwith a heat bath that is described by the canonical
N V T ensemble (constant particle number N , volume V , andtemperature T ). This ensemble connects microscopicquantities and thermodynamical properties via statisti-cal relations described by the canonical partition function Z . In thermal equilibrium, the probability for a discreteenergetic state is p m = g m e − βE m /Z , where g m denotesthe density of states at each energy E m ; β = 1 /k B T isthe inverse thermal energy and k B is the Boltzmann con-stant. In this work the units are chosen so that k B = 1.For a discrete ensemble of energetic states, the partitionfunction reads Z = X m g m e − βE m = e − βE X m g m e − β ( E m − E ) , (1) -5 -2.5 0 2.5 5 E l n ( g ) ε FIG. 1. Pictorial demonstration of the discretization of acontinuous density of states over an energy range which isdivided in n bins of size ε . Here the bins are labeled from 0to n −
1, thus the energy of the m th bin is E m = E + mε .All states with energy between E m and E m + ε are recordedin the m th bin g m . where we have extracted the Boltzmann factor of theground state for future convenience. All essential ther-modynamic quantities such as entropy and response func-tions like the heat capacity derive from the free energy F = − ln Z/β .For the subsequent analysis of a model with a con-tinuous energy spectrum, it is necessary to discretizethe density of states. Estimates obtained by means ofgeneralized-ensemble Monte Carlo methods such as mul-ticanonical [27, 28] and Wang-Landau sampling [29] arenaturally discrete in energy space (see Fig. 1). If the en-ergy bin size is chosen to be ε , the partition function (1)can be rewritten as Z = e − βE n − X m =0 g m e − βmε , (2)where n denotes the total number of bins.Defining x ≡ e − βε , the partition function can assumethe form of a polynomial Z = e − βE n − X m =0 g m x m = e − βE n − Y j =1 ( x − x j ) . (3)In the latter expression, the polynomial was decomposedinto linear factors ( x − x j ), where x j denotes the j th zero(or root) of the polynomial. With the polynomial definedin this way, the density of states can cover the entire spaceof energy for both positive and negative energies. Notethat x ≥
0; if T → x →
0, whereas x →
1, if T → ∞ .In Eq. (3), Z is written as a polynomial of degree n − n −
1, generally complex, roots. Since Z ∈ ℜ and for a finite system always Z > g m are nonzero positive real numbers, theroots must occur as complex conjugate pairs a j ± ib j with a, b ∈ ℜ . Real-valued roots must be negative. Once the partition function is determined thermody-namic quantities can be extracted from the the Helmholtzfree energy F . The internal energy is U = h E i = − ∂ ln Z∂β (4)and, most interesting for the following consideration, thespecific heat at constant volume reads c V = 1 N (cid:18) ∂U∂T (cid:19) V = k B β N ∂ ln Z∂β . (5)Inserting the factorization (3), these quantities can alsobe expressed by the Fisher zero components: U = E + n − X j =1 (cid:18) εxx − x j (cid:19) = E + n − X j =1 εx ( x − a j )( x − a j ) + b j ! , (6)and c V = k B x (ln x ) N n − X j =1 (cid:18) − x j ( x − x j ) (cid:19) = k B x (ln x ) N n − X j =1 − a j ( x − a j ) + b j (2 x − a j )[( x − a j ) + b j ] ! . (7)Obviously, this expression can only become singular at x = a j , if b j = 0, i.e., if the j th zero lies on the positivereal axis. According to Yang and Lee, zeros that comearbitrarily close to the real axis in the thermodynamiclimit mark the transition points. This is essential for ourstudy as we are interested here exclusively in transitionproperties of polymers of finite length. Therefore, we donot expect to find any real-valued zeros in the analysis ofthe complex-zero space of these systems. Rather, we willidentify the zeros closest to the positive real axis whichare called the leading zeros because they contribute mostto the quantity of interest, if x ≈ a j . If such zeros have arather isolated appearance in the distribution of the ze-ros in the complex map near the positive real axis, theyrepresent a signal in that quantity that might become asingularity in the infinitely large system. At least, in thefinite system, they indicate increased thermal activity.Canonical quantities such as the specific heat typicallypossess a peak or a “shoulder” in those regions in tem-perature space.Technically, apart from finite-size scaling, there are twopossibilities to define transition points for finite systemsby means of partition function zeros. Either one con-siders the zero as if it lies on a circle (in first-order liketransitions, the transition-state zeros distribute indeednear a circular line), in which case the radius defined via | x j | = a j + b j can be used to locate the intersection pointon the positive real axis: x c ≡ a ′ j = | x j | . Alternatively,since b j will be small near the positive real axis, one canalso simply choose x c = a j ≈ | x j | . Either way, by per-forming the projection upon the real axis, a specific-heatsingularity is mimicked even for a finite system . Thetransition point can then be defined by T c = − εk B ln | x j | . (8)On this basis, conclusions about the structural transi-tions of finite-length flexible polymers will be drawn inthis study, but these transitions should not be confusedwith the strictly defined thermodynamic phase transi-tions in the Yang-Lee sense.The accurate estimation of the partition function ze-ros requires two separate parts that for a complex sys-tem can only be accomplished computationally. First,generalized-ensemble Monte Carlo simulations have tobe performed to obtain the density of states. Second,all zeros of the polynomial form of the partition functionmust be identified. Since a polynomial of degree five orhigher has no algebraic solution in general, as stated bythe Abel-Ruffini theorem, the zeros can only be foundby means of numerical computation. We will review thepolymer model and the simulation and analysis methodsused in the following. B. Coarse-grained polymer model
A linear polymer of length L is formed by concatena-tion of L identical chemical units called monomers. Eachmonomer is composed of several atoms, thus the size ofthe chain suitable for simulation is limited by the compu-tational resources and methods currently available. Forthe study of generic thermodynamic properties of poly-mers, however, all-atom models can typically be replacedby a simpler coarse-grained representation with effectiveinteractions. We here consider such a generic coarse-grained model for linear, elastic, flexible polymers [16].Non-bonded monomers interact pairwise via a truncatedand shifted Lennard-Jones (LJ) potential V modLJ ( r ij ) = V LJ (min( r ij , r c )) − V LJ ( r c ) , where r ij denotes the distance between the i th and the j th monomer, r c is the cutoff distance, and V LJ ( r ) = 4 ǫ (cid:20)(cid:16) σr (cid:17) − (cid:16) σr (cid:17) (cid:21) is the standard LJ potential. In this work the LJ param-eters were chosen as ǫ = 1, σ = 2 − / r , and r c = 2 . σ .The elastic bonds between monomers adjacent alongthe chain are modeled by the finitely extensible nonlinearelastic (FENE) potential [32] V FENE ( r ii +1 ) = − K R ln " − (cid:18) r ii +1 − r R (cid:19) . This potential possesses a minimum at r and divergesfor r → r ± R . K is a spring constant and we set theparameters as R = 0 . r = 0 .
7, and K = 40. C. Numerical methods
1. Monte Carlo sampling in a generalized ensemble
Since the simulation of structural phases of polymersis challenging, even for a coarse-grained model and mod-erate system sizes, a sophisticated advanced Monte Carloupdate set [30] was applied in combination with multi-canonical sampling [27, 28, 30]. The majority of movesconsisted of attempted displacements of single monomerswithin a sphere around their original location. Depend-ing on energy E and number of monomers N the radiiof these spheres were chosen such that high acceptancerates could be achieved for all energies and system sizes.In addition, we used bond-rebridging moves, where allmonomers keep their position, but the linkage betweenthem is altered. Furthermore, a novel cut-and-pastemove was developed in which one monomer is removedand reinserted in an entirely different location within thepolymer chain.Most of the data were produced in a single simula-tion by sampling a generalized “grand-multicanonical”ensemble [30]. The main goal was to avoid free energybarriers by enabling the system to change its size. There-fore, in addition to the trial update schemes describedabove, a Monte Carlo step was introduced by means ofwhich single monomers could randomly be added or re-moved. A weight function W ( E, N ) assured that all en-ergies and sizes were visited with the same probability.It was tuned using a delayed Wang-Landau procedure,in which the modification factor of the original Wang-Landau method is made weight-dependent. If the mul-ticanonical weight function at Monte Carlo “time” t isdenoted by W t , then it is modified after the next updateto W t +1 ( E, N ) = W t ( E, N ) /f W t ( E,N ) /W t − d ( E,N ) (9)for E = E t − d , N = N t − d . For other values of E and N , the weight remains unchanged as in a conventionalmulticanonical simulation. Therefore, the effect of theWang-Landau modification factor f to smooth out thefree-energy landscape is delayed by d . This slows downthe saturation speed of Wang-Landau sampling and en-ables a better efficiency in exploring phase space regionsof low entropy at low energy, in particular in isolated re-gions that might contain hidden barriers. For the poly-mer system considered here, this is particularly relevantin the solid-solid transition regime. A sufficiently largedelay for the polymer model considered here is obtainedby the choice d = 10 .Once the weights had converged data were generatedin a grand-multicanonical production that consisted ofapproximately 2 × Monte Carlo moves and consumedabout 0.5 CPU years.
2. Zeros finder
Computing the zeros of polynomials can be posed asan eigenvalue problem [33, 34]. Consider the matrix pair( A , B ) where A = · · · − g · · · − g · · · − g · · · − g ... ... ... . . . ... ...0 0 0 · · · − g n − (10)is the Frobenius companion matrix related to a monicpolynomial [35] of degree n [36], and B = · · · · · · · · · · · · · · · g n . (11)Then a straightforward computation shows thatdet ( x B − A ) = X m g m x m = P ( x ) . (12)On the other hand, the well-known generalized eigenvalueproblem (GEP) [37] can be stated asdet ( λ B − A ) = 0 . (13)By comparing Eqns. (12) and (13) one finds that eigen-values of the matrix pencil ( A , B ) are the zeros of P ,i.e., x k = λ k . The GEP can be solved by the QZ algo-rithm [38], just after performing a balance on the matrixpair ( A , B ), which is very important for accuracy [39–41].Both of these algorithms can be found in LAPACK [42].Alternatively, as implemented in Mathematica [43], onecan write a companion matrix of P as C = − g / g − g / g − g / g · · · − g n − / g − g n / g · · · · · · · · · · · · . (14)Then the zeros of P are obtained directly by diagonal-ization of C and given by x k = 1 λ k . (15)This method is more time consuming but also more ro-bust than the previous one.We employed both methods for the estimation of thepartition function zeros (3).
3. Microcanonical inflection-point analysis
An alternative approach to unravel transition prop-erties of finite-size systems is the direct microcanonicalanalysis [19] of caloric quantities derived from the entropy S ( E ) = k B ln g ( E ). The basic idea is that the interplay ofenergy and entropy and, in particular, changes of it, sig-nal cooperative system behavior that can be interpretedas a transition (and in the thermodynamic limit as aphase transition) of the system. Then first and higherderivatives of S ( E ) reveal the transition points of thesystem in energy space. However, since the first deriva-tive is the reciprocal microcanonical temperature, β ( E ) ≡ T − ( E ) = (cid:18) ∂S ( E ) ∂E (cid:19) N,V , (16)energetic transition points can also be associated withtransition temperatures. Transitions occur, if β ( E ) re-sponds least sensitively to changes in the energy. Theslope of the corresponding inflection points can be usedto distinguish first- and second-order transitions system-atically. If γ ( E ) = (cid:18) ∂β ( E ) ∂E (cid:19) N,V = (cid:18) ∂ S ( E ) ∂E (cid:19) N,V (17)exhibits a positive-valued peak at the inflection point,the transition resembles a first-order transition, whereasa negative-valued peak indicates a second-order transi-tion. This method is called microcanonical inflection-point analysis [1]. In the following, we will compare thetransition temperatures obtained from the leading zeroswith microcanonical estimates.
III. RESULTS AND DISCUSSION
Based on the density of states estimates obtained inmulticanonical simulations, we calculated the partitionfunction zeros for the elastic flexible polymer modelfor chain lengths L ranging from 13 to 309 monomers.The structural transition behavior was investigated pre-viously by conventional canonical statistical analysis of“peaks” and “shoulders” of fluctuating energetic andstructural quantities as functions of the canonical tem-perature [16, 17]. Subsequently, the densities of statesof this set of polymers were analyzed systematicallyby means of microcanonical inflection-point analysis,with particular focus on the typically hardly accessiblelow-temperature transition behavior (freezing, solid-solidtransitions) [1]. The microcanonical analysis is based onestimates of the microcanonical entropy and its deriva-tives , and therefore requires highly accurate data. There-fore, it is not only interesting from the statistical physicspoint-of-view to study the partition function zeros, butalso for practical purposes. The major information aboutstructural transitions is already encoded in the corre-sponding leading zeros which are rather simple to iden-tify. The partition function zero method thus turns out -1 -0.5 0 0.5 1 Re(x) -1-0.500.51 I m ( x ) -1 -0.5 0 0.5 1 Re(x) -1-0.500.51 I m ( x ) (a) (b) -1 -0.5 0 0.5 1 Re(x) -1-0.500.51 I m ( x ) -1 -0.5 0 0.5 1 Re(x) -1-0.500.51 I m ( x ) (c) (d)FIG. 2. Complex plane map of the partition function zeros for chain size: (a) L = 35, (b) L = 55, (c) L = 90, and (d) L = 300.The leading zeros are highlighted as follows: From x = 0 to 1 green squares denote “solid-solid” transitions, magenta diamondsdenote “liquid-solid” transitions, and blue circles denote “gas-liquid” transitions. to be a robust method for the identification of transi-tion points. It is, therefore, highly interesting to verifywhether this method is capable of finding indications forthe same transitions that have already been identified bymeans of microcanonical inflection-point analysis.Figure 2 shows the distributions of the zeros identifiedfrom the discretized densities of states for specific chainlengths L = 35 , , , and 300 and using the energy binsizes ε = 0 . , . , .
20, and 0 .
29, respectively. It isworth noting that the zeros, and thus their distribution,do generally depend on the choice of ε , but the transitiontemperature estimates remain widely unaffected if ε is changed. Moreover, since the data series used for theestimation of the density of states are finite, differentsimulation runs yield different values of the zeros.Note that we plot the zeros differently than Ref. [14].In our case they are strictly confined within a circle withradius 1 (the boundary at 1 corresponds to infinite tem-perature). We also define the transition temperature dif-ferently for a finite system. Ref. [14] considers only thereal part of the leading zero, whereas we prefer the abso-lute value, motivated by the fact that at first-order tran-sitions the zeros lie on a circle whose radius is a uniqueestimator for the transition temperature. Re(x) -0.2-0.100.10.2 I m ( x ) FIG. 3. Zoom into the zeros map for L = 35. Black circlesand red triangles represent the zeros obtained in two differentsimulations. Whereas the positions of nonleading zeros vary,the leading zeros are very close to each other and the overalldistribution pattern is very similar. The blue squares repre-sent the average values of the leading zeros over ten differentsimulations. Error bars are shown for the leading zero thatcorresponds to the liquid-solid transition; in the other casesthe error is smaller than the symbol size. The section of the map for L = 35 shown in Fig. 3 con-tains sets of zeros obtained in two independent simula-tions (circles and triangles). By standard jackknife erroranalysis [44–48], the statistical error of the componentsof the complex zeros was estimated from ten independentsimulations and error bars are shown for the leading zeros(squares) only (if larger than symbol size). Thus, for theanalysis of transitions, the method is sufficiently robustand enables the identification of transition points.We only analyze here the zero maps for L = 35 , , , and 300, because these system sizes are representative forthe various transition behaviors that have been system-atically and uniquely identified for polymer chains withlengths in the above mentioned interval in canonical [16–18] and microcanonical analyses [1]. From these studiesit is known that in this model polymers with “magic”length L = 13 , , , , . . . possess a second-order-like collapse (“gas-liquid”) transition and a very strongfirst-order-like freezing or “liquid-solid” transition fromthe compact, globular liquid phase into an almost per-fect icosahedral Mackay structure [49], where the facetsare arranged as fcc overlayers. For intermediate chainlengths, the optimal packing in the solid phase can beMackay or anti-Mackay (hcp overlayers), depending onthe system size and the temperature. In other words,for certain groups of chain lengths, an additional “solid-solid” transition can be found, in which anti-Mackayoverlayers turn into energetically more preferred Mackayfacets at very low temperatures [1, 16–18]. This behaviorof finite particle systems is also well known from atomicclusters [50–53].For the systems explicitly discussed here, this meansthat we expect to find three transitions for L = 35 and 90, whereas the solid-solid transition is absent for L = 55.For L = 300, the liquid-solid and the solid-solid transi-tion merge and occur at about the same temperature.These transitions can be distinguished microcanonically,but not canonically. Therefore, we do not expect to findindications of separate transitions in the analysis of theleading zeros.As earlier analyses revealed [1, 16], the liquid-solidand solid-solid transitions for system sizes 31 ≤ L ≤ L = 38 that forms a truncated fcc octahedron, these poly-mers crystallize in two different ways by cooling downfrom the liquid phase [16]. With high probability, morethan one icosahedral nucleus crystallizes out of the liquidby forming anti-Mackay overlayers and by an additionalsolid-solid transition turns into a single icosahedral nu-cleus with 13 monomers and a Mackay overlayer formedby the remaining ones. Alternatively, with lower proba-bility, the anti-Mackay multi-core structure can also formout of the liquid via an intermediate unstable phase dom-inated by a single-core structure with Mackay overlayer.Therefore, the anti-Mackay solid phase is a mixed phasethat also contains Mackay morphologies. Therefore theliquid-solid transition for these system sizes does not ex-hibit the same characteristic as for larger polymers andis actually second-order-like [1]. To conclude, all threestructural transitions for L = 35 are second-order-like.The corresponding zero maps shown in Figs. 2(a) and 3indeed reveal three separate pairs of leading zeros thatrepresent these transitions.The polymer chain containing 55 monomers is“magic”. For this reason, it exhibits a particularly strongliquid-solid transition at T ≈ .
33 into a perfect icosahe-dral conformation [16] with complete Mackay overlayer.A stable anti-Mackay phase does not exist and, there-fore, no solid-solid transition occurs. Consequently, thezero map shown in Fig. 2(b) reveals only two sets of lead-ing zeros representing the Θ collapse and the nucleationtransition. The most striking feature is the observationthat there is an increased accumulation of zeros on a cir-cle that contains the pair of the leading zeros associatedwith the liquid-solid transition. The circular distributionhas to be attributed to the self-reciprocity of the parti-tion function polynomial [54] at a phase transition withcoexisting phases in which case the energetic canonicaldistribution is bimodal and virtually symmetric. There-fore, the circular pattern can be interpreted as the sig-nature of first-order-like transitions in the map of Fisherpartition function zeros.For the polymer with L = 90 monomers, the structuraltransitions can clearly be identified in the correspondingzeros map [Fig. 2(c)]. The liquid-solid transition into theanti-Mackay solid phase is represented by a circular zerosdistribution, but neither the collapse transition nor thesolid-solid crossover to icosahedral Mackay structures ex-hibit obvious features in the zero distribution other thanprominent locations of the leading zeros. In correspon-dence with the previous microcanonical analysis, these T c V From zerosDirectly from g T ∆ c V ( x10 - ) T c V From zerosDirectly from g T ∆ c V ( x10 - ) (a) (b) T c V From zerosDirectly from g T ∆ c V ( x10 - ) T c V From zerosDirectly from g T ∆ c V ( x10 - ) (c) (d)FIG. 4. Heat capacity curves for chain sizes: (a) L = 35, (b) L = 55, (c) L = 90, and (d) L = 300. Plotted are the curvesobtained from the zeros of the partition function and, for comparison, by direct calculation from the density of states. Theinset shows the relative differences between them. The small deviations make it evident that all zeros were identified correctly.The vertical lines are located at the transition temperatures calculated from the leading zeros. Dashed and solid lines representfirst- and second-order-like transitions, respectively. transitions are classified as of second order. It is worthmentioning that the chain length L = 90 is close to thethreshold length ( L ≈ L >
110 until the next “magic” limit L = 147 isreached [1, 17], i.e., the Mackay phase is the only stablesolid phase. Microcanonically speaking, the solid-solidtransition lies energetically within the latent heat intervalof the first-order liquid-solid transition and can no longerbe resolved in the canonical analysis (the specific heat ex-hibits only one sharp peak in these cases [17]). The zerosmap shown in Fig. 2(c) reveals a very pronounced circu-lar distribution and the projected intersection point withthe positive x axis corresponds indeed to the liquid-solidtransition temperature.While L = 90 is a length below the anti-Mackay– Mackay threshold, our last example L = 300, is abovethe corresponding threshold in the following segmentof chain lengths that lies between two magic lengths,147 < L ≤
309 ( L = 309 is the next “magic” chainlength). The most surprising feature is that in tempera-ture space liquid-solid and solid-solid transitions merge,whereas energetically both can be distinguished clearly asfirst-order-like transitions [1]. The trend is that the solid-solid transition will shift to higher microcanonical tem-peratures than the liquid-solid transition when increasing L towards L = 309. This microcanonical crossover be-havior has already been known in other systems and is apure finite-size effect [26]. The corresponding root mapshown in Fig. 2(d) displays only the general canonicalbehavior; therefore, only one circle represents this first-order-like double-transition.For the explicit estimation of the transition tempera-tures from the Fisher zeros according to Eq. (8), thereis the ambiguity to use either the absolute values of the TABLE I. Comparison of transition temperatures for solid-solid (ss), liquid-solid (ls), and gas-liquid (gl) transitions for L =35 , , , and 300 as obtained by the partition function zero method ( T z ) and by microcanonical inflection-point analysis ( T m ).These estimates are compared to peak positions of the heat-capacity curves ( T ss , ls c V ) and fluctuations of the radius of gyration( T gl d h R i /dT ), respectively. The maximum 1 σ tolerance of all estimates is ± L T ssz T ssm T ss C V T lsz T lsm T ls C V T glz T glm T gl d h R i /dT
35 0 .
15 0 .
14 0 .
14 0 .
39 0 .
39 0 .
38 1 .
39 1 .
39 1 . .
33 0 .
33 0 .
33 1 .
53 1 .
51 1 . .
26 0 .
26 0 .
27 0 .
33 0 .
33 0 .
33 1 .
68 1 .
65 1 . .
44 N/A 0 .
43 0 .
43 0 .
43 1 .
97 1 .
88 1 . complex zeros or their real parts only, T tr = − εk B ln ( a j + b j ) ≈ − εk B ln a j . (18)Both values differ for finite systems, but converge in thethermodynamic limit. Since we already know that dis-tributions of zeros for first-order-like transitions are cir-cular, we chose to define transition points by means ofthe absolute values (corresponding to the radius of thecircle). For the four examples that we discuss here inmore detail, the corresponding values are listed in Ta-ble I. These estimates are in very good agreement withthe transition temperatures obtained by microcanonicalanalysis. Since the Θ transition is only represented bya weak shoulder in the heat capacity curves shown inFig. 4, we consider in these cases the corresponding peakpositions of the fluctuations of the radius of gyration, d h R gyr i /dT as a more appropriate indicator of these tran-sitions. This is a general problem of the canonical analy-sis of fluctuating quantities and the major reason for theintroduction of methods that enable a unique identifica-tion of transition points even for finite systems.For this reason both the zeros method and the mi-crocanonical inflection-point analysis are more useful forthe definition of unique transition temperatures than theconventional approach of the quantitative analysis of fluc-tuating quantities. Furthermore, the analysis of zero dis-tributions or microcanonical inflection points allow thediscrimination between first- and second-order-like tran-sitions. This information is not easily accessible from or-dinary canonical statistical analysis. In Figs. 4(a)- 4(d),vertical lines are located at the positions of the transitiontransitions obtained by the analysis of the Fisher zeros.Fig. 5 summarizes our results of the Fisher zero anal-ysis for all chain lengths in the interval 13 ≤ L ≤ L ≥
55) that exhibit more stable structuralphases.
IV. SUMMARY
We calculated Fisher partition function zeros for ageneric model of flexible, elastic polymers on the basisof accurate estimates of the densities of states for chain
16 32 64 128 256 T T MicrocanonicalLeading zeros
FIG. 5. Transition temperatures of conformational transi-tions for elastic, flexible polymers with chain lengths rangingfrom L = 13 to 309. The black dots represent the transitiontemperatures obtained from the leading zeros of the parti-tion function. For comparison, the transition temperaturesobtained by microcanonical inflection-point analysis are alsoshown (red triangles). ≤ L ≤ ACKNOWLEDGMENTS
We would like to thank Shan-ho Tsai, Lucas A. S.M´ol, and Bismarck V. Costa for inspiring discussions.JCSR thanks CNPq (National Council for Scientificand Technological Development, Brazil) for the supportunder grant No. PDE 202122/2011-5. MB acknowl-edges partial support of this project by the NSF underGrant No. DMR-1207437 and by CNPq under Grant No.402091/2012-4. [1] S. Schnabel, D. T. Seaton, D. P. Landau, and M. Bach-mann, Phys. Rev. E , 011127 (2011).[2] P. Ehrenfest and T. Ehrenfest, The Conceptual Founda-tion of the Statistical Approach in Mechanics (CornellUniversity Press, Ithaca NY, 1959), reprint of the Ger-man 1912 edition.[3] L. D. Landau, Nature , 840 (1936).[4] M. E. Fisher, Rep. Prog. Phys. , 615 (1967).[5] L. P. Kadanoff, W. G¨otze, D. Hamblen, R. Hecht, E. A. S.Lewis, V. V. Palciauskas, M. Rayl, J. Swift, D. Aspnes,and J. Kane, Rev. Mod. Phys. , 395 (1967).[6] P. Heller, Rep. Prog. Phys. , 731 (1967).[7] C. N. Yang and T. D. Lee, Phys. Rev. , 404 (1952);T. D. Lee and C. N. Yang, Phys. Rev. , 410 (1952).[8] M. E. Fisher, in Lectures in Theoretical Physics vol. 7C ,ed. by W. E. Brittin (University of Colorado Press, Boul-der, 1965), Chap. 1.[9] W. Janke and R. Kenna, J. Stat. Phys. , 1211 (2001);Comp. Phys. Comm. , 443 (2002); Nucl. Phys. B(Proc. Suppl.) , 905 (2002).[10] C.-O. Hwang, Phys. Rev. E , 042103 (2009).[11] J. S. M. Fonseca, L. G. Rizzi, and N. A. Alves, Phys.Rev. E , 011103 (2012).[12] N. A. Alves and U. H. E. Hansmann, Phys. Rev. Lett. , 1836 (2000).[13] J. Wang and W. Wang, J. Chem. Phys. , 2952 (2003).[14] M. P. Taylor, P. P. Aung, and W. Paul, Phys. Rev. E ,012604 (2013).[15] C.-N. Chen, Y.-H. Hsieh, and C.-K. Hu, EPL , 20005(2013).[16] S. Schnabel, T. Vogel, M. Bachmann, and W. Janke, Chem. Phys. Lett. , 201 (2009).[17] S. Schnabel, M. Bachmann, and W. Janke, J. Chem.Phys. , 124904 (2009).[18] D. T. Seaton, T. W¨ust, and D. P. Landau, Phys. Rev. E , 011802 (2010).[19] D. H. E. Gross, Microcanonical Thermodynamics (WorldScientific, Singapore, 2001).[20] T. Chen, X. Lin, Y. Liu, H. Liang, Phys. Rev. E , 114907 (2009); Phys. Rev. E , 050801(R) (2009).[22] M. Bachmann, Phys. Proc. , 1387 (2010).[23] C. Junghans, M. Bachmann, and W. Janke, Phys. Rev.Lett. , 218103 (2006); J. Chem. Phys. , 085103(2008); Europhys. Lett. , 40002 (2009).[24] L. Wang, T. Chen, X. Lin, Y. Liu, H. Liang, J. Chem.Phys. , 244902 (2009).[25] M. M¨oddel, W. Janke, and M. Bachmann, PhysChem-ChemPhys , 11548 (2010).[26] J. Gross, T. Neuhaus, T. Vogel, and M. Bachmann,J. Chem. Phys. , 074905 (2013).[27] B. A. Berg and T. Neuhaus, Phys. Lett. B , 249(1991); Phys. Rev. Lett. , 9 (1992).[28] W. Janke, Physica A , 164 (1998); B. A. Berg, Comp.Phys. Commun. , 397 (2003); M. Bachmann, Phys.Scr. , 058504 (2013).[29] F. Wang and D. P. Landau, Phys. Rev. Lett. , 2050(2001); Phys. Rev. E , 056101 (2001); D. P. Landau,S.-H. Tsai, and M. Exler, Am. J. Phys. , 1294 (2004).[30] S. Schnabel, W. Janke, and M. Bachmann, J. Comput.Phys. , 4454 (2011). [31] T. Vogel, M. Bachmann, and W. Janke, Phys. Rev. E , 061803 (2007).[32] R. B. Bird, C. F. Curtiss, R. C. Armstrong, and O. Has-sager, Dynamics of Polymeric Liquids , 2nd ed. (Wiley,New York, 1987); K. Kremer and G. S. Grest, J. Chem.Phys. , 5057 (1990); A. Milchev, A. Bhattacharaya,and K. Binder, Macromolecules , 1881 (2001).[33] P. Van Dooren and P. Dewilde, Lin. Alg. Appl. , 545(1983).[34] W. H. Press, S. A. Teukolsky, W. T. Vetterling, andB. P. Flannery, Numerical Recipes in Fortran 77: TheArt of Scientific Computing , 2nd ed. (Cambridge Uni-versity Press, New York, 1992), Chap. 9.5.[35] In a monic polynomial of degree n , the n th coefficient isunity, i.e., g n = 1.[36] H. Linden, Lin. Alg. Appl. , 41 (1998).[37] D. S. Watkins, “The Generalized Eigenvalue Problem”,in The Matrix Eigenvalue Problem (SIAM, 2007), p. 233.[38] C. Moler and G. Stewart, SIAM J. Numer. Anal. , 241(1973).[39] B. Parlett and C. Reinsch, Numer. Math. , 293 (1969).[40] R. Ward, SIAM J. Sci. Stat. Comp. , 141 (1981).[41] D. Lemonnier and P. Van Dooren, SIAM J. Matrix Anal.Appl. MathWorld – A Wolfram Web Resource , http://mathworld.wolfram.com/PolynomialRoots.html .[44] R. G. Miller, Biometrika , 1 (1974).[45] B. Efron, The Jackknife, the Bootstrap, and Other Re-sampling Plans (SIAM, Philadelphia, 1982).[46] W. Janke, “Statistical Analysis of Simulations: DataCorrelations and Error Estimation”, in
Proceedings of theWinter School “Quantum Simulations of Complex Many-Body Systems: From Theory to Algorithms” , John vonNeumann Institute for Computing, J¨ulich, NIC Seriesvol. , ed. by J. Grotendorst, D. Marx, and A. Mura-matsu (NIC, J¨ulich, 2002), p. 423.[47] D. P. Landau and K. Binder, A Guide to Monte CarloSimulations in Statistical Physics , 4th ed. (CambridgeUniversity Press, New York, 2014), in press.[48] M. Bachmann,
Thermodynamics and Statistical Mechan-ics of Macromolecular Systems (Cambridge UniversityPress, New York, 2014).[49] A. L. Mackay, Acta Crystallogr. , 916 (1962).[50] J. A. Northby, J. Chem. Phys. , 6166 (1987).[51] D. J. Wales and J. P. K. Doye, J. Phys. Chem. A ,5111 (1997).[52] E. G. Noya and J. P. K. Doye, J. Chem. Phys. ,104503 (2006).[53] P. A. Frantsuzov and V. A. Mandelshtam, Phys. Rev. E , 037102 (2005).[54] W. Chen, J. Math. Anal. Appl.190