# 1 Introduction

## 1 Introduction

The detailed structure of bulk amorphous polymers is a topic of scientific interest because it is necessary for the microscopic understanding of their properties. However, because of the amorphous nature of polymer melts and glasses, structural information is difficult to obtain experimentally.

In particular, it is interesting to know how amorphous a polymer melt is on a local scale, i.e. how much residual order is left on a local scale, how far local order extends before it disappears into the long-range disorder of amorphous systems , and how the local order depends on the molecular architecture. The interest has recently been revived by solid-state NMR studies of Graf et al. , from which it was inferred that a melt of polybutadiene is far more ordered than hitherto expected. The alignment of polymer chains is restricted to a local scale, there is no sign of nematic ordering.

In order to get a better understanding of local packing and ordering effects computer simulations are very helpful, because the system is precisely known and because one has access to all data including positions and velocities of all particles at all times. Atomistic simulations may be useful in order to get an understanding of a specific system whereas simplified models yield the properties of generic polymer melts. Additionally, they need much less simulation time which allows to tackle relatively big systems for long times . Therefore, a simple bead-spring model may be a good starting point to elaborate generic packing effects.

There was some work done for semiflexible chains, both by Monte Carlo and molecular dynamics, mostly to study liquid crystals  or focusing on confined systems . The influence of chain stiffness on the dynamic structure factors of polymer melts was investigated also by analytical theory by Harnau et. al  who found major discrepancies to the fully flexible system for large scattering vectors, i.e. on short distances.

In a previous article , we showed that there is considerable local chain alignment even in melts of fully flexible chains (persistence length: 1 monomer diameter). This persistence length originates from excluded volume interaction. If there was no interaction at all (except for connectivity) the persistence length would be zero (e.g. polycatenans). In the present contribution, this model is extended to include some more information about the chemical architecture of the polymer. We firstly introduce bending potentials of different strength, in order to study the effect of semiflexibility of single chain structure as well as on the mutual local orientation of neighboring chains. Secondly, we study models with alternating stiffness in an attempt to mimic simplistically polymers with rigid subunits connected by more flexible links, like polybutadiene and polyisoprene with their alternating single and double bonds which are currently under investigation experimentally .

## 2 Simulated System

We performed molecular dynamics simulations (for details of the parallel program POLY, see ref. ) of melts of polymer chains at a density $`\rho ^{}=0.85`$ and temperature $`T^{}=1`$ at a timestep $`\delta t^{}=0.01`$ using a truncated and shifted Lennard-Jones potential (Weeks-Chandler-Anderson potential) for the excluded-volume interaction between all beads. Lennard Jones reduced units are used throughout this paper where the mass $`m`$, the potential well depth $`ϵ`$ and the radius of the potential minimum $`\sigma `$ define the unit system.

$$V_{LJ}(r)=4ϵ\left[\left(\frac{\sigma }{r}\right)^{12}\left(\frac{\sigma }{r}\right)^6\right]+ϵ,r<r_{cutoff}=\sqrt[6]{2}\sigma $$
(1)   
and a finitely extendable non-linear elastic (FENE) potential

$$V_{FENE}(r)=\frac{\alpha }{2}\frac{R^2}{\sigma ^2}\mathrm{ln}\left(1\frac{r^2}{R^2}\right),r<R=1.5\sigma ,\alpha =30$$
(2)   
for the connection of neighboring beads. Additionally, a bond angle potential

$$V_{angle}=x\left(1\frac{𝐫_{i1,i}𝐫_{i,i+1}}{r_{i1,i}r_{i,i+1}}\right)$$
(3)   
is used. This model system was already widely studied both for flexible  and for semiflexible or liquid crystalline polymer systems .

