Molecular Modeling Simulations to Predict Compatibility of Poly(vinyl

Molecular Modeling Division, Center of Excellence in Polymer Science, ... 128 ACE Extension, South Kensington Campus, Imperial College, United Kingdom...
0 downloads 0 Views 554KB Size
J. Phys. Chem. B 2007, 111, 2431-2439

2431

ARTICLES Molecular Modeling Simulations to Predict Compatibility of Poly(vinyl alcohol) and Chitosan Blends: A Comparison with Experiments Sheetal S. Jawalkar,† Kothapalli V.S.N. Raju,‡ Shivaraj B Halligudi,§ Malladi Sairam,⊥ and Tejraj M. Aminabhavi*,† Molecular Modeling DiVision, Center of Excellence in Polymer Science, Karnatak UniVersity, Dharwad, India 580 003, Organic Coatings and Polymers DiVision, Indian Institute of Chemical Technology, Hyderabad, India 500 007 Catalysis DiVision, National Chemical Laboratory, Pune, India 411 008, and Department of Chemical Engineering and Chemical Technology, 128 ACE Extension, South Kensington Campus, Imperial College, United Kingdom ReceiVed: October 18, 2006; In Final Form: January 17, 2007

Molecular modeling simulations are the most important tools to predict blend compatibility of polymers that are otherwise difficult to predict by experimental means. Conflicting reports have been reported on the blend compatibility of poly(vinyl alcohol), PVA, and chitosan, CS polymers. Since both the polymers are widely used in pharmaceutics as drug-loaded particulates and as separation membranes, we felt it necessary to investigate their compatibility over the practical range of compositions. In this paper, we attempt to study the compatibility of PVA and CS polymers using molecular modeling strategies to understand the interactions between CS and PVA polymers to predict their compatibility from atomistic simulations. Flory-Huggins interaction parameter, χ, was computed at 298 K to assess the blend compatibility at different ratios of the component polymers. Miscibility was observed for blends below 50% of PVA, while immiscibility was prevalent at compositions between 50 and 90% PVA. Computed results confirmed the experimental findings of dynamic mechanical thermal analysis, suggesting the validity of modeling strategies employed. Plots of Hildebrand solubility parameter and cohesive energy density calculated at 298 K supported these findings. The χ values for blends, which satisfied the criteria of miscibility of polymers computed by atomistic simulations, agreed with the solubility criteria related to order parameters calculated from mesoscopic simulations. Miscibility between PVA and CS polymers is attributed to hydrogen bond formation and to an understanding of which of the interacting groups of CS, i.e., -CH2OH or -NH2, are responsible in blend miscibility. This was further confirmed by molecular dynamics simulations of radial distribution functions for groups or atoms that are tentatively involved in interactions. These results are correlated well to obtain more realistic information about interactions involved as a function of blend composition. Computed freeenergy from the mesoscopic simulation for blends reached equilibrium, particularly when the simulation was performed at higher time step, indicating stability of the blend system at certain compositions.

Introduction Chitosan (CS), a naturally occurring polysaccharide, has been extensively studied in drug delivery applications,1-4 as a membrane in pervaporation separation5-7 and a solid polyelectrolyte in fuel cells.8 In the polymer literature, considerable efforts have been made to modify the basic structure of CS by blending9,10 to improve its properties. Poly(vinyl alcohol) (PVA) is another widely studied polymer that finds applications in drug delivery11 as well as a membrane in separations.12 Blends of CS and PVA have been investigated using a variety of experimental techniques.13-16 The experimental protocols, though * Corresponding author fax: 91-836-2771275; e-mail: aminabhavi@ yahoo.com. † Karnatak University. ‡ Indian Institute of Chemical Technology. § National Chemical Laboratory. ⊥ Imperial College.

accurate, are time-consuming and expensive; often the results are contradictory. Among the theoretical methods used to study polymer-polymer blend compatibility, mesoscopic dynamics17,18 (MesoDyn) and dissipative particle dynamics (DPD) methods19,20 have been used to treat polymeric chains in a coarse-grained (mesoscopic) level by grouping atoms together up to the persistence length of the polymers. MesoDyn is a dynamic mean field density functional theory in which the dynamics of phase separation is described by Langevin-type equations to study polymer diffusion. This method was successfully employed for concentrated aqueous solutions of triblock copolymers21-25 and employs soft interaction potentials allowing large time-scale simulations. Recent trends in the use of molecular mechanics (MM) and molecular dynamics (MD) simulations on bulk polymers have led to the calculations of important properties of polymeric blends with greater accuracy26 that are otherwise difficult to

10.1021/jp0668495 CCC: $37.00 © 2007 American Chemical Society Published on Web 02/21/2007

2432 J. Phys. Chem. B, Vol. 111, No. 10, 2007 SCHEME 1: Poly(vinyl Alcohol) and Chitosan Polymers

study experimentally or for which experimental data are not available or often when such results are conflicting. Foremost among the widely used systems are those of PVA and CS polymers (see Scheme 1), whose blends have been used in many areas,13,14 but the recent publications of Chang et al.27 and Yang et al.28 suggested that these are not fully compatible at the molecular level. However, several papers published in the literature27-30 point to the conclusion that CS and PVA blends are compatible, depending upon the method of fabrication and/ or at certain blend compositions. This has prompted us to undertake a detailed investigation on the study of blend compatibility of PVA and CS polymers over a practically useful range of compositions using computer-assisted MD simulations. Experimental study on the blend compatibility of CS with PVA indicated that PVA and CS are immiscible at higher compositions of PVA (i.e., >50%), while miscibility occurs at lower compositions of PVA. However, blending of PVA and CS poses difficulties due to the tight intramolecular structure as a result of hydrogen-bonding interactions. Experimental investigations of dynamic mechanical thermal analysis (DMTA) suggested that these polymers are not completely miscible over the entire range of compositions, but they start to mix at and beyond 50:50 wt. % ratio of PVA/CS. Thus, it is important to investigate miscibility/immiscibility characteristics of these polymers using the MD simulation protocols. It is well-known that the chemical structure of side chain of vinyl polymers plays an important role in their miscibility behaviors. CS has four hydroxyl groups, which are partially hydrolyzed and, hence, few carboxylic acid groups are present. Thus, CS can be considered as a strongly interacting polymer, which helps to form compatible blends.30 Therefore, it is of interest to understand which of the interacting groups of CS are responsible for blend compatibility based on the calculation of radial distribution function (RDF). Evidence of partial miscibility of blends of these polymers in solid film forms has been studied experimentally by measuring their Tg by the DMTA technique. Strong evidence of immiscibility of blends of PVA and CS at higher compositions of PVA has been observed, indicating that blends are not completely miscible in all proportions at all compositions. Even though the reason for this behavior has not been fully explored, our own DMTA data and the present simulation strategies suggest that the presence of bulky NH2 group in CS (MW ) 800 000; Tg ) 42.35 °C; density ) 0.670 g/cm3) restricts the

