Classical Force Field for Hydrofluorocarbon Molecular Simulations

Article pubs.acs.org/JPCA

Classical Force Field for Hydrofluorocarbon Molecular Simulations. Application to the Study of Gas Solubility in Poly(vinylidene fluoride) V. Lachet,*,† J.-M. Teuler,‡ and B. Rousseau‡ †

IFP Energies nouvelles, 1−4 Avenue de Bois-Préau, 92852 Rueil-Malmaison, France Laboratoire de Chimie Physique, Université Paris-Sud, UMR 8000 CNRS, 91405 Orsay, France

‡

S Supporting Information *

ABSTRACT: A classical all-atoms force field for molecular simulations of hydrofluorocarbons (HFCs) has been developed. Lennard-Jones force centers plus point charges are used to represent dispersion−repulsion and electrostatic interactions. Parametrization of this force field has been performed iteratively using three target properties of pentafluorobutane: the quantum energy of an isolated molecule, the dielectric constant in the liquid phase, and the compressed liquid density. The accuracy and transferability of this new force field has been demonstrated through the simulation of different thermophysical properties of several fluorinated compounds, showing significant improvements compared to existing models. This new force field has been applied to study solubilities of several gases in poly(vinylidene fluoride) (PVDF) above the melting temperature of this polymer. The solubility of CH4, CO2, H2S, H2, N2, O2, and H2O at infinite dilution has been computed using test particle insertions in the course of a NpT hybrid Monte Carlo simulation. For CH4, CO2, and their mixtures, some calculations beyond the Henry regime have also been performed using hybrid Monte Carlo simulations in the osmotic ensemble, allowing both swelling and solubility determination. An ideal mixing behavior is observed, with identical solubility coefficients in the mixtures and in pure gas systems.

1. INTRODUCTION Poly(vinylidene fluoride) (PVDF or PVF2) is a partially fluorinated polymer produced by the polymerization of 1,1difluoroethylene. It is a semicrystalline polymer with a glass transition temperature close to 238 K and a melting temperature ranging from 433 to 453 K. Due to its good barrier properties, this polymer is used in many industrial areas to protect materials from gas or liquid contamination. For instance, PVDF is widely used in the offshore oil and gas industry in order to ensure the impermeability of flexible pipes, both from the seawater outside and from the fluids conveyed inside. In order to prevent gases or liquids from damaging the structure of the pipes, it is necessary to acquire a good understanding of the permeation phenomenon and to develop reliable models to predict permeation behaviors in polymeric materials. Molecular simulation is an attractive tool to fulfill these requirements, as it relies on methods with few assumptions and offers the possibility to reveal the microscopic mechanisms that can play an important role in the barrier properties. Thanks to the increase of computing power and to the advances in simulation algorithms, several authors have investigated the modeling of gas solubility in different polymers at a molecular level using Monte Carlo or molecular dynamics approaches (see review paper by Theodorou.1) However, to the best of our knowledge, none of these studies deals with PVDF. This can be explained by the fact that, from the molecular simulation point of view, PVDF is a complex system with specific structural and © 2014 American Chemical Society

chemical features. Even small oligomers of PVDF have been little studied in the literature due to the lack of available atomistic force fields. Unlike perfluorinated alkanes that have been the subject of many detailed molecular simulation studies,2−10 we have identified for hydrofluorocarbons (HFCs), i.e., partially fluorinated hydrocarbons, only two united-atoms (UA) force fields proposed by Hariharan and Harris11 and by Shin et al.,12 respectively, and a single all-atoms (AA) force field proposed by Byutner and Smith.13 All three models are nonpolarizable models. The models of Hariharan and Shin are nonelectrostatic UA force fields; i.e., in such models, CH2, CH3, CF2 and CF3 groups are represented by single force centers without electrostatic charges. These two potentials differ only by the numerical values of the involved parameters and by the combination rules used. Hariharan and Harris use geometric rules while Shin et al. work with the Lorentz−Berthelot rules. Note that these two force fields have been developed for the study of partially fluorinated compounds with a fully fluorinated part and a fully hydrogenated part, and not for compounds with a strict alternation of CH2 and CF2 groups as it is the case for the compounds of interest in this study. The model proposed by Byutner and Smith is a AA force field with all intramolecular parameters and electrostatic charges derived from quantum Received: July 10, 2014 Revised: December 5, 2014 Published: December 5, 2014 140

dx.doi.org/10.1021/jp506895p | J. Phys. Chem. A 2015, 119, 140−151

The Journal of Physical Chemistry A

Article

calculations performed on small oligomers of poly(vinylidene fluoride). Nonbonded interactions have been parametrized from quantum calculations performed on CF4 + CF4 and CF4 + CH4 dimers. This model requires a quantum calculation for each studied molecule followed by an empirical adjustment of the obtained electrostatic charges to account implicitly for polarization effects. Hence, this force field is not directly transferable for the study of a new hydrofluorocarbon. However, as far as poly(vinylidene fluoride) is concerned, it is probably the most accurate model available for the study of the polymer melted or amorphous phase properties. For completeness, it is worth mentioning that some atomistic force fields are also available in the literature for PVDF crystalline polymorphs,14−18 but these models fail to accurately reproduce conformational energies and are thus not designed to study amorphous properties.13 The limited number of simulation studies carried out on HFCs may be linked to the limited number of available experimental data in the open literature. We indeed found in the open literature experimental data only for the smaller molecules of the series, namely, 1,1,1-trifluoroethane (TFE), 1,1,1,3,3,3-hexafluoropropane (HFP), 2,2-difluoropropane (DFP), 2,2,4,4-tetrafluoropentane (TFP), and 1,1,1,3,3-pentafluorobutane (PFB). This latter molecule is the most studied with available measurements of compressed liquid densities,19−21 saturated liquid densities,20,22 vapor pressures,19 surface tensions,22 vapor-phase thermal conductivities,19 liquidphase thermal conductivities,22 and dielectric constants in the liquid phase21 as well as viscosities.22 However, there are almost no data for the partially fluorinated hydrocarbons with longer chains. In this work, we propose a new force field for HFCs and PVDF, based on modifications and reparametrizations of the Byutner force field. The performance of this force field has been evaluated on its ability to reproduce the thermophysical properties of PFB, as well as the PVT properties of amorphous PVDF. These evaluations have been carried out using either Monte Carlo simulations, molecular dynamics simulations, or hybrid Monte Carlo simulations. The new proposed force field has also been used to study the thermodynamic behavior of several gas + PVDF mixtures. The solubility of different gases (CH4, CO2, H2S, H2, N2, O2, and H2O) in PVDF at infinite dilution has been computed using test particle insertions in the course of a NpT hybrid Monte Carlo simulation. For CH4 and CO2, some calculations beyond the Henry regime have also been performed using hybrid Monte Carlo simulations in the osmotic ensemble. These calculations concern both the solubility of pure gases as well as the solubility of the equimolar CH4 + CO2 binary mixture. This work is organized as follows: the molecular simulation methodologies employed in this work are described in section 2, with emphasis on the hybrid Monte Carlo approach used to efficiently sample the numerous internal degrees of freedom of PVDF chains. The proposed force field is presented in section 3, highlighting its accuracy and its transferability to different fluorinated compounds and to different property calculation. In section 4, we present our simulation results on gas solubility in melted PVDF. Section 5 gives our conclusions.