To first approximation, there is $`\frac{l_p}{l_b}=\frac{x}{k_BT}`$ where $`l_p`$ is the persistence length (see section 3)and $`l_b`$ the bond length. In our units, the numerical values for $`x`$ and $`l_p`$ therefore coincide. This potential is applied to every bead, to every 2nd bead, every 3rd bead etc. The latter is a useful model for polymers with alternating stiffness such as single-bond, double-bond sequences or for copolymers with different persistence lengths of the constituents. In the following, we refer to a system with angular potential strength $`x`$ and a (topological) distance of $`y`$ monomers between two applications of the bond potential as $`x`$-$`y`$ system. In this sense, a fully flexible chain is referred as 0-1 chain. A 5-2 chain for example has $`x=5`$ applied to every second bond angle.

All simulated systems contained 500 chains of 50, 100 or 200 monomers each, so the overall number of particles was between 25.000 and 100.000 in a cubic periodic box.

The short-chain systems (50 monomers) could be observed until the auto-correlation function of the end-to-end vector $`𝐑_s`$ was decayed. Figure 1a shows the reorientation in the case of the 5-1 system which has the longest relaxation time. The second Legendre polynomial $`P_2(z)=(3z^21)/2`$ is used for consistency with analyses further below.

At this time, approximately, the mean square displacement of a single monomer begins to coincide with that of the center of mass. For the longer chains, we first waited until the radius of gyration and the end-to-end distance did not change systematically any more but fluctuated only around their mean values. The loss of local orientation of shorter chain segments is shown in figure 1b. One sees that there are two regimes. On short time scales ($`t^{}<5000`$), there is a fast decay due to local processes. On long time scales, however, there is a long tail which is determined by the overall motion of the whole chain. For the investigations in the following this overall motion is not important.

All systems were simulated at least for $`t^{}=20000`$. We trust that the static chain properties were well equilibrated, because the overall properties like gyration radius settled and at least local orientation decayed. Moreover, the error estimation for $`𝐑_s`$ and $`𝐑_g`$ was performed according to a binning analysis . The correlation times for the observed properties resulting from this analysis were also exceeded substantially. In the case of 5-1 with 200 monomers, which has the longest equilibration times, this “binning time” is about $`t_b^{}=8000`$. Such systems were then simulated for $`t^{}=40000`$ to 80000.

## 3 Chain Structure in the Melt

In this section, we investigate the effect of the melt environment on single chains. The presence of the other chains screens out the excluded volume interaction and the chain statistics of a self-avoiding walk appropriate for chains in good solvent turns into a simple random walk . This is, for example, evident in the single chain structure functions (see figure 8 in section 5). In the semiflexible case, one expects that, at large scales, the Gaussian statistics (random walk) is fulfilled, whereas on short scales the local stiffness is relevant. Two concepts can be used for analysis: One is the idea of a Kuhn length $`l_K`$ which is defined via

$$l_K=\frac{R_s^2}{l_b(N1)}.$$
(4)   
This assumes that the melt consists of “blobs” of length $`l_K`$ which contain inside all the local information which is not relevant on the long scales.

The second idea is the persistence length $`l_p`$ which derives from the worm-like chain model . It corresponds to the decay length of the correlation of bond orientations (the tangent vector) along the chain.

$$\mathrm{cos}\alpha (s)=𝐮(s)𝐮(0),s:\text{ monomer index}$$
(5)   
which can be shown to decay exponentially in this model

$$\mathrm{cos}\alpha (s)=e^{sl_b/l_p}.$$
(6)   
Since we do not only apply the bond angle potential to every bond, but also investigated systems with alternating stiff and flexible bonds (e.g. 5-2), the persistence lengths of these systems are not a priori known (at least the $`x`$-$`y`$ case $`y1`$). In order to determine $`l_p`$ the bond correlation function (eq. 5) was determined, in 100 configurations after the equilibration and the initial decay was fitted with an exponential $`e^{l/l_p}`$ (see table 1). If the bending potential was applied to every monomer the decay was well approximated by an exponential and the decay length $`l_p`$ was not too far from the expected value $`x`$ from the bond angle potential. This is in agreement with Monte Carlo results for stronger stiffness . In the case of alternating stiffness, minor deviations from exponential decay were observed (see figure 2a). The error in the bond correlation is about 0.03. Hence, the systems with very short persistence lengths, were difficult to determine because only very few points for fitting the decay were available and, therefore, the resulting error bars are not negligible. However, in all cases a fit over more than one order of magnitude was possible.