Jawalkar et al. free rotation around -C-C- bond of its backbone, thereby hindering the chain mobility, resulting in an increase of Tg. On the other hand, hydroxyl moiety of PVA (MW ) 125 000; Tg ) 64.1 °C; density ) 1.269 g/cm3) is small compared to NH2 group in CS and hence, it provides an ability to crystallize PVA much faster than CS. Sandoval et al.29 as well as Castro et al.30 studied specific interactions in blends containing CS and functionalized polymers through MD simulations. In the present study, MD simulations of the oligomers of PVA and CS have been performed at the ambient temperature over wide range of blend compositions. Cohesive energy density, CED of the pristine polymers in the blends was obtained to compute solubility parameter, δ of the blends as a function of blend composition. Glass transition temperatures of the individual polymers and their blend films were determined by DMTA experiments as a function of compositions of CS. FloryHuggins interaction parameter, χ, was calculated to understand the energetics in the mixing of polymers, which indicated favorable interactions. Interaction parameters obtained from atomistic simulations along with other structure-dependent (monomer number and length, characteristic ratio, etc.,) parameters were subsequently supplied in the mesoscopic simulations. Short-scale characteristics were then absorbed into large structureless beads. Coarse-grained representation of the systems was used to study the blends of high molecular weight and for extension of time scales that allowed the observation of phase separation. The degree of order (order parameter) of the phases formed was investigated. These results indicated a close agreement with the reported literature data27,28 as well as with our own experimental results, appropriately predicting the miscibility/immiscibility trends of the blends as a function of blend compositions.31 Experimental Section Blend Preparations. Chitosan solutions were prepared by dissolving the required amount of CS with 85% deacetylation in 2% acetic acid solution at ambient temperature with constant stirring overnight, and the solution was filtered before use. PVA solutions were prepared by dissolving the required amount of PVA in distilled water heated at 80 °C under constant stirring for 4 h. Then, required amounts of CS and PVA solutions were mixed, the mixture was stirred well for 24 h, films were cast in a petri dish at ambient temperature and kept for 48 h. Using the above method, PVA-CS films containing 10-90 wt. % of CS were prepared. These are designated, respectively as PVACS-10 to PVA-CS-90. Individual films were prepared by casting the respective polymer solutions in petri dishes and subsequently, drying them at ambient temperature for 24 h. Dynamic Mechanical Thermal Analysis (DMTA). Tan δ curves of PVA, CS, and their blends were recorded on a Rheometric Scientific DMTA IV instrument in the tensile mode operating at 1 Hz frequency. Samples were scanned from -50° to 200 °C at the heating rate of 3 °C/min. These measurements were done at the Indian Institute of Chemical Technology, Hyderabad. Simulation Strategies. MD simulations have been performed on the blends of PVA and CS at ambient temperature (298 K) using the software packages of MesoDyn and MesoDyn Interfaces (Accelrys, Inc., San Diego, CA) with the Material Studio Modeling (version 3.2) installed on Windows-2000. Using the rotational isomeric state (RIS) theory that describes the conformations of unperturbed chains, initially the structures were built. For the good blending of polymers, amorphous phases were checked for space filling regularly after the