(MC) simulations in the canonical ensemble (constant number of particles, volume, and temperature) or in the isothermal− isobaric ensemble (constant number of particles, pressure, and temperature). During MD simulations, equations of motion were integrated using the velocity-Verlet algorithm23 with a time step of 1 fs, in the canonical ensemble using Nosé− Hoover thermostat24 or in the NpT ensemble based on the external bath coupling algorithm from Berendsen et al.25 During MC simulations, the following Monte Carlo moves, with their corresponding attempt probabilities, were used: molecular rigid translations (20%), molecular rigid rotations (20%), atomic translations (30%), internal rotations (29%), and volume changes (1% in the NpT ensemble). The internal rotation or flip move is the rotation of one force center around the axis formed by two nearest neighbors. This move, originally designed to sample linear parts of molecules, has been extended to branched or AA molecules.26 In such cases, not only the chosen force center is rotated around a selected axis but also all substituents or ramifications of this force center. The amplitude of molecular translations, molecular rotations, atomic translations, internal rotations, and volume changes was adjusted during the simulation to achieve an acceptance ratio of 40%. In most cases, the MC simulations involved an equilibration run of around 30 million Monte Carlo steps (one step corresponding to a single Monte Carlo move) and a production run of 70 million Monte Carlo steps. A few MC simulations in the NVT Gibbs ensemble27 (GEMC) were also performed in order to assess liquid−vapor equilibrium properties of pure compounds. In such simulations, the two phases are simulated in individual boxes without an explicit interface, and particle transfers between boxes are performed in order to ensure phase equilibrium. The different MC moves that were used during these simulations are molecular rigid translations (15%), molecular rigid rotations (15%), atomic translations (30%), internal rotations (19%), transfers with configurational bias (20%), and volume changes (1%). In MD and MC simulations, a spherical cutoff equal to half of the simulation box was used to calculate Lennard-Jones interactions, and classical tail corrections were employed. For long-range electrostatic energy, the Ewald sum technique was used with the Gaussian width and number of reciprocal vectors parameters chosen such that the relative error in the reciprocal and real space energy be of the order of 10−5. Depending on the studied system, a total number of 150−300 molecules was used for the simulations of HFCs. For pure PVDF simulations, a box with 15 chains having each 70 carbon atoms was considered. 2.2. Hybrid Monte Carlo Simulations. The study of gas + PVDF mixtures was carried out using Monte Carlo simulations performed at fixed pressure and temperature. To evaluate the Henry constant of the different studied gases, we used the Widom method28 where insertion tests of a phantom gas molecule are carried out in a pure PVDF matrix during a NpT Monte Carlo simulation. The Henry constant of the gas is obtained by the following expression:

2. METHODOLOGY 2.1. Standard Molecular Dynamics and Monte Carlo Simulations. The study of HFCs and pure PVDF was carried out using either molecular dynamics (MD) or Monte Carlo

where ρN is the number density of the system, ΔU is the potential energy change associated with the insertion test, V is the simulation box volume, kB is the Boltzmann constant, and β = 1/kBT.

⎡ ⟨V exp( −β ΔU )⟩ ⎤ H = ρN kBT ⎢ ⎥ ⎣ ⎦ ⟨V ⟩

141

(1)



Article

an average acceptance ratio of the MD move of 40%, and (iii) the acceptance probability of this MD move is given by

To study gas solubilities in PVDF beyond the Henry regime, we used Monte Carlo simulations in the osmotic ensemble. In this ensemble, a single phase (the polymer-rich phase in our case) is considered with a fixed amount of polymer chains, Np, in equilibrium with a gas or a gas mixture with imposed chemical potential μi for each gas species. Temperature T and pressure P are also imposed, whereas the number of molecules of the different gas species and the volume of the simulation box are allowed to vary. The configurational part of the partition function in this ensemble is given by the following expression:

Pacc(i→j) = min{1, exp{−β Δ/)}

where / is the Hamiltonian of the system. A hybrid Monte Carlo simulation of a gas + PVDF mixture involved the following MC moves with the following attempt probabilities: gas insertions/deletions for osmotic simulations or equivalently gas insertion tests for NpT simulations (35%), volume changes (5%), rigid molecular translations (5%), rigid molecular rotations (5%), and MD moves (50%). A total number of 1− 2 million MC steps was performed for each studied system. For these simulations, the PVDF matrix was modeled using 15 chains having each 70 carbon atoms. All of the simulations performed in this work were carried out with locally developed simulation packages Newton41 and Gibbs.42

n

Θ=

∑ ∑ ∑ exp[−βU + V

Ni

U

∑ βμj Nj − βPV ] j=1

(3)

(2)

where Ni is the number of molecules of gas species i and U is the potential energy of the system. This expression is consistent with those found in previous works.29,30 A preliminary simulation of the gas mixture at T, P, and gas composition yi is performed using Widom test particle insertions in order to calculate the chemical potential μi of each species in the gas phase. Whatever the statistical ensemble that is used, an efficient sampling of the configurational space of polymer chains is required in these Monte Carlo simulations. Efficient algorithms need to be employed in order to explore the numerous conformations that can be adopted by these long-chain molecules. Most of the Monte Carlo simulations of polymeric systems published so far make use of Monte Carlo moves specifically designed for long-chain molecules, such as concerted rotation,31 reptation with configurational bias,32 double bridging, or internal double bridging.33−36 Such moves are very efficient for linear chains, but their extension to more complex systems, such as branched or AA chains among others, is not trivial. In the present work, to face sampling problems in PVDF simulations, we decided to use hybrid Monte Carlo.37,38 Hybrid Monte Carlo consists of using a short molecular dynamics run as a Monte Carlo move, in addition to other usual moves. Such an approach presents several advantages: (i) it promotes the sampling of the configurational space by performing cooperative moves, (ii) it allows one to sample all intramolecular degrees of freedom of complex molecules, for which there is not always suitable Monte Carlo moves, (iii) the MD moves performed in a hybrid simulation are force biased, i.e., driven by the force field, which can improve the sampling efficiency, and (iv) by using a Monte Carlo algorithm, hybrid Monte Carlo prevents truncation errors that are inherent to integration methods using finite differences. In order to respect microreversibility, the molecular dynamic steps performed in the course of an hybrid Monte Carlo simulation must respect two conditions. First, time reversibility must be observed. Thus, the model used should not contain any holonomic constraint and care must be taken when using thermostat or barostat.39,40 Second, the MD run must be allowed to run forward and backward in time; i.e., the sign of the time step must be randomly chosen before the MD move to ensure the possibility of performing one move and the corresponding reverse move as required in any MC simulation. Our implementation of hybrid Monte Carlo is the following: (i) random generation of initial moments from a Gaussian distribution at temperature T, (ii) use of 50 microcanonical molecular dynamics cycles with an adjustable time step (with randomly chosen sign) so as to reach

