Analysis of the Miscibility of Polymer Blends Through Molecular Dynamics Simulations

1. Introduction

Molecular Dynamics Simulations are important tools in the prediction of the properties of polymer materials and have therefore become invaluable aids inthe design of new materialssuited to particular applications.The combination of Molecular Simulations and the Flory-Huggins theory allows the study of the compatibility of polymer blends. MD Simulations constitute the simplest type of simulation for complex systems containing different types of interactions between bonded and non bonded atoms. These interactions are summarized in the forcefield files. The systems may be very sensitive to some details of the model; particularly the choice of the forcefield. Among different properties, the interaction parameters can be calculated from the numerical trajectories of polymer blends, which constitute statistical samples.

Polylactides (PLAs) have been studied extensively for a number of applications due to their potential utility in a number of growing technologies. PLAs are biocompatible and biodegradable aliphatic polyesters. They have been used in a variety of applications such as biomedical materials for tissue regeneration, orthopedics, drug delivery matrices, sutures and scaffolds (Gupta A.et al., 2007). In addition, PLAs synthesized from renewable resources are not harmful for the environment so they are desirable materials for packaging applications to reduce the impact of plastic packaging residues on the environment (Tsuji H. et al., 2005). PLAs include two optically active and crystallizable isomers, namely poly(L-lactide) (PLLA) and poly(D-lactide) (PDLA), and the optically inactive and amorphous mixture of isomers termed poly(DL-lactide), (PDLA). Semicrystalline PLA is a brittle material with good mechanical properties such as high tensile strength and modulus. The compatibility of PLAs in polymer mixtures is of great scientific and technological interest (Uras R., et al. 2010). However, only a few miscible counterparts have been reported for PLAs: poly(methylmethacrylate) (PMMA) (Eguiburu J.L. et al, 1998; Zhang G. et al., 2003) poly(methylacrylate) (PMA) (Eguiburu J.L. et al, 1998), poly(vinyl acetate) (PVAc), (Gajria A. M. et al., 1997; Yoon J.S. et al., 1990; Ogata N. et al, 1997, Park J.W., 2003) and poly(vinyl phenol) (PVPh) (Meaurio E. et al, 2005a, 2005b). In addition, the investigation of the phase behaviour of PLA/PVPh blends was controversial, since partial miscibility was initially reported (Zhang et al, 1998a, 1998b),but complete miscibility was later demonstrated for the PLA/PVPh blends(Meaurio E. et al, 2005a, 2005b). In addition, PLAs are also miscible with styrene-vinyl phenol copolymers (STVPh) (Zuza et al, 2008) within a wide range of compositions. In the PLA/STVPh blends, immiscibility was observed with the styrene rich copolymers and with pure Poly(styrene) (PS), (Zuza et al, 2008). The latter result was indeed expected considering the difference of solubility parameters (δ = 9.5 and 10.1 (cal/cm³)^1/2 for PS and PLA respectively; Coleman M.M. et al, 1990), and considering also the lack of specific interactions in the PLA/PS system.

The present workreportsthe MD Simulation results obtained with the Discover and Amorphous Cell modules of the Materials Studio software suite for blends of Poly(DL-Lactide) (PDLLA) withstyrene-co-vinyl phenolcopolymers (STVPh). The introduction of the energetic results obtained using Molecular Dynamics Simulation techniques in a modified Flory-Hugging Theory makes possible foretelling the phase behaviour of PDLLA/STVPh blends of different composition.

The immiscibility of the PDLLA/PS blends was confirmedby our research group in a recent work. At the same time PDLLA is completely miscible with PVPh due to the presence of attractive hydrogen bonding interactions between the hydroxyl groups (-OH) of PVPh and the carbonyls groups (-C=O) of PDLLA (Meaurio E. et al, 2005a, 2005b). The present workdeals with the effect of the insertion of hydroxyl groups in the styrene repeat units (rendering vinyl phenol (VPh) repeat units) of Poly(styrene) (PS) on the miscibility of the resulting styrene-vinyl phenol copolymers (STVPh) with PDLLA. The results provide a detailed description of the factors governing the miscibility of the system and allow assessing the number of –OH groups necessary in the copolymer to attain miscibility with PDLLA as a consequence of the establishment of intermolecular –OH O=C hydrogen bonds. Our research work in this system covers both the experimental and the modelling approaches, allowing the comparison of the results obtained from both methods. Certainly, MD Simulation results are in good agreement with the experimental measurements.