Compatibility of Poly(vinyl Alcohol) and Chitosan Blends TABLE 1: Solubility Parameters of PVA and CS Polymers solubility parameter (cal/cm3)0.5 polymer repeat units molecular weight COMPASS expt (ref 46) PVA CS

80 20

3520 3580

11.9 10.2

12.6 10.3

construction of the initial amorphous cell. If the two component chains were not well mixed (sufficient intermolecular contacts) in the initial configuration, then it was discarded and a new one was attempted. This was attempted for more than 12 independently generated structures for each blend; out of these, one structure was selected to estimate solubility parameter, δ. The repetitions of these calculations were accurate to (0.5 units for PVA and (1.8 units for CS. Bulk phases were first generated using the amorphous cell after minimizing the initial structures using the conjugate gradient method (CGM), which utilizes the Polak-Ribiere algorithm with a convergence level of 0.001 kcals/mol/Å. The COMPASS (condensed-phase optimized molecular potentials for atomistic simulation studies) force field,32 which is based on the PCFF (polymer consistent force field) model was used for computing interatomic interactions. This force field has been widely used to optimize and predict the structural, conformational, and thermophysical condensed phase properties of the molecules including polymers. The combination of valence terms, including diagonal and cross-coupling terms for bond stretching, bond angle bending, and dihedral angle distortions, Columbic term for electrostatic interactions and Lennard-Jones 6-12 potential for van der Waals interactions, comprised the energy expression. The analysis of results was carried out after construction of the amorphous cell and dynamics of 300 ps was simulated. The cell multipole method33-35 was used to calculate nonbonded interactions. This method has been proved to be efficient to simulate big systems, since it scales linearly with the number of atoms, N (compared to cutoffs that scale as N2), which requires modest hardware memory. The periodic box is divided into M cubic cells with M = N/4. In each cell, the cells in the nearest neighborhood contribute to the near-field potential, while others to the far-field potential (short and long-range interactions). Interactions between atoms in the near-field cells are calculated directly for each pair of atoms, but interactions for atoms in the far-field cells are computed via expansions of multipole moments (charges, dipoles, quadruples and octapoles) around the center of each

J. Phys. Chem. B, Vol. 111, No. 10, 2007 2433 cell. Since it is assumed that interactions of atoms in the central cell with far-field atoms vary little from one atom position to another (compared to interactions between the near-field atoms), Taylor coefficients of the expansions are calculated only for each cell and not for each atom. Calculated solubility parameter values for individual polymers are given in Table 1. In the COMPASS force field approach, total energy, ET of the system is represented by the sum of bonding and nonbonding interactions given as the following:

ET ) Eb + E0 + Eφ + Eoop + Epe + Evdw + Eq

(1)

Here, the first five terms represent the bonded interactions, which correspond to energies associated with the bond, Eb; bond angle bending, E0; torsion angle rotations, Eφ; out of loop, Eoop; and potential energy, Epe. The last two terms represent nonbonded interactions, which consist of van der Waals term, Evdw, and electrostatic force, Eq. In COMPASS, Evdw is invariably described by the Lennard-Jones 6-12 potential, while the electrostatic energy is calculated from the partial charges of atoms in the system as estimated by the charge-equilibration method.36 Electrostatic interaction was calculated by the Ewald summation method,37 since it calculates long-range interactions more accurately. The potential and nonbonded energy fluctuations observed during molecular dynamics of 300 ps are displayed in Figure 1. To choose the polymer chain length for MD simulation of both CS and PVA for blending, individual polymers are first subjected to minimization and then dynamics was performed to find the solubility parameter, which is plotted vs the number of repeating units of the polymer. When a stable value is obtained, it confirms that the number of repeating units are sufficient for simulation. In the present work, PVA consisting of 80 repeating units and CS of 20 units have been used to construct the amorphous cell of 1:1 blend. Density of the chosen polymers were taken from the literature38 as F(PVA) ) 1.269 g/cm3 and F(CS) ) 0.670 g/cm3. The system investigated consists of 1024 atoms for 1:1 blend combination. Further, properties of interest such as CED, δ and χ have been calculated. MD simulations under constant temperature and density (NVT ensemble) conditions have been performed for each configuration using the Discover program of Accelrys. It uses the Anderson thermostat for controlling the temperature fluctuation. For dynamics, Andersen’s stochastic method was used to control the temperature. A group-based cutoff of 8.5 Å with a spline

Figure 1. Potential and nonbonded energy fluctuations observed during molecular dynamics of 300 ps.

2434 J. Phys. Chem. B, Vol. 111, No. 10, 2007

Jawalkar et al.

χ)

z.∆Emix RT

(2)

where z is coordination number, whose value for cubic lattice model was taken to be six, R is molar gas constant (cal/mol), and T is temperature in Kelvin, at which the simulation was performed. Energy of mixing, ∆Emix needed to compute χ was calculated as

∆Emix ) φA

Figure 2. Snapshot of unit cell for 1:1 (PVA:CS) blend. The display style of the blend is shown in CPK model (Corey-Pauling-Koltun models) where radius is dependent on van der Waals radius of the element. Carbon atom, gray; hydrogen atom, white; oxygen, red.