3. FORCE FIELD DEVELOPMENT A new AA force field is proposed for HFCs, based on modifications and reparametrizations of the Byutner force field.13 The Byutner model is a quantum chemistry based force field developed in 2000 to study PVDF bulk properties such as PVT properties and chain dimensions. In this model, intramolecular parameters and electrostatic charges are derived from ab initio calculations performed on four small oligomers of poly(vinylidene fluoride): 1,1,1,3,3-pentafluorobutane (PFB), 1,1,1,3,3,5,5,5-octofluoropentane (OFP), 2,2,4,4-tetrafluoropentane (TFP), and 2,2,4,4,6,6-hexafluoroheptane (HFH). Different sets of partial charges have been determined for the different studied molecules, and an empirical adjustment of these charges is proposed to account implicitly for polarization effects in real systems. All of the atoms have a nonzero electrostatic charge. Nonbonded interactions have been parametrized from quantum calculations performed on CF4 + CF4 and CF4 + CH4 dimers. These nonbonded interactions are modeled through a Buckingham exponential/6 potential plus a Coulombic contribution. In this force field, intramolecular interactions are composed of bonding and bending contributions, both modeled using harmonic potentials, and a torsional contribution modeled with a Fourier series. Within a given molecule, van der Waals and electrostatic interactions between atoms separated by exactly three bonds, i.e., 1−4 intramolecular interactions, are systematically taken into account without weighting factors. Finally, van der Waals cross-interactions, i.e., between unlike atoms, are not calculated using a combination rule but are given individually for each pair of atom types considered. This model is probably the most accurate one to study bulk properties of HFCs and PVDF. However, an empirical adjustment of electrostatic charges is required for each studied molecule. Therefore, we propose in the present work a new force field for HFCs and PVDF, with special attention paid to the predictive ability and transferability of the model. This new force field is based on the standard Lennard-Jones dispersion− repulsion potential in order to allow a widespread use with other Lennard-Jones models if fluid mixtures are considered. 3.1. Functional Forms and Associated Parameters. In this new force field, the dispersive−repulsive interactions between two force centers i and j belonging either to two different molecules (intermolecular interactions) or to the same molecule but separated by three bonds or more (intramolecular 142



Article

parameters involved in the different energy contributions of the force field. For instance, the 1−4 intramolecular electrostatic and Lennard-Jones interactions lead to a dependency of the torsional parameters to the values of electrostatic charges and Lennard-Jones parameters. Similarly, partial electrostatic charges cannot be determined on the basis of dipole moments, since these dipole moments are very sensitive to the conformation of the molecule, the conformation of the molecule being dependent on the whole set of intramolecular parameters. In order to minimize the number of parameters to be adjusted, bonded parameters corresponding to bond and bend vibrations as well as torsional parameters associated with C− C−C−C, H−C−C−C, and H−C−C−F dihedral angles have been reused from Byutner parametrization. Similarly, the Lennard-Jones parameters of the OPLS-AA force field have been reused for carbon and hydrogen atoms.49 The determination of remaining parameters (atom charges, fluorine LJ parameters, and F−C−C −F torsional parameters) has been performed iteratively until getting a set of parameters that allows reproducing accurately three target properties: the energy of an isolated molecule computed using density functional theory (DFT), the dielectric constants in the liquid phase, and the compressed liquid densities. This iterative procedure has been initiated from Byutner force field for the atom charges and F−C−C−F torsional parameters. The initial values for the LJ 12−6 parameters for fluorine were obtained by fitting the required LJ parameters to the Buckingham exponential/6 potential proposed by Byutner. Using the aforementioned values for the fixed parameters and the aforementioned initial values for the parameters to be optimized, we have conducted the optimization procedure as follows: (i) adjustment of fluorine σ and ϵ LJ parameters with current values for atom charges and F−C−C−F torsional parameters using the optimization procedure proposed by Bourasseau et al.,51 (ii) adjustment of the F−C−C−F torsional potential parameters using the current values for fluorine LJ parameters and atom charges [this adjustment was performed by fitting the F−C−C−F torsion parameters so that the total energy of an isolated molecule calculated with the force field matches the DFT energy (see section 3.2.1 for more details)], and (iii) empirical scaling of the atomic electrostatic charges using the current values for fluorine LJ parameters and F−C− C−F torsional parameters in order to best reproduce the dielectric constants in the liquid phase (see section 3.2.4 for more details). Only PFB was considered in this optimization procedure since experimental data in sufficient quantities are not available for other HFCs. The bonded and nonbonded parameters involved in this new force field are given in Table 1 and Table 2, respectively. It is worth specifying that the form of the torsional potential given in Byutner work differs from the one used in this work, although simple analytical expressions allow one to convert from one parameter set to another (see Supporting Information (SI)). As proposed by Byutner, H−C− C−C and H−C−C−F torsional energies are not explicitly accounted for in our force field; the corresponding parameters are thus all set to zero. 3.2. Force Field Validation. The performance of our new force field has been evaluated on its ability to reproduce several thermophysical properties of PFB (see Figure 1). The PFB molecule has been chosen given the availability of experimental data. For compressed liquid densities, a comparison with the

interactions) are described through a 12−6 Lennard-Jones (LJ) potential: ⎡⎛ ⎞12 ⎛ ⎞6 ⎤ σij σij ULJ(rij) = 4ϵij⎢⎢⎜⎜ ⎟⎟ − ⎜⎜ ⎟⎟ ⎥⎥ r ⎝ rij ⎠ ⎦ ⎣⎝ ij ⎠

(4)

where rij, ϵij, σij are the interparticle distance, the LJ well depth, and the LJ size, respectively. The Kong43 combining rules are used to compute the cross-LJ parameters: ϵijσij12

13 ⎡ ⎛ ϵ σ 12 ⎞1/13⎤ ϵiiσii12 ⎢ jj jj ⎥ ⎟ = 13 ⎢1 + ⎜⎜ , 2 ϵiiσii12 ⎟⎠ ⎥ ⎝ ⎣ ⎦

ϵijσij 6 = (ϵiiσii 6 ϵjjσjj 6)1/2 (5)

The choice of the Kong rules is due to the fact that they seem to offer better results than those of Lorentz−Berthelot for simulation of liquid−vapor phase equilibrium of mixtures.44 Furthermore, the use of a Lennard-Jones potential associated with these combining rules allows us to be consistent with our previous works on gas permeation in polyethylene.45−47 The electrostatic interactions between molecules (intermolecular electrostatic interactions) or between two parts of the same molecule (intramolecular electrostatic interactions) are calculated as the sum of the interactions between pairs of point charges i and j located on the atoms using the Coulomb potential: qiqj Uelec(rij) = 4π ϵ0rij (6) where qi is the magnitude of charge i, ϵ0 the vacuum permittivity, and rij the distance between charges i and j. The Ewald sum technique is used to calculate these electrostatic interactions. For intramolecular contributions, interactions between charges separated by three bonds and more are taken into account. We decided to take into account explicitly the 1−4 interactions, both for dispersive−repulsive and electrostatic contributions. This approach, sometimes coupled with the use of a scaling factor, is commonly used in the literature when multifunctional molecules are considered.48−50 The bonded (bonding, bending, and torsional) interactions are given by 1 k bond(r − r0)2 2

Ubond(r )/kB =

(7)

where kB is the Boltzmann constant, kbond is the bonding constant, r is the distance between the two bonded atoms, r0 is the equilibrium bonded distance, and Ubend(θ )/kB =

1 k bend(θ − θ0)2 2

(8)

where kbend is the bending constant, θ is the bending angle associated with atoms separated by two bonds, θ0 is the equilibrium bending angle; 8

Utors(χ )/kB =

∑ an[cos χ ]n n=0

(9)

where the χ angle has been defined such that χ = ϕ + π, with ϕ being the dihedral angle associated with atoms separated by three bonds and an the corresponding torsional coefficients. The main difficulty when dealing with molecules bearing several functional groups is the interdependency of the 143



Article

range of 0−180°, with increments of 10°. For each angle value, a geometry optimization was performed and the corresponding DFT energy was calculated and used for the fitting of the C− C−C−C torsional energy with the analytical expression given in eq 9. The accuracy of this adjustment is shown in Figure 2

Table 1. Bonded Force Field Parameters for HFCs and PVDF kbond (K/Å2)

r0 (Å)

bonding

ref

C−C C−H C−F bending

1.534 1.085 1.357 θ0 (deg)

310889 329709 502514 kbend (K/rad2)

13 13 13 ref

C−C−F F−C−F C−C−C C−C−H C−C−F torsion

107.74 105.27 118.24 108.45 109.27

90579 120772 80817 43176 38748

13 13 13 13 13 ref

an (K)

C−C−C−C

F−C−C−C

a0 a2 a4 a6 a8 a0 a2 a4 a6 a8

= = = = = = = = = =

920.89 −1323.47 221.42 402.57 0.00 1069.96 −1764.99 1676.33 −638.83 0.00

a1 a3 a5 a7

= = = =

1841.78 −5042.25 3421.89 0.00

13

a1 a3 a5 a7

= = = =

1035.05 −3636.48 2943.89 0.00

this work

Figure 2. Energy of a PFB molecule in vacuum as a function of the C− C−C−C dihedral angle. An angle equal to 0° corresponds to a cis conformation, whereas an angle equal to 180° corresponds to a trans conformation. Comparison between values derived from DFT calculations (red open circles, red line, cubic spline fit of the DFT data) and from force field calculations (blue squares).