2. Computational details

MD Simulations were performed with the commercial software suite Materials Studio (v 4.1) supplied by Accelrys (San Diego, CA, USA). The Discover (molecular mechanics and dynamics simulation) and Amorphous Cell program modules were used to perform the computational work. Discover module is a molecular simulation program that provides a broad range of simulation methods, giving the ability to study molecular systems and a variety of material types (http://accelrys.com). It also enables to perform structural characterization and property prediction for the modelled systems. Amorphous Cell module is a suite of computational tools that allows to construct representative models of complex amorphous systems and to predict key properties. Amorphous Cell module employs the combined use of the arc algorithm developed by Theodorou and Suter and the scanning method of Meirovitch. This algorithm is suitable for longer chains and builds molecules in two separate steps, involving the creation of an initial guess structure, followed by relaxation of the structure to a state of minimum potential energy (Theodorou et al, 1985; Meirovitch, 1983).

Choosing the appropriate force field is one of the most important factors in obtaining real and reproducible results. For this study, the COMPASS force field has been selected because it enables accurate prediction of structural, conformational, vibrational, and thermo physical properties for a broad range of molecules in the isolated and condensed phases, and under a wide range of conditions of temperature and pressure. (http://www.scripps.edu/rc/softwaredocs/msi/cerius45/compass/COMPASSTOC.doc.html).

COMPASS (Condensed-phase Optimized Molecular Potentials for Atomistic Simulation Studies) is based on PCFF (Polymer Consistent Force Field) (Rigby D. et al., 1999). COMPASS is the first ab initio force field used for modelling interatomic interactions. The potential energy of a system can be expressed as a sum of valence (or bond), crossterm, and non-bond interaction energies given as:

E_TOTAL= E_valence+ E_crossterm + E_nonbond

The energy of valence interactions (E_valence), corresponds to energies associated with bond stretching (E_bond), valence angle bending (E_angle), dihedral angle torsion (E_torsion), and inversion (E_inversion) or also called out-of plane interaction energy (E_oop). Modern forcefields like COMPASS include a new term named Urey–Bradley (E_UB), that considers interactions between atom pairs in 1–3 configurations (i.e., atoms bound to a common atom):

E_valence= E_bond + E_angle + E_torsion + E_oop + E_UB

Cross terms included in E_crossterm increase the accuracy of the force-field by introducing correction factors to the valence energy to account for the interdependence existing between different valence terms. For example, the term named E_bond-bond considers stretch-stretch interactions between two adjacent bonds. Similarly, the COMPASS force-field includes stretch–bend, bend–bend, stretch–torsion, bend–torsion and bend–bend–torsion terms.

The non-bond interaction term (E_non-bond), represents the secondary interactions between non-bonded atoms. It is defined as the sum of the van der Waals energy (E_vdW), the Coulomb electrostatic energy (E_Coulomb), and the hydrogen bond energy (E_Hbond), given as:

E_non-bond = E_vdW + E_Coulomb + E_Hbond

In COMPASS, E_vdW is described by the Lennard-Jones 6-12 potential and the electrostatic energy is obtained from the partial charges of atoms in the system as estimated by the charge-equilibration method (Rappe A.K. and Goddard W.A., 1991). Electrostatic interactions were calculated by Ewald summation method, highly accurate in the calculation of long range interactions (Jawalkar S.S. et al, 2008). An accuracy of 0.0001kcal/mol with an update width of 5Ǻ was applied to evaluate electrostatic interactions (Meirovitch, 1983). Moreover, cell multi-pole method was used to calculate the non-bonded interactions due to its efficient capacity to simulate big systems.

2.1. Determination of parameters for the prediction of miscibility

The Hildebrand solubility parameter (δ), describes the attractive strength between the molecules of the material,exerting a huge influence on blend miscibility. The last few hundred picoseconds of the trajectory file were used to calculate the solubility parameter, described as the root of the cohesive energy density (CED), equation (4):

$$\delta =\sqrt{CED}$$E1

The molecular weight of the polymer is an important factor for the MD Simulation. In many cases they can not be performed using voluminous systems due to data storage space limitations of the computers. Therefore, calculating the minimum chain length representing the high molecular weight polymer is crucial.This minimum length can be determined from the solubility parameter values calculated for chains of the pure polymers with different lengths. The molecular size at which the solubility parameter becomes nearly independent from the length and reaches a nearly constant value can be assumed to represent the real polymer chain (Zhang M. at al, 2003).

The cohesive energy densities (CED) are known to vary considerably with the molecular weight of the polymer. If all intramolecular forces are eliminated, CED is defined as the increase in energy per mole of a material (Gestoso P. and Brisson J., 2001), given as:

$$CED=\left(\raisebox{1ex}{$E_{coh}$}\!\left/ \!\raisebox{-1ex}{$V$}\right.\right)$$E2

The computational results of CED allow the calculation of the energy of mixing, ΔE_mix, according to equation 6:

$$\Delta E_{mix}={\Phi }_A{\left(\frac{E{}_{coh}}{V}\right)}_A+{\Phi }_B{\left(\frac{E{}_{coh}}{V}\right)}_B-{\left(\frac{E{}_{coh}}{V}\right)}_{mix}$$E3

Where the terms in parenthesis represent the cohesive energies of the pure polymers (A and B) and the blend (mix), and $${\Phi }_A$$ and $${\Phi }_B$$ represent the volume fractions of the polymers in the blend.

Finally, the Flory-Huggins interaction parameter,χ, can be calculated from ΔE_mix according to (Case F.H. et al, 1994):

$$\chi =\left(\frac{\Delta E_{mix}}{RT{\varphi }_A{\varphi }_B}\right)V_m$$E4

where V_m is the molar volume of the repeat unit chosen as reference (PLA). R is the molar gas constant and T is the temperature of the simulation in Kelvin. In general, a positive value of the Flory-Huggins interaction parameter is considered indicative of the immiscibility of high molecular weight polymer blends, but actually the critical valueof χ obeys equation 8:

$${\left({\chi }_{AB}\right)}_{critical}=\frac{1}{2}{\left(\frac{1}{\sqrt{m_A}}+\frac{1}{\sqrt{m_B}}\right)}^2$$E5

where m_A and m_B represent the degree of polymerization of the pure polymers. If the interaction parameter of the blend is smaller than χ_critical the system is miscible in the whole composition range. If χ is slightly larger than the critical value, the blend exhibits partial miscibility: two phases coexist containing both components. For larger values of χ, the components are completely immiscible. By comparison of the values of χ calculated by atomistic simulation with the critical value (equation 8) the miscibility behaviour of the system can be predicted (Jawalkar S.S. et al., 2005).

3. Results and discussion

3.1. Construction of the models

Initial models were generated using the Flory´s rotational isomeric state theory (RIS), (Flory P.J., 1989). This method can be used as a basis for generating the conformations of unperturbed chains. The RIS approach is a powerful tool to predict conformational properties of polymer chains when statistical weight matrices are known (Blomqvist J., et al, 2001). The major advantage of this method is that it allows the quick evaluation of the miscibility while simultaneously includes the specific interactions between dissimilar molecules in a mixture (Zeng F.L. et al, 2009). Minimization of the system was carried out with 5000 steps using the conjugate gradient method (CGM). CGM method improves the line search direction by storing information from the previous iteration, utilizing the Polak-Ribiere algorithm. The minimized structures were examined to attain a suitable distribution of the chains in the blend and a correct number of contacts between the molecules. Furthermore, to avoid excessive overlaps between the chains, modified conditional probabilities are used to account for the nonbonded interactions between the atoms to be placed and the rest of the system.

Three-dimensional cubic unit cells with periodic boundary conditions were constructed using the Amorphous Cell module. Then, cubic cells were minimized and refined by molecular dynamics calculations. MD simulation runs were equilibrated for 200ps at 298K with time step of 1fs in the NVT ensemble (constant number of particles, N, volume, V, and temperature,T). This protocol was followed for all simulated structures.

3.2. Analysis of amorphous models

As indicatedbefore, the miscibility of poly(lactide)/poly(styrene) (PLA/PS) and poly(lactide)/poly(vinylphenol) (PLA/PVPh) blends was investigated by our research group in a recent paper using Molecular Modelling (Martínez de Arenaza et al., 2010). First, the solubility parameters of the pure polymers (PLA, PS and PVPh) were calculated. Polymer chains of different length were generated using the Amorphous cell module with the aim of determining the minimum representative length for the polymer chains. The minimum length at which the solubility parameter adopts a constant value can be consideredthe representative length formodelization purposes.As can be seen in Figure 1, the solubility parameter decreasesas chain length increases. In case of PDLLA, the solubility parameter becomes nearly constant above 20 repeat units (M_w=1440g/mol). However, only 10 repeat units of PS (M_w=1040g/mol) and 10 units of PVPh (M_w=1200g/mol) are necessary.

Figure 1.

Calculated solubility parameter versus number of repeat units for (Δ) PVPh, (◊) PS and (○) PDLLA.

Table 1 summarizes naming conventions, molecular weights (M_W), and solubility parameter values obtained from the MD Simulation (δ_MD) and from experimental measurements (δ_exp), corresponding to the pure polymers and to the copolymers. In the STVPh-10 and STVPh-20 copolymers 10 and 20% of ST units have been replaced by VPh units respectively. As can be observed, the calculated solubility parameters for al the polymers investigated in this paper show negative deviations, about 10-20%,relative to the experimental values. This behaviour has been reportedby other authors (Jawalkar S.S. and Aminabhavi T.M., 2006; Gestoso P., 2001; Mu D. et al., 2008).

The MD Simulations performed using polymer chains of representative length provide us the cohesive energy density (equation 5), necessary to calculate the energy of mixing (equation 6), from which the Flory-Huggins interaction parameter can be estimated (equation 7). In case of the PDLLA/PS system, calculated χ values were above the χ_critical line (see Figure 2), indicating immiscibility for these blends. On the contrary, negative interaction parameter values were calculated in the whole range of compositions for the PDLLA/PVPh system, indicating the miscibility of these blends (Figure 2). Figure 3 shows the snapshots corresponding to amorphous unit cells of 1:1 (mol:mol) compositions for PDLLA/PS and PDLLA/PVPh blends. In the snapshots carbon atoms are coloured grey, hydrogen atoms white and oxygen atoms are red coloured. The presence of -OH groups in PVPh allows the formation of hydrogen bonds with the carbonyl groups (-C=O) of the PDLLA as can be observed in Figure 4. These new strong interactions are the responsible for the compatibility of the PDLLA/PVPh blends, favouring the miscibility of the system. Numerous authors have studied the effect of hydrogen bonding in PVPh (Gestoso, 2003a and 2003b).

System	Acronyms	Repeat units	Mw (g/mol)	δMD (cal/cm^3)^0.5	δexp (cal/cm^3)^0.5
Poly(DL-lactide)	PDLLA	20	1440	8.6	10.6
Poly(styrene)	PS	10	1040	8.4	9.4
Poly(vinylphenol)	PVPh	10	1200	10.4	12.0
Poly(styrene-co-vinyl phenol) (10 mol% VPh)	STVPh-10	10	1056	8.8	-
Poly(styrene-co-vinyl phenol) (20 mol% VPh)	STVPh-20	10	1072	8.21	-

Table 1.

Calculated and experimental solubility parameters for the polymers studied in this paper.

Figure 2.

Difference between Flory-Huggins interaction parameters of PDLLA versus weight fraction of PDLLA

Figure 3.

Snapshots of amorphous unit cells of 1:1 (mol:mol) composition for (a) PDLLA/PS and (b) PDLLA/PVPh blends.

Figure 4.

Snapshot showing the formation of hydrogen bonds between the -OH groups of PDLLA and the -C=O groups of PVPh (O-H O=C).

In this work the analysis of the miscibility between PDLLA and PS or PVPh is extended to blends of Poly(DL-lactide) withPoly(styrene-co-vinyl phenol) copolymers (STVPh) in order to establish the number of vinylphenol units (VPh) necessary in the copolymer to achieve complete miscibility. Moreover, a comparison between the experimental and the modelling analyses is also provided here. The properties of the pure isolated polymers have been investigated from the modelization of single chains of the pure homopolymers (PDLLA, PS and PVPh) or copolymers (STVPh-10 and STVPh-20). Table 2 presents the experimental densities (Sarazin et al, 2003), molar volumes (V_M = M₀/ρ) and CED obtained from the MD simulations of the single-chain cells. The densities and molar volumes of the copolymers have been calculated assuming molar volume additivities. The model for the STVPh-10 copolymer consistson a single chain with10 repeat units containing 9 styrene repeat units and 1 VPh repeat unit (10% styrene units replaced). Likewise, STVPh-20 symbolizes a chain of 10 units containing 8 styrene repeat units and 2 VPh repeat units(20% styrene units replaced; see Figure 5). Cell sizes depend on the system, being 26.51 Å for PDLLA, 12.12 Å for PS, 11.69 Å for PVPh, 14.82 Å for STVPh-10 and 14.77 Å for STVPh-20.

System Number	Polymer	Density(g/cm^3)	MolarVolume (cm^3/mol)	CED; (10^-7cal/m^3)
1	PDLLA	1.247	57.7	7.388
2	PS	1.05	99.0	6.967
3	PVPh	1.25	96.0	10.89
4	STVPh-10	1.07	98.7	7.801
5	STVPh-20	1.09	98.4	5.425

Table 2.

Experimental densities and calculated CED for the pure polymer and copolymers.

Figure 5.

Snapshot of the STVPh-20 copolymer chain.

The next step in the miscibility investigation was the simulation of one 1:1 amorphous cell between PDLLA and two different STVPh copolymers: STVPh-10 and STVPh-20. The system was simulated using the specific parameters calculated previously from the individual polymers as the composition, molar volume and density. Table 3 summarizesthe MD Simulation results for the two blends proposed. PDLLA/STVPh blends behave quite different when the percentage of ST units replaced by VPh units increases from 10 to 20%. As can be seen in table 3, according to MD simulations the PDLLA/STVPh blend obtained after replacing one styrene unit in PS by one VPh unit shows a positive interaction parameter of moderate magnitude, χ = 0.60. This value is above the critical data line at 0.1457, so the PDLLA/STVPh-10 blend is immiscible. On the contrary, the existence of two VPh units in the PS chain (20% VPh units) results in a negative value of χ, indicating the miscibility of the PDLLA/STVPh-20 blends. Figure 6 displays 1:1 (mol:mol) amorphous unit cells of PDLLA/STVPh-10 and PDLLA/STVPh-20 blends.

Cubic Cell; 1:1	Density; (g/cm^3)	ΧPDLLA/STVPh	COMPATIBILITY
PDLLA/STVPh-10	1,158	0.60	IMMISCIBLE
PDLLA/STVPh-20	1.148	-0.29	MISCIBLE

Table 3.

Simulation details of two 1:1 (mol:mol) cubic unit cells of PDLLA/STVPh.

Figure 6.

Snapshots of 1:1 amorphous unit cell of (a) PDLLA/STVPh-10 blend (immiscible) and (b) PDLLA/STVPh-20 blend (miscible).

The accurate calculation of the minimum content of VPh units in the copolymer necessary to achieve miscibility requires the construction of additional models with intermediate compositions. With this aim, the polymer and copolymers chains lengths were doubled in order to increase the percentage of styrene units susceptible of substitution. Consequently, blends composed by PDLLA chains containing 40 repeat units and copolymer chains containing 20 repeat units were built to perform MD Simulations. Thus, the miscibility study was repeated with the new polymers and copolymers. As described before, the STVPh-5, STVPh-10,STVPh-15 and STVPh-20 acronyms represent the percentage of VPh units (5, 10, 15 and 20% respectively)insertedin the copolymer chain. Table 4 shows the CED obtained by MD Simulations for the new polymer chains, along with other selected properties.

Number	Systems	Repeat; units	Mw; (g/mol)	Density; (g/cm^3)	CED; (10^-7cal/m^3)
1	PDLLA	40	2880	1.247	6.766
2	STVPh-5	20	2400	1.06	6.932
3	STVPh-10	20	2720	1.07	6.970
4	STVPh-15	20	3040	1.08	7.544

Table 4.

CED and other selected properties for the double-length polymers and copolymers.

Then, the mixtures were investigated showing again completely different behaviours with the degree of VPh in the copolymer chain. Table 5 displays the results of 1:1 cubic unit cells of blends. Recalling that, as discussed before (equation 8), χ_critical depends on the degrees of polymerization of the pure polymers (m_A and m_B), the χ_critical corresponding to the double-length chains decreases to 0.0728 (m_PDLLA and m_PVPh-S are 40 and 20 respectively). Molecular Modelling results for the PDLLA/STVPh-10 blends provide a χ value of 0.22 suggesting the immiscibility of the system according to the Flory-Huggins theory. However, χ for PDLLA/STVPh-15 blends was negative (-1.18), indicating the miscibility of these blends.

Cubic Unit Cell; 1/1	Density; (g/cm^3)	ΧPDLLA/STVPh	COMPATIBILITY
PDLLA/STVPh-5	1.153	0.39	IMMISCIBLE
PDLLA/STVPh-10	1.158	0.22	IMMISCIBLE
PDLLA/STVPh-15	1.163	-1.18	MISCIBLE
PDLLA/STVPh-20	1.163	-1.33	MISCIBLE

Table 5.

MD Simulation results for 1:1 cubit unit cell

Table 5 actually proves the miscibility of the 1:1 blends but does not tell about the rest of compositions. Therefore, we have carried out MD simulations for PDLLA/STVPh-15 blends of different composition, and Figure 7 shows the snapshots of two cubic amorphous unit cells for the PDLLA/STVPh-15 blends of 2:3 and 3:2 molar compositions. The range of compositions studied for this system actually includes the whole range of compositions, and the results are summarized in table 6.

Figure 7.

Snapshot of the cubic amorphous unit cells for PDLLA/STVPh-15 miscible blends of different molar composition: a) 2:3 and b) 3:2.