TABLE 2: Simulation Details for PVA/CS Blendsa system no.

number of chains per unit cell

composition (wt. % PVA)

density (g/cm3)

χ (eq 2)

1 2 3 4 5 6 7 8 9 10 11

1CS chain 10:90 (PVA:CS) 20:80 (PVA:CS) 30:70 (PVA:CS) 40:60 (PVA:CS) 50:50 (PVA:CS) 60:40 (PVA:CS) 70:30 (PVA:CS) 80:20 (PVA:CS) 90:10 (PVA:CS) 1PVA chain

0 10 20 30 40 50 60 70 80 90 100

0.6700 0.7299 0.7898 0.8497 0.9096 0.9695 1.0294 1.0893 1.1492 1.2091 1.2690

N/A 0.1710 0.0110 0.0402 0.0650 0.5330 0.2180 0.2050 0.2600 0.2340 N/A

a

N/A, not applicable.

width of 1 Å was applied to evaluate nonbonded interactions, which owed to the separation between atoms. In the present study, for pristine polymers (PVA or CS) and their blends, bulk amorphous states were built using the cubic unit cells subjected to periodic boundary conditions. Simulations were performed with the Discover MD simulation modules. Systems were built with a 3D periodicity and equilibrated in the NVT ensemble at 298 K. This equilibration was usually done for 5 ps with the dynamics that was followed by the data accumulation run lasting at least 100 ps with configurations saved every 5 ps. The snapshot unit cells (CPK model) for 1:1 blend of PVA/CS is shown in Figure 2. Carbon atoms are gray, hydrogen white, and oxygen red in color. The detailed model construction procedure is described by different ratios of number of chains of PVA to number of chains of CS in the unit cells. The number of chains per unit cell, chain length, and density values are summarized in Table 2. Density of the blends was calculated from the density of individual polymers and volume fraction of each polymer. In general, the initial amorphous structure is in a relatively high-energy state. Therefore, before performing MD calculations, we had to perform the energy minimization after the system comes close to minimum. Simulation time depends upon the number of atoms in the system, but simulation was carried out until total energy of the system was stabilized. The last few hundred ps of the trajectory files were used to calculate δ defined as the square root of CED as well as Flory-Huggins39 interaction parameter, χ given by

( ) Ecoh V

A

( ) ( )

+ φB

Ecoh V

-

B

Ecoh V

mix

(3)

In the above equation, subscripts A, B, and mix represent CED values of PVA, CS, and PVA/CS blends, respectively, by considering the identity: CED ≡ (Ecoh/V). Symbols φA and φB represent the volume fractions of PVA and CS, respectively. Using the calculated value of energy of mixing for CS and PVA, χ was calculated from eq 2. Pair Correlation Function. The pair correlation function is a measure of probability that an atom at the origin of an arbitrary reference frame located in a spherical shell of infinitesimal thickness at a distance, r from the reference atom. This concept also embraces the idea that the atom at the origin and the atom at distance, r may be of different chemical species, A and B. The resulting function is given by the symbol, gAB(r), referred to as the radial distribution function (RDF), which was calculated by the average of static relationship of every given pair of particles, AB as follows:

〈nAB(r)〉 gAB(r) ) 2 4πr ∆FAB

(4)

where < nAB(r) > is average number of atom pairs between r and r + ∆r and FAB is density of atom pairs of type AB. The total pair distribution function, g(r), which gives a measure of spatial organization of atoms about the central atom, can be used to demonstrate long-range order of the structure. Therefore, g(r) was used to distinguish between amorphous and crystalline structures of the polymers. RDF was calculated for various pairs of atoms of molecules, since it gives an insight as to how the atoms pack in a amorphous structure. To perform RDF calculations, we have considered the interactions between oxygen atom of hydroxymethyl group and nitrogen atom of amino groups of CS with hydrogen atoms of hydroxyl groups of PVA. Results and Discussion Dynamic Mechanical Thermal Analysis (DMTA). One of the most common methods to estimate polymer-polymer compatibility is to determine Tg of the blend and compare it with the Tg of the individual polymers. This method is used here to study the compatibility of CS and PVA polymers. The tan δ versus T curves of these blends is given in Figure 3. This figure gives the overlay of tan δ data of CS and with various blends of PVA. This shows that there is a peak at 40.8 °C which corresponds to pure CS and another peak at 64.1 °C corresponds to pure PVA, whereas the Tg of blends are between these two values. There is a single and gradual increase in Tg value containing lower amounts of PVA, i.e., up to 50%. Then the broad peaks were observed afterward, i.e., 60% and above. This shows that CS and PVA polymers are miscible up to 50% and the peaks became broad for blends containing more PVA, which suggests that blends containing 20-40% PVA and 80-60% CS are miscible.

Compatibility of Poly(vinyl Alcohol) and Chitosan Blends

J. Phys. Chem. B, Vol. 111, No. 10, 2007 2435

Figure 5. Solubility parameter of CS vs repeating units of CS.

Figure 6. Cohesive energy density vs weight % PVA.