where the energy obtained by summing all of the intramolecular contributions of our force field (bonding, bending, torsional, intramolecular Lennard-Jones, and electrostatic energies) is compared to the DFT energy minus the energy of the minimum energy conformer. 3.2.2. PFB Compressed Liquid Densities. Compressed liquid densities under different pressure and temperature conditions have been calculated using Monte Carlo simulations in the NpT ensemble. Two temperatures have been investigated, namely, 289.15 and 393.15 K, for pressures in the range of 0.1−10 MPa. The calculated densities at 289.15 K are shown in Figure 3,

Table 2. Nonbonded Force Field Parameters for HFCs and PVDF atom

ϵ (K)

σ (Å)

ref

q (e)

C/CF3 C/CF2 C/CH3 C/CH2 H F

33.212 33.212 33.212 33.212 15.097 30.109

3.500 3.500 3.500 3.500 2.500 2.983

49 49 49 49 49 this work

+0.949000 +0.765000 −0.776875 −0.650250 +0.225875 −0.283250

ref this this this this this this

work work work work work work

Figure 3. Liquid densities of PFB at 289.15 K. Comparison between experimental data of Marrucho et al.19 (black filled circles) and simulated data obtained with different force fields: this work (blue open circles), Byutner and Smith13 (green open squares), Hariharan and Harris11 (red open diamonds), and Shin et al.12 (orange open triangles).

Figure 1. PFB (1,1,1,3,3-pentafluorobutane) molecule.

predictions obtained with the other available force fields, namely, Byutner, Hariharan, and Shin, is also proposed. 3.2.1. Energy of an Isolated PFB Molecule. In order to validate the intramolecular contributions of our force field, a series of DFT calculations was performed on an isolated PFB molecule using Gaussian with the B3LYP functional and the 6311G(d,p) basis set. We have chosen this method because it proved to reproduce energy values of a reference and timeconsuming ab initio method, CCSD(T)/aug-cc-pvdz, for a limited number of PFB conformers. The DFT energy was calculated for different C−C−C−C torsional angles in the

together with the experimental data of Marrucho et al.19 and the simulated data obtained with other force fields (Shin, Hariharan, and Byutner). A very good agreement between experimental and calculated values is achieved when using our new force field, with an average deviation of 0.7%. On the contrary, UA force fields of Hariharan and Shin lead to a drastic underestimate of densities in the liquid phase with an average deviation of 12% with the Hariharan model and an average deviation of 23% with the Shin model. The results obtained 144



Article

with the Byutner force field are more satisfying, but still significant deviations with experiments are observed with an average overestimate of the densities of 6.1%. The same behavior is observed at the other studied temperature (see the SI). 3.2.3. PFB Liquid−Vapor Equilibrium Properties. Figure 4 shows the saturated liquid and vapor densities of PFB. In this

ρl − ρv = B(T − Tc)β

(11)

where β is a universal exponent, taken to be equal to 0.325, ρl and ρv are the densities of the liquid and vapor phases, respectively, and T is the temperature. A and B are adjustable parameters. The critical point we have obtained (ρc = 459.79 kg/m3 and Tc = 469.35 K) is in reasonable agreement with the one given in the NIST webbook53 (ρc = 474 kg/m3 and Tc = 460 K). 3.2.4. Dielectric Constants in PFB Liquid Phase. The validation of the partial electrostatic charges proposed in our force field has been done through molecular dynamics calculations of the dielectric constant ϵr of liquid PFB at three temperatures, namely, 263, 283, and 303 K, under a pressure of 1 MPa. These calculations are based on the fluctuations of the total dipole moment of the simulation box, using the following expression (in SI units):54,55 N

ϵr − 1 = Figure 4. Comparison between experimental and simulated densities of PFB in the liquid−vapor coexistence region. NIST correlation52 (black line), Marrucho et al.19 experimental data (red open squares), Fröba et al.22 experimental data (green open diamonds), and GEMC simulations (blue filled circles). Notice that the critical point from Marrucho is an estimate based on Dohrn correlation; it is not measured (see details in the Marrucho et al. work19).

N

⟨|∑i = 1 μi |2 ⟩ − ⟨|∑i = 1 μi |⟩2 3ϵ0VkBT

(12)

where μi is the dipole moment of molecule i. As can be seen in Figure 6, the molecular dipole moment distribution exhibits

figure, predictions obtained by GEMC simulations with the new force field are compared to experimental data of Marrucho et al.19 and Fröba et al.,22 as well as to the NIST correlation.52 The numerical values corresponding to our GEMC simulations are given in the SI. The PFB saturated vapor pressures, computed using the Virial equation, are shown in Figure 5, where they are Figure 6. Distribution of molecular dipole moments in a molecular dynamics simulation of PFB at 300 K using the new force field.

two peaks, the first one around 2.5 D and the second one around 4.5 D. μi varies significantly as molecules undergo conformational changes governed by the two energy minima in the intramolecular energy curve (see Figure 2). Therefore, during the optimization procedure of the electrostatic charges, we avoided using the configurational average of the molecular dipole moment. We preferred to rely on the dielectric constant, since it reflects changes in molecular dipole moments. The dielectric constants obtained in the NpT ensemble are presented in Table 3 where they are compared to experimental

Figure 5. Comparison between experimental and simulated vapor pressures of PFB. NIST correlation52 (black line), Marrucho et al.19 experimental data (red open squares), and GEMC simulations (blue filled circles).

Table 3. Calculated (This Work) and Experimental (Ribeiro et al.21) Densities and Dielectric Constants of PFB at 1 MPa

compared to experimental data of Marrucho et al.19 and to the NIST correlation.52 For both studied properties, a good agreement is found (average absolute deviations of 1.7% for saturated liquid densities and of 15.1% for vapor pressures). The accuracy of our new force field in the vicinity of the PFB critical point was also investigated. The critical coordinates Tc and ρc of PFB were estimated by a least-squares fit of the law of rectilinear diameters (eq 10) and a scaling law (eq 11): ρl + ρv = ρc + A(T − Tc) (10) 2

T (K)

ρsim‑NpT (kg·m−3)

ρexp (kg·m−3)

ϵsim‑NpT r

ϵexp r

263 283 303

1347.4 ± 1.0 1307.7 ± 0.9 1263.3 ± 1.1

1334.90 1294.23 1252.34

14.7 ± 0.8 12.4 ± 0.6 11.7 ± 0.3

14.90 12.90 11.23

data of Ribeiro et al.21 The calculated densities are also included in this table and compared to the Ribeiro data. An accurate evolution with temperature of both densities and dielectric constants can be observed. 3.2.5. PFB Viscosity along the Saturation Curve. Although the new force field was optimized to reproduce thermodynamic 145



Article

A good transferability of our force field to these three new molecules is found, with an average deviation between simulated and measured densities of 5.9% for TFE, 2.8% for HFP, and 3.2% for TFP. The ability of the force field to predict saturated properties of 2,2-difluoropropane has also been investigated using NVT GEMC. The obtained values are presented in the SI, where they are compared to experimental data from Henne et al.62 Finally, the accuracy of the force field toward dynamical properties has been evaluated (using the same methodology as described in section 3.2.5) by computing the saturated liquid viscosity of HFP for three temperatures in the range of 273−333 K. Experimental data from Laesecke and Defibaugh63 and simulated values are reported in Table 5.

properties, it is interesting to see if it can be used for transport properties predictions. We computed PFB shear viscosity at three different temperatures along the saturation curve and compared our values with Fröba et al.22 experimental data. The shear viscosity was computed from 30 ns molecular dynamics runs in the NVT ensemble and using the Green−Kubo formalism. More details on the procedure can be found elsewhere.41 Experimental and simulated viscosities are reported in Table 4 along with liquid experimental density Table 4. Experimental (Fröba et al.22) and Calculated (This Work) PFB Viscosities in the Liquid Phase at Three Temperatures along the Saturation Curve T (K)