It is not clear if for the $`x`$-$`2`$ case the bond correlation has to follow an exponential law. However, we found this always to be the case. From a simple argument, the effective persistence length in the case of persistence lengths $`l_{p1}`$ and $`l_{p2}`$ for alternating angles is

$$\frac{1}{l_p}=\frac{1}{2}\left(\frac{1}{l_{p1}}+\frac{1}{l_{p2}}\right),$$
(7)   
because

$$e^{2l/l_p}=e^{l/l_{p1}}e^{l/l_{p2}}.$$
(8)   
This is exactly true for all points with even monomer distances in the bond correlation function. This result may be generalized to a repetitive sequence of $`n`$ different bond angle potentials

$$\frac{1}{l_p}=\frac{1}{n}\underset{j=1}{\overset{n}{}}\frac{1}{l_{pj}}.$$
(9)   
A more elaborate calculation in the framework of a generalized wormlike chain model with varying stiffness yields the same result.

The persistence length of the fully flexible model is found to be exactly one monomer distance which is on average $`l_b=0.97`$. The bond correlation functions (of inner monomers) show in the very beginning a decay with a persistence length which is close to the expected value (on the length scale of about 5 monomers). This suggests that the very local orientation correlation is determined by the “true” potential strength whereas on longer scales finite size or many chain effects contribute considerably. This is especially reflected in the persistence length values for the 5-1 chains (see figure 2b), where the bond correlation function for the shorter chain shows at distances $`s>5`$ substantial differences to the longer chain. They may be attributed to finite chain length effects. All bond correlation functions were determined starting from the innermost monomers in order to avoid end effects as much as possible.

Also the end-to-end distances and radii of gyration of the corresponding chains were calculated. They are also presented in table 1. Upon increasing $`x`$, the $`x`$-1 systems stretch the chains considerably (see figure 3). A much larger bending force constant is needed for the $`x`$-2 systems than for $`x`$-1 if one wants the same $`R_g`$. Therefore the systems $`x`$-2 behave more like fully flexible chains with a bigger monomer which is most strongly seen in the persistence length. Note that even in the 100-2 case where an almost rigid and a fully flexible bond alternate, the chain stretching is not as strong as in the 2-1 case. So there is a fundamental difference between these two scenarios.

At least in the $`x`$-1 cases we find $`2l_pl_K`$ as expected from the wormlike chain model . The relation $`R_s^26R_g^2`$ for the Gaussian chain is well fulfilled in most of our cases. The larger deviations, e.g. in the 5-1 case with 50 monomers, may be attributed to finite chain length effects. Therefore, we do not see substantial deviations from Gaussian behavior.

## 4 Melt structure

Local orientation of neighboring chains may be measured by the spatial orientation correlation function OCF. To this end, we define unit vectors between adjacent monomers $`𝐮=\frac{𝐫_i𝐫_{i1}}{|𝐫_i𝐫_{i1}|}`$. The scalar products between two such unit vectors describe the angle between chain tangent vectors

$$\mathrm{cos}\alpha (r)=𝐮_{chain1}𝐮_{chain2}.$$
(10)   
The distance $`r`$ denotes the distance between the centers of mass of the respective chain segments. In order to compare better to NMR experiments, as well as to avoid the distinction between head and tail of the chain, we use the second Legendre polynomial $`P_2(r)=\frac{1}{2}(3\mathrm{cos}^2\alpha (r)1)`$.