Figure 3. (a) Tan δ versus T plot of CS/PVA blends from PVA 0-50% (b) Tan δ versus T plot of CS/PVA blends from PVA 60-100%.

monomer units of CS. Almost a similar dependency was observed for the blends of poly(ether imide) and poly(carbonate) studied by Zhang et al.,41 wherein it was found that for polycarbonate, 10 repeating units were sufficient to perform the simulation. Similarly, according to Fan et al.,42 21 repeating units of polycarbonate were enough to achieve the simulation. In the present work, the tendency of PVA and CS polymers to mix at a specific composition was examined via the calculations of cohesive energy density of the blends and of pure components. It is, therefore, realized that solubility parameters of individual polymers (PVA or CS) exert an influence on their blend miscibility. Figure 6 displays the typical variation of CED with varying weight fraction of CS over the entire composition of the blend. An increase in CED curve is observed with increasing composition of PVA of the blend, indicating lesser interaction at higher compositions of PVA. As a further test of the compatibility/incompatibility behaviors of PVA and CS polymers, we have computed the critical value of χ using the following equation:39

(χAB)critical ) Figure 4. Solubility parameter of PVA vs repeating units of PVA.

Atomistic Simulation. Size of the macromolecule used in model calculations is important to calculate thermodynamic parameters of interest. Thus, it is of interest to understand what minimum molecular size is sufficient to represent the real polymer chain. To determine this minimum size, solubility parameters of CS and PVA at the carefully chosen different molecular weights have been computed up to a point even if the increase of polymer molecular weight did not change the solubility parameter. As displayed in Figure 4, solubility parameter curve of PVA varies within a narrow range and then levels off beyond 100 repeating units of PVA. Computed values of δ are plotted vs repeating units of CS (up to 100 monomer units) in Figure 5. It is observed that, beyond 20 monomer units, the δ values did not vary much and leveled off beyond 60

(x

1 2

1 mA

+

1

)

xmB

2

(5)

where mA and mB represent the degree of polymerization (actual number of repeating units) of A and B. It may be noted that even though the χ parameter has been discussed at various levels40,43 to be a dependent parameter on mixture compositions, in the present research it was realized that composition dependence of χ led to further computational problems. Hence, the calculated values of χ parameter are good enough to obtain the reproducible data.26 Notice that blends are miscible if χ is < (χAB)critical. If χ is considerably greater than the critical value, then blends are immiscible (i.e., they form two separate phases). If the value of χ is slightly greater than the critical value, then blends are partially miscible in which case, both the phases are present, one for each component of the blend. So long as the measured χ and the corresponding critical values are calculated based on the same reference volume, comparing the measured χ with the critical values would provide a good indication of

2436 J. Phys. Chem. B, Vol. 111, No. 10, 2007

Jawalkar et al.

Figure 7. Flory-Huggins interaction parameter vs weight % of CS.

Figure 9. Plot of free energy density vs time step for (4) 40:60(PVA/ CS), (O) 1:1(PVA/CS), and (0) 80:20 (PVA/CS) blend systems.

Figure 8. Molar enthalpy and vs weight fraction of CS.

the degree of blend miscibility. Therefore, one can use χAB to characterize blend miscibility.40 Results of χ vs weight fraction of CS calculated from eqs 2 and 5 are displayed in Figure 7. For the blends of 10:90 (PVA/CS) to 40:60 (PVA/CS), the value of χ lies below the χcritical line. For other compositions, i.e., 1:1 and above (with increasing composition of PVA), the points lie above the critical value of χ indicating blend immiscibility. The value of χcritical is around 0.18, which represents the limiting line for the blend to be miscible. The χ value is higher than χcritical for 50:50 PVA/CS blends and above, indicating that blends are immiscible at higher compositions of PVA, whereas χ value for PVA/CS blends containing less than 40% PVA are below the χcritical, indicating the miscibility of PVA and CS polymers at that composition26 (see Figure 7). Enthalpy of mixing, ∆Hm of the blend system is another important parameter used to understand their miscible/immiscible trends; this was computed as

∆Hm ) Vm{(CEDA)0.5 - (CEDB)0.5}2φAφB

(6)

Here, Vm is total volume of the blend and CED values used in computing ∆Hm were taken from MD simulations. The plot of ∆Hmix vs weight fraction of CS is displayed in Figure 8. Notice that for the blend to be miscible, larger contributions are seen from Etorsion, Ebond angle, and Eoop terms. Once the atoms for 1:1 blend ratio are fixed, the blends of 10:90, 20:80, 30:70, 40:60, 60:40, 70:30, 80:20, and 90:10 of PVA:CS have been studied to investigate whether they are miscible or not. Here, the time required to attain equilibrium depends upon the size of the system, which was determined by running the simulation until complete minimization. Once the average total energies were calculated for the system, then ∆Hm was calculated for each blend system using eq 6. Mesoscopic Simulations. As a further proof, the MesoDyn program was also used to simulate the phase separation

Figure 10. Mesophase order parameter for blends of (4) 40:60(PVA/ CS), (O) 1:1(PVA/CS), and (0) 80:20 (PVA/CS)] as a function of time step.