ρexp (kg·m−3)

ηexp (mPa·s)

ηsim (mPa·s)

263 283 303

1328.3 1287.7 1245.6

0.822 0.613 0.468

1.22 ± 0.12 0.76 ± 0.05 0.58 ± 0.03

Table 5. Experimental (Laesecke and Defibaugh63) and Calculated (This Work) HFP Viscosities in the Liquid Phase at Three Temperatures along the Saturation Curve

values at saturation. Computed viscosities are overestimated with relative deviation ranging from 25% at high temperature up to 50% at 263 K. The poor result at low temperature may be due to the fact that the system is close to its melting point56 of 238.14 K. At such high liquid densities, viscosity prediction is always difficult.57,58 3.2.6. Properties of Other HFCs. We end the force field validation by a fully predictive use of this new model for determining the liquid densities of three other partially fluorinated hydrocarbons: 1,1,1-trifluoroethane (TFE), 1,1,1,3,3,3-hexafluoropropane (HFP) and 2,2,4,4-tetrafluoropentane (TFP). TFE was studied at 283.15 K for pressures in the range of 2−15 MPa. The simulated data for this compound are compared to experimental data of Nakamura et al.59 HFP was studied at two temperatures (293 and 373 K) for pressures up to 20 MPa, and simulated data are compared with measurements of Yin and Wu.60 TFP was studied at a single temperature (293 K) for pressures in the range of 0.1−8 MPa. A single experimental measurement from Stepanov et al.61 could be found for this species. The obtained results are shown in Figure 7.

T (K)

ρexp (kg·m−3)

ηexp (mPa·s)

ηsim (mPa·s)

273.38 303.10 333.15

1434.31 1341.36 1230.20

0.394 0.269 0.185

0.54 ± 0.10 0.33 ± 0.05 0.22 ± 0.03

Due to the lack of experimental data available, it was not possible to test this transferability to longer molecules such as small PVDF oligomers.

4. SIMULATION RESULTS 4.1. Pure Amorphous PVDF PVT Properties. In this section, we check that volumetric properties of amorphous PVDF are correctly predicted by our new force field. We realized molecular dynamics simulations using a NpT reversible algorithm24 and computed density at a pressure of 0.1 MPa and temperatures in the range of 373−493 K, i.e., below and above the PVDF melting temperature. The initial configuration for the whole set of simulations was built at high temperature and low density (0.1 g/cm3) followed by a rapid isothermal compression and an equilibration run of 5 ns at 493 K and up to 20 ns at 373 K. Production runs lasted for 8 ns at 493 K and up to 20 ns at 373 K. Computed densities are presented in Figure 8 along with statistical errors. As can be seen, no crystallization is obtained during the course of the simulation (this is expected as homogeneous nucleation is a rare event on MD time and length scales, and it requires large-scale

Figure 7. Liquid densities of 1,1,1,3,3-pentafluorobutane (PFB), 1,1,1,3,3,3-hexafluoropropane (HFP), 2,2,4,4-tetrafluoropentane (TFP), and 1,1,1-trifluoroethane (TFE). Symbols correspond to the simulated data obtained with our new force field, and lines represent experimental data (except for TFP where a single experimental point is shown). PFB at 289.15 K (black filled circles and black solid line),19 PFB at 395.15 K (black open circles and black dashed line),19 HFP at 293.38 K (blue filled squares and blue dotted−dashed line),60 HFP at 373.17 K (blue open squares and blue dotted line),60 TFP at 293.15 K (green open triangles and green solid line),61 and TFE at 283.15 K (green open diamonds and green dashed−double-dotted line).59

Figure 8. Pure PVDF density versus temperature at atmospheric pressure. Experimental data taken from Zöller and Walsh67 (red filled circles) and molecular dynamics results with our new force field (blue filled diamonds). Dashed line represents experimental densities in the liquid phase extrapolated to temperatures below the PVDF melting temperature. 146



Article

Table 6. Intermolecular Potential Parameters for the Different Studied Penetrantsa

a

molecule

force center

ϵ (K)

σ (Å)

q (e)

CO2 CO2 O2 O2 O2 N2 N2 H2 H2 CH4 H2S H2S H2S H2O H2O H2O

C O O q− CoM N CoM CoM H CoM S H q− O H q−

28.129 80.507 43.183

2.7570 3.0330 3.1062

+0.6520 −0.3260

36.000

3.3000

36.700

2.9580

149.920 250.000

3.7330 3.7300

93.200

3.1589

−2.1000 +4.2000 −0.5075 +1.0150 −0.9328 +0.4664 +0.4000 +0.2500 −0.9000 +0.5564 −1.1128

geometry lC−O = 1.1490 Å lO−O = 0.9700 Å lq−−q− = 0.2000 Å lN−N = 1.0980 Å lH−H = 0.7900 Å

lS−H = 1.3400 Å θH−S−H = 92.00° lS−q− = 0.1862 Å lO−H = 0.9572 Å θH−O−H = 104.52° lO−q− = 0.1546 Å

ref 68 68 73 73 73 71 71 69 69 74 70 70 70 75 75 75

CoM means the center of mass of the molecule.

intermediate axial position). The interaction parameters for these different force fields are given in Table 6. 4.2.2. Henry Constants of Different Penetrants in Melted PVDF. Infinite dilution calculations have been performed at atmospheric pressure and 493 K, i.e., above the melting point of PVDF. For some of the studied penetrants, experimental solubility data76−78 are available in the literature for temperatures below the melting temperature, but no data have been found for melted PVDF. The logarithm of the reciprocal of the Henry constant of the different studied penetrants is plotted in Figure 9 as a function