Figure 4a shows inter-chain orientation correlation functions of different systems. The first minimum $`(r<1)`$ is close to $`P_2=\frac{1}{2}`$ which would indicate a perfect perpendicular ordering. Two chains which come so close can only pack perpendicular because of the excluded volume interaction. The radial distribution function (RDF, see below) shows that there are very few such contacts. The first peak $`(r1.2)`$ shows a preferred parallel alignment at the distance of the first neighbor. A second parallel peak follows at $`r2`$. The intervening minima $`(r1.6)`$ get weaker for stronger orientation which indicates a stronger local parallel ordering. The OCF decays to zero with $`r`$ because the system is globally isotropic, not nematic.

The local ordering is only slightly different for the systems 0-1, 2-1, and 5-2, whereas the 5-1 chain shows a more pronounced local parallel orientation. For the 5-1 chains, there is residual parallel ordering even at the intermediate minimum $`(r1.6)`$, where the other systems show some perpendicular ordering. Except for the very few direct contacts, there is parallel orientation between neighboring chains. This ordering is visible up to about three monomer diameters. The more flexible systems (0-1, 2-1, 5-2) show qualitatively a similar ordering, but it is less pronounced and there is a intermediate preferred perpendicular orientation at the distance of about the first minimum in the radial distribution function $`(r1.6)`$.

The orientation depends only weakly on the chain length (figure 4b). This demonstrates that the effect is strictly local. This holds even though the global dynamics of the chains of different lengths is strongly different: The chains of length 50 are not yet entangled whereas the longer chains are already influenced strongly by entanglements (entanglement length in the 0-1 case $``$ 60 monomers ).

Also orientation correlation functions of longer chain segments are investigated. In this case, not only vectors connecting nearest neighbors but vectors connecting next-to-nearest neighbor beads or beads farther apart are taken into account (see figure 5).

$$𝐮_d:=\frac{𝐫_i𝐫_{id}}{|𝐫_i𝐫_{id}|}$$
(11)   
It is clear (figure 6a) that the effect of local parallel chain orientation is not restricted to segments of 2 monomers only. It persists when larger chain fragments are analyzed. On the other hand, the degree of ordering decreases with the segment size considered. Figure 6b shows again the more pronounced local ordering in the 5-1 case compared to the more flexible chains. The 2-1 and the 13-2 systems coincide. Their persistence lengths are quite similar, and for the bigger segment sizes the exact local realization of this persistence length seems to average out.

Our results are also in qualitative agreement with an early lattice Monte Carlo investigation of shorter chains . Lattice models, however, are biased in favor of orientation correlation.

The inter-chain radial distribution function $`g(r)`$ (RDF) on large scales does not change much in the case of added stiffness (figure 7a). However, there are some differences on very local scales. Both the second and third neighbor peaks are farther apart for stronger stiffness. Furthermore, the minimum between the first and second neighbor shell is not as pronounced as in the more flexible cases. The local stretching allows a closer approach of chains. This leads to a reduction of the expected correlation hole. In fully flexible systems the number of neighbors of one monomer being on the same chain increase with increasing chain length. No effect of chain length is seen here (figure 7b) which again reflects the strict locality of the structure formation in the melt.

## 5 Structure functions

The structure of single chains and of the overall melt may be additionally characterized by static structure functions. Figure 8 shows the single-chain and melt structure functions of our systems. The (isotropically averaged) melt structure function is defined as

$$S_{melt}(k)=\frac{1}{N}|\underset{m=1}{\overset{N_C}{}}\underset{j=1}{\overset{n_b}{}}\mathrm{exp}(ikr_j^m)|^2=S_{SC}(k)S_{inter}(k)$$
(12)   
where $`S_{SC}`$ denotes the single chain structure function

$$S_{SC}(k)=\frac{1}{N}\underset{m=1}{\overset{N_C}{}}|\underset{j=1}{\overset{n_b}{}}\mathrm{exp}(ikr_j^m)|^2.$$
(13)   
The first sums run over all chains ($`N_C`$: number of chains, $`m`$: chain index), the second along the chains ($`n_b=N/N_C`$: number of beads along the chain, $`j`$: monomer index along the chain).