dynamics of the blends at the mesoscopic level. This approach is based on the dynamic variant of mean-field density functional theory,18 which is similar to classical dynamic random phase approximation (RPA).42 Here, the polymer chains are modeled as ideal Guassian chains consisting of beads, each representing the monomer of the polymer (Kuhn statistical segments). Free energy of the system was then calculated in terms of bead distribution functions. Positions of the beads were then correlated to each other by converting it to a many-body problem. To overcome this problem, interchain correlations were neglected and the system is approximated by a set of independent Gaussian chains embedded in the mean-field. A key rudiment of the dynamic density functional theory is that on a coarsegrained time scale, the free energy function is minimized. Details of these derivations were given earlier.26 Parametrization. Polymer chain lengths were determined from the degree of polymerization and characteristic ratios of the polymers. The expression for MesoDyn chain length (Nmeso) is given by

NMeso )

Mp MmCn

(7)

where Mp is polymer molecular weight, Mm is monomer molecular weight, and Cn is characteristic ratio.44 MesoDyn input parameter is related to Flory-Huggins interaction parameter through the following relation:

ν-1ij ) χijRT

(8)

Compatibility of Poly(vinyl Alcohol) and Chitosan Blends

J. Phys. Chem. B, Vol. 111, No. 10, 2007 2437

Figure 11. Density profile slices for 40:60 (miscibile) and 80:20 (immiscible) blend of PVA/CS showing the extent of mixing; green color represents pristine PVA, red color for pristine CS and gray shaded areas correspond to interfaces between them.

Figure 12. Radial distribution function calculation for PVA/CS blend representing oxygen atom of hydroxyl methyl group of CS and nitrogen atom of CS relative to the distance of hydrogen atom of hydroxyl group of PVA.

where χij is taken from atomistic simulations carried out for each blend at each composition (as discussed before), R is molar gas constant (8.314 J/mol.K) and T is 298 K. The Nmeso and ν-1IJ (interacting energies between pair of interacting bead types) are the input parameters for MesoDyn calculations. For 50:50 (PVA/CS) blend, Nmeso taken for PVA and CS were 9.6 and 1.99 units, respectively, while ν-1IJ was 22.7 kJ/mol. For 3:1 (PVA/CS) blend, the values of Nmeso taken for PVA and CS were 0.902 and 0.48, respectively, while ν-1IJ was chosen to be 21.9 kJ/mol. For 1:3 blend, Nmeso was 0.869 for PVA and 0.488 for CS with ν-1IJ ) 25.84 kJ/mol. To account for numerical stability, the time step for performing the simulation was chosen in such a way that dimensionless time step, τ, used by the program was 0.5 (i.e., between 0 and 1) and bond length was 1.154 nm throughout. Thus, 200 µs time step was used for PVA/CS mixtures at 40:60 and 80:20 compositions showing miscibility and immiscibility, respectively. A constant noise parameter of 75.002 was maintained for the entire simulation (since too high or too low value will lead to system unstability42), which was applied to shorter chains with longer statistical units.18 The grid dimensions were 32 × 32 × 32 nm and size of the mesh over which density variations are to be plotted in MesoDyn length units using the Grid spacing field were 1 nm. Bond length was 1.1543 Å times the cell length to guarantee isotropy of all the grid-restricted operators. The temperature 298 K corresponds to the temperature at which χ parameter was calculated for a total simulation of 200 µs. Nonideality in the system was incorporated through the effective external potentials, which were obtained from the χ parameter

for the pair of species chosen. Hence, the chemistry of the system is determined by the values included in the external potentials. After performing the calculation of energy of mixing using the atomistic simulation, it was related to Flory-Huggins parameter as

χ)

( )

∆Emix Vmon RT

(9)

However, there are other alternative routes to χ parameter, but not at least from the experimental data, e.g., vapor pressure or interfacial tension data. In the mesoscopic simulation, after setting up the initial configurations, systems were let to evolve toward equilibrium (phase separation or mixing). More specifically, the initially homogeneous phases were cooled in the twophase region by tuning on interactions at time, t ) 0. Then, blend compatibility was deduced from the final equilibrium morphologies. Order Parameter. During the simulation step, free energy should asymptotically approach a stable value as the system attains dynamic equilibrium. However, the free-energy density is not routinely calculated for real systems, and hence, its direct comparison with experimental data is not possible. However, the evolution of free energy is a good measure of the stability of a system. Plots of free-energy density vs time step given in Figure 9 shows that the system reached equilibrium particularly when it was simulated for higher time interval, indicating its stability. The order parameter, Pi, defined as the volume average of the difference between local density squared and the overall

2438 J. Phys. Chem. B, Vol. 111, No. 10, 2007

Jawalkar et al.

density squared, is given by the following integral:15

Pi )

1 V

∫V [ηi2(r) - ηi2]dr

(10)