System; number	PLA; molar ratio; (%)	Number of; chains in the; cubic cell; PLA/STVPh-15	Density; (g/cm^3)	ΧPDLLA/STVPh
1	0	0/1	1.07	-
2	14.88	1/5	1.099	-0.87
3	19.15	1/4	1.105	-1.06
4	24.00	1/3	1.114	-1.39
5	38.71	2/3	1.141	-1.25
6	44.83	1/1	1.158	-1.43
7	65.45	3/2	1.176	-0.73
8	65.29	2/1	1.888	-0.55
9	79.12	4/1	1.212	-1.62
10	82.56	5/1	1.217	-1.49
11	100	1/0	1.247	-

Table 6.

Simulation details for the PDLLA/STVPh-15 blends.

As can be seen, in table 6, negative interaction parameters have been obtained across the whole range of compositions. Therefore, MD simulations predict complete miscibility for the PDLLA/STVPh-15 system. Figure 8 displays the dependence on composition for the interaction parameter. Considering the typical errors associated to these calculations, χ does not show a clear dependence with composition, it rather fluctuates about the value χ = -1.0.

These results can be compared with the experimental investigation of the PDLLA/STVPh system carried out recently by our research group (Zuza E. et al., 2008). High molecular weight STVPh copolymers of different composition were blended with high molecular weight PDLLA, and the analysis of the blends by Differential Scanning Calorimetry (DSC) showed that at least 16% of VPh units were necessary in the copolymer to achieve miscibility. Polystyrene and the STVPh copolymers containing less than 16% of VPh units were immiscible with PLA-s (Zuza E. et al., 2008). As can be seen, MD simulations show good agreement with the experimental results, and are a valuable tool in the prediction of the miscibility and phase behaviour of polymer blends.

Figure 8.

Dependence of the interaction parameter with composition for the PDLLA/STVPh-15 blends according to MD simulations.

4. Conclusions

The miscibility of Poly(DL-Lactide) (PDLLA) withstyrene-co-vinyl phenolcopolymers (STVPh) has been investigated using MD simulations and the results have been compared with the experimental information available. The MD simulations indicate the formation of hydrogen bonds between the –OH groups in the VPh repeat units and the C=O groups in PDLLA, in agreement with the results reported recently for these systems using FTIR spectroscopy (Meaurio E. et al, 2005a, 2005b; Zuza E. et al., 2008). According to the MD simulations, 10-15% VPh units must be introduced in the PS chain to achieve miscibility with PDLLA. This result is in very good agreement with the experimental results obtained in high molecular weight polymers, in which miscibility was only achieved for blends of PDLLA with STVPh copolymers containing at least 16 mol% of VPh repeat units. The good agreement between the modeling and the experimental results indicates that MD simulations are a valuable tool in the prediction of the miscibility and phase behaviour of polymer blends; particularly in cases where the polymers under consideration are not available.