In the limit for $`k0`$ the chain structure is no more visible but we just see a massive object which is related to the first plateau in $`S_{SC}`$. The next “scaling” regime is connected to the fractal nature of the chains. The self similarity yields a decay with $`k^{\frac{1}{\nu }}`$, where $`d=\frac{1}{\nu }`$ is the fractal dimension of the chain. For a Gaussian chain we have $`\nu =\frac{1}{2}`$. Stretching leads to a smaller fractal dimension resulting in a less steep decay. In the large-$`k`$ range we deal with structure of the size of one monomer. A bead spring model has no structure on a shorter scale. The melt structure function is the Fourier transform of the density-density correlation function. It shows therefore some additional peaks which correspond to peaks in the RDF. Hence, it contains not only information about the single-chain structure but also about the overall structure of the whole system.

The single-chain structure functions $`S_{SC}`$ of the $`x`$-1 systems look very similar, whereas there are strong differences between the 3-1 and the 3-2 system, the latter behaving very much like a fully flexible system. The crossover to the scaling regime is shifted to smaller $`k`$-vectors for stiffer chains. There is additionally a crossover ($`k=0.8`$) between two regimes in the decay which means that there is different fractal chain structure on different length scales. At larger scales (small-$`k`$ regime) the slope does not differ much from the fully flexible case, whereas on intermediate scales larger deviations occur, which indicate local chain stretching. The 3-2 system, however, is close to the fully flexible (Gaussian) system. But its slope is not as steep, which hints at a slight stretching of the chains compared to the Gaussian chain. This minor difference between the alternating stiffness and the fully flexible chains supports our earlier suspicion that the alternating chains behave like renormalized bead-spring chains with larger effective monomers. On the other hand, the stiffened chains with true semiflexibility are strongly different on intermediate scales.

The melt structure functions $`S_{melt}(k)`$ differ also between the $`x`$-1 case and the $`x`$-2 case. The latter is again very similar to the 0-1 system. The slopes in the regimes around $`k1`$ are clearly different, whereas the fine structure revealing peaks connected to the neighboring shells is quite similar. The exact positions of these peaks differ, however, which shows again that the distance to the nearest neighbors is slightly altered with stiffness. Figure 8b shows that all differences are on local scales $`(k>0.1)`$. The static structure functions coincide for $`k0`$, where the overall structure on large scales gets important. So the melt structure function shows that all systems behave similar on large scales but the systems with homogeneous and alternating stiffness differ on local scales.

## 6 Conclusions

The static structure of semiflexible polymers in the melt was determined. The stiffness strongly affects the persistence length, end-to-end distance and radius of gyration. Stiff chains are more stretched than chains only interacting via excluded volume. This also affects the local mutual ordering of the chains. Stiff chains pack more parallel on local scales whereas the overall structure remains isotropic. The chain length does not influence this strictly local phenomenon. Systems consisting of alternating stiff and flexible links behave similar to systems with much weaker overall persistence lengths. Their structure is similar to the structure of fully flexible chains with larger monomers. Their overall persistence length is smaller than expected by a analytical calculation using the persistence lengths of the respective potentials. Finally, the overall local structure of the melt differs considerably for alternating and homogeneous stiffness.

The static data presented here shows already ordering effects but its effect on dynamical properties relevant to NMR experiments can not be inferred. In order to compare directly to NMR experiments on real polymer melts, dynamical investigations are needed. Such simulations and analyses are presently being performed with the mesoscopic model of this article. Moreover, detailed atomistic simulations are underway for melts of specific polymers. They provide directly the time evolution of the atom-atom vectors monitored in the experiments, at least for short times. These simulations will be mapped onto more coarse grained simulations like the one presented here. This mapping will allow us to make the connection between the $`(x,y)`$ parameters of our model and real polymers. The presented static properties are a first step on the way of understanding ordering phenomena as examined in NMR experiments.

## Acknowledgements

Valuable discussions with Andreas Heuer, Mathias Pütz, and Heiko Schmitz are gratefully acknowledged.