where ηi is dimensionless density (volume fraction) for species, i. Order parameters with large values indicate strong phase segregation, which is observed in case of 80:20 blend of PVA/ CS; conversely, very small values indicate the miscible nature of 1:3 PVA/CS blend. Notice that the plot of order parameter vs time steps in ns is constant over the time scale from 0 to 4000 ns as shown in Figure 10. It is observed that for 1:1 blend, order parameter is >0.1, indicating the immiscible nature of the blends, while for 1:3 blend of PVA/CS, it is almost close to 0.1, indicating blend miscibility (see distinct phase miscibility with well defined interfaces shown in Figure 11). The density profiles presented in Figure 11 show the phase separation for different time scale simulations with miscible and immiscible blends. Notice that phase separation proceeds via the diffusion of components through the interfaces. Hydrodynamics becomes important, especially at later stages of phase separation barrier to further diffusion.45 Results of this study are in accordance with the internal experimental findings (see DMTA curves in Figure 3). Figure 12 displays the total pair correlation function as a function of r for bulk simulated structures where one could see the variation of g(r) for a polymer backbone containing carboncarbon atoms of the simulated structures of PVA and CS. The peaks of g(r) indicate the presence of definite correlation between atoms within that radius, while the absence of any peaks beyond 4 Å distance indicates that there is no long-range interaction in the systems. However, at long distances, g(r) approaches unity, which is quite probable for a purely amorphous system. The peaks observed at distances less than 4 Å are assigned to specific distance of the closely coupled atoms. Pair correlation function is useful in the structural investigations of both solid and liquid packing (local structure) for studying specific interactions such as hydrogen bonding. Chitosan interacts with PVA through -OHCH2 and/or -NH2 of CS because of the presence of hydroxyl group in PVA. The g(r) function of about 2.5 for PVA/CS (1:1) indicates a strong interaction between these groups at lower concentration of PVA in the blend (Figure 12). This can be explained on the basis of free rotation of hydroxymethyl groups of CS, which may interact with the -OH groups of PVA. The probability of interaction between -NH2 of CS and -OH of PVA would be small at these compositions. However, at higher concentrations of PVA (PVA/CS, 80:20), the interactions between -OH groups of CS and PVA are at minimum, indicating their immiscibility. Conclusions The MD simulations using COMPASS force field could effectively allow us to understand the miscibility or immiscibility nature of PVA and CS blends. Atomistic and mesoscopic simulations indicated blend miscibility above 40 wt. % of PVA in the blend. MD calculations were supported by the experimental DMTA data. Plot of cohesive energy density vs weight fraction of PVA showed that values of CED increased as the PVA content of the blend increased. The χ values for blends below 0.18 (i.e., below χcritical), which was observed for 40:60 PVA/CS blend further supported the simulation results; this also agreed quantitatively with the mesoscopic simulation. The order parameter, being close to 0.1, supported the miscible nature of CS and PVA above 40 wt. % of PVA in the blend. Mesoscopic density slices clearly displayed the phase separations between