simulations and/or coarse-grained models to be observed.64−66) The density linearly increases with decreasing temperature, in good agreement with the experimental data67 for temperatures above the PVDF melting point. Simulated densities are in good agreement with the extrapolated values of the experimental liquid densities below the melting temperature. If we account for chain length effect on density, the predicted values of our new model for large polymer chains are overestimated by 1−2%, a result consistent with the one obtained for short HFCs. 4.2. Infinite Dilution Calculations. 4.2.1. Penetrant Force Fields. The EPM2 model from Harris and Yung68 has been used for CO2. The molecule is treated as a rigid linear body, with C−O distance equal to 1.149 Å. Three LennardJones sites and three partial charges located on each atom are used to represent dispersion−repulsion and electrostatic interactions. We used the Darkrim model69 to represent the hydrogen molecule. In this model, a single Lennard-Jones center is placed at the center of mass of the molecule. A quadrupolar moment is created by placing two positive charges on each hydrogen atom and one negative charge on the molecular center of mass. The hydrogen sulfide molecule has been described by the Kristóf and Liszi UA model,70 which involves a single Lennard-Jones site and four electrostatic charges. The force field chosen for nitrogen was developed by Delhommelle.71 Using this definition, the N2 molecule consists of two Lennard-Jones centers separated by 1.098 Å. It involves two negative electrostatic charges on nitrogen atoms and a positive one located on the center of mass of the molecule. For oxygen, we have used an adaptation of the force field developed by Vrabec et al.72 proposed by Boutard et al.73 This model, which avoids the explicit consideration of a quadrupole point, consists of two Lennard-Jones centers spaced by 0.970 Å and three electrostatic charges: a positive charge in the center of mass and two negative charges 0.200 Å apart from the center of mass on the O−O axis. Methane is described as a single Lennard-Jones site.74 Finally, the TIP4P/2005 force field75 was chosen to model water molecules. This force field involves one Lennard-Jones center and three electrostatic charges (two positive charges located on the hydrogen atoms and one negative charge that is shifted from the oxygen nucleus to an

Figure 9. Logarithm of the reciprocal of Henry constant of different penetrants in PVDF at 493 K as a function of the penetrant critical temperature. Dashed line is a guide for the eyes.

of the penetrant critical temperature. A linear behavior is observed, in agreement with the observations reported several years ago by Durrill and Griskey79 for different polymeric materials at 461 K. This finding demonstrates that the penetrant solubilities in the Henry regime are governed primarily by the thermodynamic properties of the penetrant. 4.3. CO2 and CH4 Solubilities in Melted PVDF beyond the Henry Regime. In this section, we present solubility data of pure CH4, pure CO2, and CO2 + CH4 mixtures in melted PVDF calculated at 493 K using hydrid Monte Carlo simulations in the osmotic ensemble. Results for pure gases, expressed as weight concentrations, are presented in Figure 10 for pressures ranging from 1 to 10 MPa. As can be seen, the gas concentration in PVDF is linearly proportional to the pressure in this pressure range, following the Henry straight line. CO2 is 147



Article

Figure 12. Gas solubility coefficients versus gas fugacity for pure CH4 (red filled triangles) and CO2 (blue filled circles) and for CH4 (red open triangles) and CO2 (blue open circles) in equimolar mixtures at 493 K.

Figure 10. CO2 (blue open circles) and CH4 (red open triangles) concentration in melted PVDF at 493 K versus gas pressure. The solid lines represent infinite dilution calculations whereas symbols correspond to calculated values in the osmotic ensemble.

much more soluble in PVDF than CH4 with weight concentrations 10 times higher. No experimental data have been found in the literature to directly validate our results on melted PVDF, but the behavior observed in our simulations is consistent with experimental results reported for semicrystalline PVDF.76−78 An interesting feature of the osmotic ensemble is its ability to give, in addition to the amount of absorbed penetrant, the associated swelling degree of the polymer phase, Q defined as Q=

Vswell − Vpure Vpure

Si(P , yi ) =

Ci fi

(14)

where Ci is the penetrant concentration of species i, f i is the fugacity of species i, P is the total gas pressure, and yi is the molar composition of the gas phase. In this plot, partial gas solubility coefficients in the mixture are compared to pure gas solubility coefficients. An ideal behavior is observed, with identical gas solubility coefficients in the mixture or in pure systems. The solubility of the CH4 + CO2 mixture can thus be expressed from pure gas solubilities, at least under the thermodynamic conditions involved here.

× 100 (13)

where Vpure is the volume of the pure polymer at the studied temperature and atmospheric pressure and Vswell is the volume of the swollen polymer. Vswell is directly obtained from simulations in the osmotic ensemble, whereas Vpure is obtained from simple NpT simulations of the pure polymer. The swelling degree of the studied systems are given in Figure 11. As

5. CONCLUSION We have proposed a new force field to model hydrofluorocarbons. This force field is based on the standard 12−6 Lennard-Jones potential which makes it easy to use in combination with common existing force fields for organic compounds. This force field has been validated using the few available experimental data for HFCs in the literature. We believe that the good agreement between our simulation results and these experimental data makes this force field suitable for quantitative predictions of thermophysical properties of other HFC compounds with different chain lengths. This force field has been applied to a long HFC chain, taken as a model of PVDF polymer. The volumetric properties of pure PVDF are recovered at temperatures above the polymer melting point. Below this melting temperature, simulated densities are in good agreement with extrapolated values of the experimental liquid density. Gas + PVDF mixtures have also been studied in order to assess infinite dilution solubilities. From the study of seven different gases, we demonstrate that the Henry constant can be related to the penetrant critical temperature, in agreement with the behavior commonly observed in other polymers. Work is now in progress to study gas solubility in PVDF below its melting temperature. The objective of this ongoing work is to try to recover experimental solubilities of gases in semicrystalline PVDF by performing MC simulations of the polymer amorphous phase in the osmotic ensemble using an ad hoc constraint to mimic the overall effect of the crystalline regions. Such an approach has already been successfully used to quantitatively calculate solubility data of many different gases in semicrystalline polyethylene in previous works.46,47,80

Figure 11. PVDF swelling associated with CO2 (blue open circles) and CH4 (red open triangles) sorption at 493 K and different gas pressures, calculated using hybrid MC simulations in the osmotic ensemble.

expected from the solubility values shown previously, a high swelling degree is observed for CO2 (Q = 5.1% at 10 MPa), while almost no swelling is recorded for CH4 (Q = 0.2% at 10 MPa). We end this study by the presentation of the solubility of an equimolar CH4 + CO2 mixture in PVDF at an imposed temperature of 493 K and a total gas pressure ranging from 1 to 20 MPa. The obtained results, shown in Figure 12, are expressed as partial solubility coefficients, Si(P,yi): 148



■

Article

(15) Karasawa, N.; Goddard, W. I. Dielectric Properties of Poly(vinylidene fluoride) from Molecular Dynamics Simulations. Macromolecules 1995, 28, 6765−6772. (16) Carbeck, J.; Rutledge, G. A Simple Approach to Electronic Charge Induction in Atomistic Simulations. I. Model Development and Parametrization. J. Comput.-Aided Mater. Des. 1993, 1, 97−110. (17) Carbeck, J.; Lacks, D.; Rutledge, G. A Model of Crystal Polarization in β-Poly(vinylidene fluoride). J. Chem. Phys. 1995, 103, 10347−10355. (18) Carbeck, J.; Rutledge, G. Temperature Dependent Elastic, Piezoelectric and Pyroelectric Properties of β-Poly(vinylidene fluoride) from Molecular Simulation. Polymer 1996, 37, 5089−5097. (19) Marrucho, I.; Oliveira, N.; Dohrn, R. Vapor-Phase Thermal Conductivity, Vapor Pressure, and Liquid Density of R365mfc. J. Chem. Eng. Data 2002, 47, 554−558. (20) Bobbo, S.; Scattolini, M.; Fedele, L.; Camporese, R. Compressed Liquid Densities and Saturated Liquid Densities of HFC365mfc. Fluid Phase Equilib. 2004, 222−223, 291−296. (21) Ribeiro, A. P. C.; Nieto De Castro, C. A.; Pai-Panandiker, R. S.; Mardolcar, U. V. Relative Permittivities of 1,1,1,2,3,3,3-Heptafluoropropane (HFC-227ea), 1,1,1,2,3,3-Hexafluoropropane (HFC-236ea), and 1,1,1,3,3-Pentafluorobutane (HFC-365mfc) in the Liquid Phase. J. Chem. Eng. Data 2007, 52, 2041−2049. (22) Fröba, A. P.; Krzeminski, K.; Leipertz, A. Thermophysical Properties of 1,1,1,3,3-Pentafluorobutane (R365mfc). Int. J. Thermophys. 2004, 25, 987−1004. (23) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Oxford Science: Oxford, U.K., 1989. (24) Martyna, G. J.; Tuckerman, M.; Tobias, D.; Klein, M. Explicit Reversible Integrators for Extended Systems Dynamics. Mol. Phys. 1996, 87, 1117−1157. (25) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Dinola, A.; Haak, J. R. Molecular-Dynamics with Coupling to an External Bath. J. Chem. Phys. 1984, 81, 3684−3690. (26) Bourasseau, E.; Ungerer, P.; Boutin, A.; Fuchs, A. H. Monte Carlo Simulation of Branched Alkanes and Long Chain n-Alkanes with Anisotropic United Atoms Intermolecular Potential. Mol. Simul. 2002, 28, 317−336. (27) Panagiotopoulos, A. Z. Direct Determination of Phase Coexistence Properties of Fluids by Monte Carlo Simulation in a New Ensemble. Mol. Phys. 1987, 61, 813−826. (28) Widom, B. Some Topics in the Theory of Fluids. J. Chem. Phys. 1963, 39, 2808−2812. (29) Escobedo, F. A. Novel Pseudoensembles for Simulation of Multicomponent Phase Equilibria. J. Chem. Phys. 1998, 108, 8761− 8772. (30) Banaszak, B. J.; Faller, R.; de Pablo, J. J. Simulation of the Effects of Chain Architecture on the Sorption of Ethylene in Polyethylene. J. Chem. Phys. 2004, 120, 11304−11315. (31) Dodd, L. R.; Boone, T. D.; Theodorou, D. N. A Concerted Rotation Algorithm for Atomistic Monte Carlo Simulation of Polymer Melts and Glasses. Mol. Phys. 1993, 78, 961−996. (32) Boyd, R. H. An Off-Lattice Constant-Pressure Simulation of Liquid Polymethylene. Macromolecules 1989, 22, 2477−2481. (33) Mavrantzas, V. G.; Boone, T. D.; Zervopoulou, E.; Theodorou, D. N. End-Bridging Monte Carlo: A Fast Algorithm for Atomistic Simulation of Condensed Phases of Long Polymer Chains. Macromolecules 1999, 32, 5072−5096. (34) Karayiannis, N.; Mavrantzas, V. G.; Theodorou, D. N. A Novel Monte Carlo Scheme for the Rapid Equilibration of Atomistic Model Polymer Systems of Precisely Defined Molecular Architecture. Phys. Rev. Lett. 2002, 88, 1−4. (35) Karayiannis, N. C.; Giannousaki, A. E.; Mavrantzas, V. G.; Theodorou, D. N. Atomistic Monte Carlo Simulation of Strictly Monodisperse Long Polyethylene Melts through a Generalized Chain Bridging Algorithm. J. Chem. Phys. 2002, 117, 5465−5479. (36) Theodorou, D. N. Understanding and Predicting StructureProperty Relations in Polymeric Materials through Molecular Simulations. Mol. Phys. 2004, 102, 147−166.

ASSOCIATED CONTENT

S Supporting Information *

Text describing conversion between different analytical expressions for torsion potential energy, figures showing PFB liquid densities at 393.15 K and DFP liquid−vapor equilibrium properties, and table listing PFB liquid−vapor equilibrium properties. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We gratefully acknowledge P. Archirel for assistance performing quantum calculations. B.R. thanks the PREDIMOL project (Grant ANR-10-CDII-007), financed by ANR and the French ministry in charge of environment.

■

REFERENCES

(1) Theodorou, D. N. In Materials Science of Membranes; Yampolskii, Y., Pinnau, I., Freeman, B. D., Eds.; John Wiley & Sons: Chichester, U.K., 2006; pp 49−94. (2) Cui, S.; Siepmann, J.; Cochran, H.; Cummings, P. Intermolecular Potentials and Vapor-Liquid Phase Equilibria of Perfluorinated Alkanes. Fluid Phase Equilib. 1998, 146, 51−61. (3) Cui, S.; Cochran, H.; Cummings, P. Vapor-Liquid Phase Coexistence of Alkane-Carbon Dioxide and Perfluoroalkane-Carbon Dioxide Mixtures. J. Phys. Chem. B 1999, 103, 4485−4491. (4) Chen, K.-H.; Walker, G.; Allinger, N. A Molecular Mechanics (MM3) Study of Fluorinated Hydrocarbons. J. Mol. Struct. (THEOCHEM) 1999, 490, 87−107. (5) Chen, K.-H.; Lii, J.-H.; Walker, G.; Y, X.; Schaefer, H., III; Allinger, N. Molecular Mechanics (MM4) of Fluorinated Hydrocarbons. J. Phys. Chem. A 2006, 110, 7202−7227. (6) Song, W.; Rossky, P.; Maroncelli, M. Modeling Alkane +Perfluoroalkane Interactions using All-Atom Potentials: Failure of the Usual Combining Rules. J. Chem. Phys. 2003, 119, 9145−9162. (7) Watkins, E.; Jorgensen, W. Perfluoroalkanes: Conformational Analysis and Liquid-State Properties from ab-Initio and Monte Carlo Calculations. J. Phys. Chem. A 2001, 105, 4118−4125. (8) Potoff, J. J.; Bernard-Brunel, D. A. Mie Potentials for Phase Equilibria Calculations: Application to Alkanes and Perfluoroalkanes. J. Phys. Chem. B 2009, 113, 14725−14731. (9) Du, Q.; Yang, Z.; Yang, N.; Yang, X. Coarse-Grained Model for Perfluorocarbons and Phase Equilibrium Simulation of Perfluorocarbons/CO2 Mixtures. Ind. Eng. Chem. Res. 2010, 49, 8271−8278. (10) Paulechka, E.; Kroenlein, K.; Kazakov, A.; Frenkel, M. A Systematic Approach for Development of an OPLS-Like Force Field and Its Application to Hydrofluorocarbons. J. Phys. Chem. B 2012, 116, 14389−14397. (11) Hariharan, A.; Harris, J. Structure and Thermodynamics of the Liquid−Vapor Interface of Fluorocarbons and Semifluorinated Alkane Diblocks: A Molecular Dynamics Study. J. Chem. Phys. 1994, 101, 4156−4165. (12) Shin, S.; Collazo, N.; Rice, S. A Molecular Dynamics Study of the Packing Structures in Monolayers of Partially Fluorinated Amphiphiles. J. Chem. Phys. 1992, 96, 1352−1366. (13) Byutner, O.; Smith, G. Quantum Chemistry Based Force Field for Simulations of Poly(vinylidene fluoride). Macromolecules 2000, 33, 4264−4270. (14) Karasawa, N.; Goddard, W. I. Force Fields, Structures, and Properties of Poly(vinylidene fluoride) Crystals. Macromolecules 1992, 25, 7268−7281. 149



Article

(37) Duane, S.; Kennedy, A. D.; Pendleton, B. J.; Roweth, D. Hybrid Monte Carlo. Phys. Lett. B 1987, 195, 216−222. (38) Mehlig, B.; Heermann, D.; Forrest, B. Hybrid Monte Carlo Method for Condensed-Matter Systems. Phys. Rev. B 1992, 45, 679− 685. (39) Faller, R.; de Pablo, J. J. Constant Pressure Hybrid Molecular DynamicsMonte Carlo Simulations. J. Chem. Phys. 2002, 116, 55− 59. (40) Faller, R.; de Pablo, J. J. Constant Pressure Hybrid Molecular Dynamics-Monte Carlo Simulations Erratum. J. Chem. Phys. 2003, 119, 7605−7605. (41) Nguyen, T. V.-O.; Houriez, C.; Rousseau, B. Viscosity of the 1Ethyl-3-methylimidazolium Bis(trifluoromethylsulfonyl)imide Ionic Liquid from Equilibrium and Nonequilibrium Molecular Dynamics. Phys. Chem. Chem. Phys. 2010, 12, 930−936. (42) Ungerer, P.; Tavitian, B.; Boutin, A. Applications of Molecular Simulation in the Oil and Gas Industry: Monte Carlo Methods, Editions Technip; IFP: Paris, France, 2005. (43) Kong, C. Combining Rules for Intermolecular Potential Parameters. II. Rules for the Lennard-Jones (12−6) Potential and the Morse Potential. J. Chem. Phys. 1973, 59, 2464−2467. (44) Delhommelle, J.; Evans, D. J. Configurational Temperature Profile in Confined Fluids. I. Atomic Fluid. J. Chem. Phys. 2001, 114, 6229−6235. (45) Faure, F.; Rousseau, B.; Lachet, V.; Ungerer, P. Molecular Simulation of the Solubility and Diffusion of Carbon Dioxide and Hydrogensulfide in Polyethylene Melts. Fluid Phase Equilib. 2007, 261, 168−175. (46) Memari, P.; Lachet, V.; Rousseau, B. Molecular Simulations of the Solubility of Gases in Polyethylene below Its Melting Temperature. Polymer 2010, 51, 4978−4984. (47) Memari, P.; Lachet, V.; Klopffer, M.-H.; Flaconnèche, B.; Rousseau, B. Gas Mixture Solubilities in Polyethylene below Its Melting Temperature: Experimental and Molecular Simulation Studies. J. Membr. Sci. 2012, 390−391, 194−200. (48) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. A 2nd Generation Force-Field for the Simulation of Proteins, Nucleic-Acids, and Organic-Molecules. J. Am. Chem. Soc. 1995, 117, 5179−5197. (49) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (50) Scott, W. R. P.; Hünenberger, P. H.; Tironi, I. G.; Mark, A. E.; Billeter, S. R.; Fennen, J.; Torda, A. E.; Huber, T.; Krüger, P.; van Gunsteren, W. F. The GROMOS Biomolecular Simulation Program Package. J. Phys. Chem. A 1999, 103, 3596−3607. (51) Bourasseau, E.; Haboudou, M.; Boutin, A.; Fuchs, A. H.; Ungerer, P. New Optimization Method for Intermolecular Potentials: Optimization of a New Anisotropic United Atoms Potential for Olefins: Prediction of Equilibrium Properties. J. Chem. Phys. 2003, 118, No. 3020. (52) Lemmon, E. W.; Huber, M.; McLinden, M. NIST Standard Reference Database 23: Reference Fluid Thermodynamic and Transport Properties-REFPROP, Version 8.0; National Institute of Standards and Technology, Standard Reference Data Program: Gaithersburg, MD, USA, 2007. (53) NIST, NIST/TRC Web Thermo Tables (WTT), NIST Standard Reference Subscription Database 3−Professional Edition, Version 2-2012-1-Pro. (54) Neumann, M. Dipole Moment Fluctuation Formulas in Computer Simulations of Polar Systems. Mol. Phys. 1983, 50, 841− 858. (55) Smith, P. E. Constant Dielectric Properties of the Simple Point Charge and Extended Simple Point Charge Water Models at 277 and 300 K. J. Chem. Phys. 1994, 100, 3169−3174.

(56) Henne, A.; Hinkamp, J. The Synthesis and Directed Chlorination of 2,2-Difluorobutane. J. Am. Chem. Soc. 1945, 67, 1194−1197. (57) Mondello, M.; Grest, G. S. Viscosity Calculations of n-Alkanes by Equilibrium Molecular Dynamics. J. Chem. Phys. 1997, 106, 9327− 9336. (58) Dysthe, D. K.; Fuchs, A. H.; Rousseau, B. Fluid Transport Properties by Equilibrium Molecular Dynamics III. Evaluation of United Atom Interaction Potential Models for Pure Alkanes. J. Chem. Phys. 2000, 112, 7581−7590. (59) Nakamura, S.; Fujiwara, K.; Noguchi, M. PVT Properties for 1,1,1-Trifluoroethane (R-143a). J. Chem. Eng. Data 1997, 9568, 334− 338. (60) Yin, J.; Wu, J. Compressed Liquid Densities of 1,1,1,3,3Pentafluoropropane (HFC-245fa) and 1,1,1,3,3,3-Hexafluoropropane (HFC-236fa). Fluid Phase Equilib. 2011, 307, 1−5. (61) Stepanov, I.; Burmakov, A.; Kunshenko, B.; Alekseeva, L.; Yagupol’skii, L. Reaction of Hydroxy and Carbonyl Compounds with Sulfur-Tetrafluoride. IX. Reactions of Aliphatic β-Diketones with SF4. Zh. Org. Khim. 1983, 19, 273−279. (62) Henne, A. L.; Renoll, M. W.; Leicester, H. M. Aliphatic Difluorides. J. Am. Chem. Soc. 1939, 61, 938−940. (63) Laesecke, A.; Defibaugh, D. R. Viscosity of 1,1,1,2,3,3Hexafluoropropane and 1,1,1,3,3,3-Hexafluoropropane at SaturatedLiquid Conditions from 262 to 353 K. J. Chem. Eng. Data 1996, 41, 59−62. (64) Meyer, H.; Müller-Plathe, F. Formation of Chain-Folded Structures in Supercooled Polymer Melts. J. Chem. Phys. 2001, 115, 7807. (65) Gee, R. H.; Lacevic, N.; Fried, L. E. Atomistic Simulations of Spinodal Phase Separation Preceding Polymer Crystallization. Nat. Mater. 2006, 5, 39−43. (66) Sommer, J.-U.; Luo, C. Molecular Dynamics Simulations of Semicrystalline Polymers: Crystallization, Melting, and Reorganization. J. Polym. Sci., Part B: Polym. Phys. 2010, 48, 2222−2232. (67) Zoller, P.; Walsh, D. Standard Pressure−Volume−Temperature Data for Polymers, 1st ed.; Technomic: Lancaster, PA, USA, 1995. (68) Harris, J. G.; Yung, K. H. Carbon Dioxide’s Liquid-Vapor Coexistence Curve and Critical Properties as Predicted by a Simple Molecular Model. J. Phys. Chem. 1995, 99, 12021−12024. (69) Darkrim, F.; Vermesse, J.; Malbrunot, P.; Levesque, D. Monte Carlo Simulations of Nitrogen and Hydrogen Physisorption at High Pressures and Room Temperature. Comparison with Experiments. J. Chem. Phys. 1999, 110, 4020−4027. (70) Kristóf, T.; Liszi, J. Effective Intermolecular Potential for Fluid Hydrogen Sulfide. J. Phys. Chem. B 1997, 101, 5480−5483. (71) Delhommelle, J. Etablissement de Potentiels d’Interactions pour la Simulation Moléculaire. Application à la Prédiction des Equilibres Liquide-Vapeur de Mélanges Binaires Alcane-Molécule Multipolaire. Ph.D. thesis, Université Paris-Sud 11, Orsay, France, 2000. (72) Vrabec, J.; Stoll, J.; Hasse, H. A Set of Molecular Models for Symmetric Quadrupolar Fluids. J. Phys. Chem. B 2001, 105, 12126− 12133. (73) Boutard, B.; Ungerer, P.; Teuler, J.-M.; Ahunbay, M. G.; Sabater, S. F.; Perez-Pellitero, J.; Mackie, A. D.; Bourasseau, E. Extension of the Anisotropic United Atoms Intermolecular Potential to Amines, Amides and Alkanols. Application to the Problems of the 2004 Fluid Simulation Challenge. Fluid Phase Equilib. 2005, 236, 25−41. (74) Möller, D.; Ó przynski, J.; Müller, A.; Fischer, J. Prediction of Thermodynamic Properties of Fluid Mixtures by Molecular Dynamics Simulations: Methane-Ethane. Mol. Phys. 1992, 75, 363−378. (75) Abascal, J. L. F.; Vega, C. A General Purpose Model for the Condensed Phases of Water: TIP4P/2005. J. Chem. Phys. 2005, 123, 234505−234512. (76) Solms, N. V.; Zecchin, N.; Rubin, A.; Andersen, S. I.; Stenby, E. H. Direct Measurement of Gas Solubility and Diffusivity in Poly (Vinylidene Fluoride) with a High-Pressure Microbalance. Eur. Polym. J. 2005, 41, 341−348. 150



Article

(77) Flaconnèche, B.; Martin, J.; Klopffer, M.-H. Permeability, Diffusion and Solubility of Gases in Polyethylene, Polyamide 11 and Poly(Vinylidene Fluoride). Oil Gas Sci. Technol. 2001, 56, 261−278. (78) El-Hibri, M.; Paul, D. Gas Transport in Poly (Vinylidene Fluoride): Effects of Uniaxial Drawing and Processing Temperature. J. Appl. Polym. Sci. 1986, 31, 2533−2560. (79) Durrill, P.; Griskey, R. Diffusion and Solution of Gases into Thermally Soften or Molten Polymers. AIChE J. 1969, 15, 106−110. (80) Memari, P.; Lachet, V.; Rousseau, B. Gas Permeation in Semicrystalline Polyethylene as Studied by Molecular Simulation and Elastic Model. Oil Gas Sci. Technol. 2013, 1−9.

151


Classical Force Field for Hydrofluorocarbon Molecular Simulations

Recommend Documents