PVA and CS interfaces. From free-energy density calculations, it was observed that the system approached stability with time. MD simulations indicated that 50:50-90:10 blends of PVA: CS are immiscible, since χ is higher than χcritical, for these systems, whereas 10:90-40:60 blends of PVA:CS are miscible. The present calculations show that at low compositions of PVA, interactions are established with hydroxymethyl group of CS. MesoDyn computations as well as radial distribution function calculations show consistency with the experimental DMTA data as well as χ calculations. Acknowledgment. We thank University Grants Commission (UGC), New Delhi, India for a major funding (F1-41/2001/ CPP-II) to establish Center of Excellence in Polymer Science at Karnatak University, Dharwad. We also thank Mr. Ashok Betraj (Apsara Innovations, Bangalore, India), Mr. Anand Gupta of Accelyrs Inc., (San Diego CA), and Dr. Mori (Japan) for assistance. This paper is Center of Excellence in Polymer Science communication no. 111. References and Notes (1) van der Lubben, I. M.; Verhoef, J. C.; van Aelst, A. C.; Borchard, G.; Junginger, H. E. Biomaterials 2001, 22, 687-694. (2) Genta, I.; Perugini, P.; Pavanetto, F. Drug DeV. Ind. Pharm. 1998, 24, 779-784. (3) Agnihotri, S. A.; Aminabhavi, T. M. J. Controlled Release 2004, 96, 245-259. (4) Agnihotri, S. A.; Mallikarjuna, N. N.; Aminabhavi, T. M. J. Controlled. Release 2004, 100, 5-28. (5) Dubey, V.; Pandey, L. K.; Saxena, C. J. Membr. Sci. 2005, 251, 131-136. (6) Anjali Devi, D.; Smitha, B.; Sridhar, S.; Aminabhavi, T. M. J. Membr. Sci. 2005, 262, 91-99. (7) Mochizuki, A.; Amiya, S.; Sato, Y.; Ogawara, H.; Yamashita, S. J. Appl. Polym. Sci. 1989, 37, 3385-398. (8) Smitha, B.; Sridhar, S.; Khan, A. A. Macromolecules 2004, 37, 2234-2239. (9) Smitha, B.; Sridhar S.; Khan, A. A. Eur. Polym. J. 2005, 41, 18591866. (10) Huang, Y.; Onyeri, S.; Siewe, M.; Moshfeghian, A.; Madihally, S. V. Biomaterials 2005, 26, 7616-7627. (11) Kurkuri, M. D.; Aminabhavi, T. M. J. Controlled Release 2004, 96, 9-20. (12) Yu, J.; Lee, C. H.; Hong, W. H. Chem. Eng. Process. 2002, 41, 693-698. (13) Srinivasa, P. C.; Ramesh, M. N.; Kumar, K. R.; Tharanathana, R. N. Carbohydr. Polym. 2003, 53, 431-438. (14) Kim, S. J.; Park, S. J.; Kim, S. I. React. Funct. Polym. 2003, 55, 53-59. (15) Naidu, B. V. K.; Sairam, M.; Raju, K. V. S. N.; Aminabhavi, T. M. Carbohydr. Polym. 2005, 61, 52-60. (16) Miura, K.; Kimura, N.; Suzuki, H.; Miyashita, Y.; Nishio, Y. Carbohydr. Polym. 1999, 39, 139-144. (17) Fraaije, J. E. M. J. Chem. Phys. 1993, 99, 9202. (18) Fraaije, J. G. E. M.; van Vlimmeren, B. A. C.; Maurits, N. M.; Postma, M.; Evers, O. A.; Hoffman, C.; Altevogt, P.; Goldbeck, W. G. J. Chem. Phys. 1997, 106, 4260. (19) Groot, R. D.; Warren, P. B. J. Chem. Phys. 1997, 107, 4423. (20) Groot, R. D.; Madden, T. J. J. Chem. Phys. 1998, 108, 8713. (21) Helfand, E.; Wasserman, Z. R. Macromolecules 1976, 9, 879. (22) Leibler, L. Macromolecules 1980, 13, 1602. (23) Melenkevitz, J.; Muthukumar, M. Macromolecules 1991, 24, 4199. (24) Matsen, M. W.; Schick, M. Phys. ReV. Lett. 1994, 72, 2660; Macromolecules 1994, 27, 6761. (25) Matsen, M. W.; Bates, F. S. Macromolecules 1996, 29, 1091. (26) Jawalkar, S. S.; Adoor, S. G.; Sairam, M.; Mallikarjuna, N. N.; Aminabhavi, T. M. J. Phys. Chem. B. 2005, 109, 15611-15620. (27) Chuang, W. Y.; Young, T. H.; Yao, C. H.; Chiu, W. Y. Biomaterials 1999, 20, 1479-1487. (28) Yang, J. M.; Sua, W. Y.; Leu, T. L.; Yang, M. C. J. Membr. Sci. 2004, 236, 39. (29) Sandoval, C.; Castro, C.; Gargallo, L.; Radic, D.; Freire, J. Polymer 2005, 46, 10437. (30) Castro, C.; Gargallo, L.; Leiva, A.; Radic, D. J. Appl. Polym. Sci. 2005, 97, 1953.

Compatibility of Poly(vinyl Alcohol) and Chitosan Blends (31) Adoor, S. G.; Manjeshwar, L. S. Krishna Rao, K. S. V.; Naidu, B. V. K.; Aminabhavi, T. M. J. Appl. Polym. Sci. 2006, 100, 24152421. (32) Rigby, D.; Sun, H.; Eichinger, B. E. Polym. Int. 1997, 44, 311. (33) Greengard, L.; Rokhlin, V. I. J. Comput. Phys. 1987, 73, 325. (34) Schmidt, K. E.; Lee, M. A. J. Stat. Phys. 1991, 63, 1223. (35) Ding, H. Q.; Karasawa, N.; Goddard, W. A., III. J. Chem. Phys. 1992, 97, 4309. (36) Rappe, A. K.; Goddard, W. A. J. Phys. Chem. 1991, 95, 3358. (37) Ewald, P. P. Ann. Physik 1921, 64, 253. (38) Grulke, E. A. Solubility Parameter Values. In Polymer Handbook, 4th Ed.; Brandrup, J., Immergut, E. H., Eds.; Wiley: New York, 1999; pp 705.

J. Phys. Chem. B, Vol. 111, No. 10, 2007 2439 (39) Flory, P. J. Principles of Polymer Chemistry; Cornell University Press: Ithaca, New York, 1953. (40) Krause, S. Polymer-Polymer Compatibility in Polymer Blends; Paul, D. R., Newman, S., Eds; Academic Press: New York, 1978; Vol. 1. (41) Zhang, M.; Choi, P.; Sundararaj, U. Polymer 2003, 44, 1979. (42) Fan, C. F.; Cagin, T.; Chen, Z. M.; Smith, K. A. Macromolecules 1994, 27, 2383. (43) Aminabhavi, T. M.; Munk, P. Macromolecules 1979, 12, 607. (44) Flory, P. J. Statistical Mechanics of Chain Molecules; Hanser: Munich, Germany, 1989. (45) Groot, R. D.; Madden, T. J.; Tildesley, D. J. J. Chem. Phys. 1999, 110, 9739. (46) Ravindra, R. Kameshwara, R. K.; Khan, A. A. Carbohydr. Polym. 1998, 36, 121.