| Property Prediction |
, binding energy component analysis, and the identification of favorable binding configurations between molecular pairs, molecules and surfaces, additives and bulk materials, and liquids and crystals.
Sections in this chapter
Introduction
Blends provides a way to shorten the discovery process by estimating the miscibility behavior of binary mixtures. These include solvent-solvent, polymer-solvent, and polymer-polymer mixtures. Blends predicts the thermodynamics of mixing directly from the chemical structures of the two components and, therefore, requires only their molecular structures and a reasonable choice of nonbond forcefield parameters as input.
Hmix and TSmix).
Blends is targeted to the following industrial areas:
, calculate phase diagrams, and plot thermodynamic functions all start with the same sequence of steps:Calculate functions
1. Specify the blends molecules.
2. Specify the packing variables.
3. Calculate the pair interaction energies.
4. Calculate the coordination numbers.
5. Fit the energy of mixing data with an analytical model.
This data can then be used to:
(T).
Analyze functions
After performing the appropriate calculations, the following analytical functions can be performed:
Extracting pairs
|
Please see the online help for more detailed information on the control panels.
|
Both theoretical models and molecular simulation techniques have been employed in an effort to predict the thermodynamic behavior of binary mixtures. For more discussion of theoretical models, please see the "Theory" section on page 54.
General methodology
To calculate pair interaction energies
1. Place the two blend molecules to be used in two model spaces
(these can either be loaded from a file or built using one or more
of the builders).
| The blend molecules used to generate the molecular pairs should be individually and fully minimized prior to performing the calculations. |
2. Select Calculate on the BLENDS card to bring up the Run
Blends control panel.
3. Enter the model number for molecules 1 and 2 in the Model
entry boxes.
4. See "To specify the packing variables" on page 46.
Specify the calculation controls
6. Enter the number of molecular pairs to be generated.
8. Enter the file prefix for the energy files to be created.
9. Select the text output level from the popup (None, Default, or
Verbose).
11. To plot an energy distribution histogram as the pairs are generated:
a. Check the Update Graph box.
b. Enter the frequency (number of pairs between updates).
c. Enter the binwidth (kcal/mol).
12. Click the Calculate button for the Eij to be calculated (E11, E12,
E21, or E22).
To calculate all Eijs (E11, E12, E21, and E22)
13. Return to the Run Blends control panel and click the Calculate
Interaction Energies button.
See "Calculating pair interaction energies" on page 58 for theoretical background on this method.
To calculate the coordination numbers
1. Place the two blend molecules to be used in two model spaces
(these can either be loaded from a file or built using one of the
builders).
2. Select Calculate on the BLENDS card to bring up the Run
Blends control panel.
3. Enter the model number for molecules 1 and 2 in the entry
boxes under Model.
4. See "To specify the packing variables" on page 46.
Specify the calculation controls
6. Enter a value for the number of clusters to be generated.
7. Enter the number of trials per cluster.
8. To obtain a different set of clusters, change the number used as
the random seed initializer.
9. Select the text output level from the popup (None, Default, or
Verbose).
12. Click the Calculate button for the coordination number to be
calculated (Z11, Z12, Z21, or Z22).
To calculate all Zs (Z11, Z12, Z21, and Z22)
Return to the Run Blends control panel and click the Calculate Coordination Numbers button.
To specify the packing variables
1. Select the Calculate item on the BLENDS card to bring up the
Run Blends control panel.
2. Click the top Packing... button to bring up the Molecule 1 Packing
control panel.
4. Select Isotropic or Axial Distribution from the popup.
5. If using axial distribution, enter the allowed angular spread
(degrees).
6. If using excluded atom constraints:
a. Check the Use Non-contact Atom List box (otherwise, leave
it unchecked).
b. Select the atoms to be excluded.
c. Click the Add Selected Atoms to List action button.
d. To clear the noncontact atoms list, click the Clear Contact
List action button.
e. To show the noncontact atoms in the model window, click
the Show Contact List action button.
8. Repeat Steps 3 through 6 for molecule 2.
Please see "Specifying packing variables" on page 63 for the theoretical background on this method.
| The interaction energies for the two blend molecules must already have been calculated and saved in .enr files (see "Calculating pair interaction energies" on page 58). Also, if you are using calculated values for Z, these calculations must already have been performed (see "Calculating coordination numbers" on page 61). |
1. Bring up the Run Blends control panel.
3. To list additional data in the text window, check the Print Data
in Text Window box.
4. Select the algorithm to be used from those listed in the Emix (T)
Model popup.
Please see "Fitting the mixing energy and calculating Chi" on page 64 for more theoretical background on this method.
To plot the interaction parameter Chi(T)
1. Calculate the interaction energies for the two blend molecules
(see "Calculating pair interaction energies" on page 58).
2. Calculate the coordination numbers (see "Calculating coordination
numbers" on page 61), unless you are using specified
values for Z.
3. Fit the mixing energy data with an analytical model (see "To fit
the mixing energy data" on page 47).
5. Click the Plot Interaction Parameter Chi(T) button.
Please see "Fitting the mixing energy and calculating Chi" on page 64 for more theoretical background on this method.
To calculate a phase diagram
1. Calculate the pairwise interaction energies for the two components
of the binary mixture (see "Calculating pair interaction
energies" on page 58).
2. Calculate the coordination number Zij for each of the four pair
combinations (see "Calculating coordination numbers" on
page 61). This step can be skipped if using your own Z values.
3. Do an analytical fit of the energy of mixing data (see "To fit the
mixing energy data" on page 47).
6. Select the Emix (T) model to be used from the popup.
8. Enter the degree of polymerization to be used for component 1
(X1) and component 2 (X2).
Please see "Calculating phase diagrams" on page 65 for more theoretical information on this method.
Plotting pair energy distributions
The generation of different orientations using the Pairs Method leads to configurations of varying energetics. Often, very favorable energies are found; many actually locate a local minimum on the potential energy surface. The distribution of the pair energies obtained can be plotted using the options on the Blends Energy Analysis control panel (see the online help for more control panel information). The total energy of interaction can be plotted, or just the electrostatic, van der Waals, or hydrogen bond terms. Inspection of these distributions can be helpful in noting the different modes and frequencies with which two molecules interact.
To plot a pair energy distribution
The interaction energies for the pair must already have been calculated and saved in a .enr file (see "Calculating pair interaction energies" on page 58).
4. To normalize the histogram, check the Normalize Histogram
box.
6. Enter the binwidth (kcal/mol).
A plot of the pair energy distribution appears in the graph window.
Plotting thermodynamic functions
Blends analysis options can be used both to calculate thermodynamic functions (enthalpy, entropy, free energy of mixing) for a binary system, and to create plots of these functions versus composition at a specified temperature. The plots generated reflect the current choice of the interaction parameter model (T) and the degree of polymerization (X1, X2) of the two components. These thermodynamic analysis options are provided on the Thermodynamic Analysis control panel (see the online help for more control panel information).
Several temperatures (isotherms) for a given thermodynamic function can also be calculated and plotted using the options on the Blends Isotherms control panel (see "Plotting thermodynamic isotherms" on page 52).
To calculate and plot thermodynamic functions
1. Select Analyze on the BLENDS card, then choose Thermodynamics
to bring up the Thermodynamic Analysis control panel.
Specify the functions to be computed
2. To plot enthalpy, check the Plot Enthalpy of Mixing (dH) box.
3. To plot entropy, check the Plot Entropy of Mixing (-TdS) box.
4. To plot the Gibbs free energy of mixing, check the Plot Free
Energy of Mixing (dG) box.
5. Enter the Temperature (K) value to be used for the plots.
6. Click the Thermodynamic Model... button to bring up the
Thermodynamic Model control panel.
7. Select the Emix (T) model to be used from the popup.
8. Enter values for the model parameters (A, B, and C) or use
those calculated from a Blends energy of mixing fit (see "Fitting
the mixing energy and calculating Chi" on page 64).
Specify the degree of polymerization
9. Enter the degree of polymerization to be used for component 1
(X1) and component 2 (X2).
Plotting thermodynamic isotherms
Several temperatures (isotherms) for a given thermodynamic function can be calculated and plotted using the options on the Blends Isotherms control panel (see the online help). The functions that can be plotted are free energy of mixing, enthalpy, and entropy. You can specify the number of curves plotted and the range of temperatures used. The plots generated reflect the current choice of the interaction parameter model (T) and the degree of polymerization (X1, X2) of the two components. You can specify the model parameters directly, or Blends can calculate these based on a molecular simulation. For more details, see "Fitting the mixing energy and calculating Chi" on page 64.
To plot more than one thermodynamic mixing function on the same plot at a fixed temperature, use the options on the Thermodynamics Analysis control panel (see "Plotting thermodynamic functions" on page 51).
To plot thermodynamic isotherms
1. Select Analyze on the BLENDS card, then choose Isotherms to
bring up the Blends Isotherms control panel.
2. Select the thermodynamic function to be plotted (Free Energy,
Entropy, or Enthalpy).
3. Enter the number of isotherms to be plotted.
4. Enter the temperatures for the first and last isotherms in the __
K to __K entry boxes.
5. Click the Thermodynamic Model... button to bring up the
Thermodynamic Model control panel.
6. Select the Emix (T) model to be used from the popup.
7. Enter values for the model parameters (A, B, and C) or use
those calculated from a Blends energy of mixing fit (see "Fitting
the mixing energy and calculating Chi" on page 64).
Specify the degree of polymerization
8. Enter the degree of polymerization to be used for component 1
(X1) and component 2 (X2).
9. Return to the Blends Isotherms control panel, and click the Plot
Mixing Isotherms button.
Extracting molecular pairs
Molecular pair configurations can be extracted from blends interaction energy files and saved in .msi files so that they can be viewed later. The configurations are sorted according to energy (lowest to highest) using the energy term specified (total energy, van der Waals, electrostatic, or hydrogen bond energy). The extracted configurations can be restricted to a maximum number having energies within a specified range. Appropriate values can be determined by viewing the Eij distribution for the corresponding pair interaction type (see "Plotting pair energy distributions" on page 49). The options used to extract molecular pairs are found on the Extract Pairs control panel (see the online help).
To extract molecular pair configurations
1. Choose Extract Pairs on the BLENDS card to open the Extract
Pairs control panel.
3. Select the energy term to be used from the popup (Total, vdW,
Coulomb, or H-Bond).
4. Enter the lowest and highest energy values for the pairs to be
extracted (in kcal/mol).
5. Enter the maximum number of pairs to be extracted.
6. Click the Extract Molecular Pairs button.
Theory
Theoretical models
Flory-Huggins model
Undoubtedly, the simplest and best-known theory of the thermodynamics of mixing and phase separation in binary systems is the Flory-Huggins lattice theory (Flory 1953). The general expression for calculating the free energy of mixing G for a binary system is:
Eq. 1
Where:
G = Free energy of mixing per mole.
= Volume fraction. is defined as:
Eq. 2
Where:
E12 = Differential energy of interaction of an unlike pair.
Eq. 3
Here Eij is the energy of a particular ij pair.
Certain deficiencies have been noted in the Flory-Huggins model, and more sophisticated theories, such as the reference interaction site model (RISM) and lattice cluster theory (LCT), have recently been developed to treat the problems beyond Flory-Huggins original approximations (Schweizer and Curro 1989, Freed 1985). The applicability of these theories is limited, however, because detailed information on each component required by these theories is often absent. Many parameters characteristic of a binary mixture are thus obtained by fitting a theoretical model with some experimental data. Prediction of the thermodynamic behavior for a system that is not well known remains a difficult task.
Molecular simulations
Recent advances in molecular simulation techniques have improved this situation. Accurate forcefields can be obtained by defining parameters using structural and spectral data. Molecular simulations can be performed on well-characterized systems leading to a better fundamental understanding of atomic level interactions. This information can then be used to predict useful physical properties of systems that are less well characterized.
Blends approach
Blends combines a modified Flory-Huggins model and molecular simulation techniques to calculate the compatibility of binary mixtures. Two important extensions to the Flory-Huggins model are employed:
. This is accomplished by first carrying out a molecular simulation of the differential energy of binding between similar and dissimilar pairs, then temperature averaging the results using the Boltzmann method (see "Calculating pair interaction energies" on page 58).
The temperature dependence of
is fitted with one of several analytical models, (T), then used in the expression for the free energy of mixing:
Eq. 4
A popular temperature model is the Kamide expression:
Eq. 5
Here, the mixture-dependent parameters A, B, and C are fully determined through a Blends molecular simulation (see "Fitting the mixing energy and calculating Chi" on page 64).
Analytical expressions of this type allow the accurate determination of first, second, and higher-order derivatives of the free energy of mixing with respect to volume fraction as a function of temperature. The derivatives are used to locate the coexisting curve (binodal) and the stability curve (spinodal) in a two-phase diagram. Critical mixing temperatures and volume fractions can easily be read from such a phase diagram (see "Calculating phase diagrams" on page 65). Thermodynamic functions, such as enthalpy, entropy, and free energy of mixing, can also be plotted (see "Plotting thermodynamic functions" on page 51).
As with any molecular simulation, the results obtained depend on the accuracy of the forcefield. By default, the Universal forcefield is used. However, you can either load a different forcefield using the Open Force Field (OFF) module or specify different parameters using the Force Field Editor. For information about the OFF and the Force Field Editor, see the Cerius2 Simulation Tools book.
Calculating pair interaction energies
Sampling problems
The differential energy of interaction of an unlike pair (Eij) can be obtained in a straightforward manner simply by calculating the energies of the four different pairs as defined by Eq. 3 of the Flory-Huggins model (see page 55). However, a question arises as to whether the Eij calculated from simple pairwise interactions is sufficiently accurate to represent interactions in the actual condensed state. It is important to take into account and to properly weigh a large number of relative orientations of the two molecules. Simple energy minimization using several selected configurations is not representative of the interaction energy for binary mixtures exhibiting normal Boltzmann distributions.
Molecular dynamics represents an improved sampling technique. However, molecular dynamics has difficulties in sampling significantly different regions of orientation space. If the temperature used is low, the two molecules remain trapped in a local energy minimum. If the temperature is too high, the molecules drift apart because sufficient kinetic energy is present for them to escape local energy minima. The use of periodic boundary conditions reduces some of these problems at the cost of longer simulation times. In general, the use of molecular dynamics for molecular sampling of binding energies is rather time consuming because the equations of motion need to be integrated with short time steps. Many time steps are needed to generate significantly different molecular orientations. These problems are illustrated by the fact that simple molecular dynamics simulations have shown that an averaged
E12 is strongly dependent on the initial configurations chosen.
Blends Monte Carlo sampling technique
Blends uses Monte Carlo atomistic simulations both to generate thousands of different molecular orientations and to calculate their pair-interaction energies; this is called the pairs method. This approach generates energetically favorable configurations by employing a Monte Carlo technique that includes excluded-volume constraints and takes temperature effects into account.
The pairs method
The pairs method results in four Boltzmann-averaged Eij values. These Eij values can be used to calculate both Eij using Eq. 3 (see page 55) and Emix(T) using Eq. 7 (see page 65).
Proper structures for each of the molecules are constructed and optimized. The overall shape of each molecule is represented by its van der Waals surface.
7. The centers of mass of both molecules are positioned near the
origin of the Cartesian coordinate frame; the coordinates of
molecule 1 (white in the above illustration) are unchanged
throughout the calculations.
9. A vector (n) that points from the origin to the surface of a unit
sphere is randomly chosen.
11. The pairwise interaction energy (Eij) of this specific configuration
is calculated and stored.
Eq. 6
14. Steps 2 through 8 are repeated for the other pair interaction
types, and four Boltzmann-averaged Eij(T) values are obtained:
Value
Molecule 1
Molecule 2
E11(T)
Component 1
Component 1
E12(T)
Component 1
Component 2
E21(T)
Component 2
Component 1
E22(T)
Component 2
Component 2
In the pair-generation method illustrated here, the orientation of molecule 2 is determined randomly. However, this may not be appropriate for oriented polymers or rigid-rod molecules such as liquid crystals. Connectivity between polymer segments should also be considered. As a result, Blends provides options that place restrictions on both molecule alignment and atom contacts during packing and, thus, allow you to obtain more representative Eij values (see "Specifying packing variables" on page 63).
By default, the mean interaction energies for every 100 pairs generated are reported in the text window along with the total binding energy, the total nonbond energy, and the energies due to van der Waals, electrostatic, and hydrogen bond interactions. Each of the pairs generated is displayed in the model window, and an energy distribution histogram is plotted in the graph window. The results of the calculations are saved in interaction energy (.enr) files. More detailed analysis of the pair-energy distributions can be obtained using the energy analysis functions (see "Plotting pair energy distributions" on page 49).
The number of pairs calculated, the interaction energy file name, and other output variables are specified on the Interaction Energies control panel. By default, 10,000 pairs are calculated, but a higher number is recommended for large molecules, molecules capable of hydrogen bonding, and polar molecules. The packing variables are specified on the Molecule 1 Packing and Molecule 2 Packing control panels (see the online help for more control panel information).
Calculating coordination numbers
The coordination number Zij is the number of molecules of type j that can be packed around a single molecule of type i. A single coordination number has a definite physical significance only when the two components of the binary mixture have similar volumes or surface areas. The difficulty in defining a coordination number for a system in which two components are not similar in size somewhat limits the applicability of the pairs method. For a binary system, at least four different combinations are possible; that is, the central molecule, as well as the surrounding molecule, can be either component 1 or component 2. This leads to four coordination numbers.
| Number | Central molecule | Surrounding molecule |
|---|---|---|
| Z11 | Component 1 | Component 1 |
| Z12 | Component 1 | Component 2 |
| Z21 | Component 2 | Component 1 |
| Z22 | Component 2 | Component 2 |
The Blends module differs from the Flory-Huggins model in that it uses an off-lattice calculation; that is, molecules are not arranged on a regular lattice as in the original Flory-Huggins theory. The coordination number Z is explicitly calculated for each of the possible molecular pairs using molecular simulations.
packing
The technique involves generating clusters in which nearest neighbors are packed around the central molecule until no more will fit.
The Z obtained depends on the number of packing trials allowed. In general, after a certain point, the number of nearest neighbors increases very slowly with an increase in the number of packing trials. The default maximum number of trials is 100.
Packing variables
In the method illustrated here, the orientation of the surrounding molecules is determined randomly. However, this may not be appropriate for oriented polymers or rigid-rod molecules such as liquid crystals. Connectivity between polymer segments should also be considered. Blends, therefore, provides options that allow you to place restrictions on molecule alignment and atom contacts during packing and, thus, obtain more representative values for Z (see "Specifying packing variables" on page 63).
The Z values are calculated for several clusters and an average Zij is obtained for each of the four pair combinations. Individual and average values are reported in the text window, and the running averages are plotted in the graph window. The plots show how well each Zij converged. By default, 100 clusters are generated for each Zij calculation.
The number of clusters generated, the number of trials per cluster, and the output variables are specified using the options on the
Z Number Calculation control panel. The packing variables are specified on the Molecule 1 Packing and Molecule 2 Packing control panels (see the online help for more control panel information).
calculate
Averaged coordination numbers are employed in the expression for the temperature dependence of the interaction parameter
:
Eq. 7
The Zij values used to calculate 
are shown on the Fit Mixing Energy control panel (see the online help). These values are updated when a coordination number calculation is performed. If you do not run a coordination number calculation, default values are employed in the (T) expression. Alternatively, you can specify the Z values to be used.
Specifying packing variables
Packing can be affected by placing restrictions on molecule alignment and atom contacts. The packing variables are set on the Molecule 1 Packing and Molecule 2 Packing control panels (see the online help).
Isotropic versus axial packing
In the pairs method and the Z calculation method described above, pair or cluster generation was illustrated using isotropic packing; that is, the orientation of the surrounding molecules with respect to the central molecule was determined randomly; all orientations were possible. However, isotropic packing may not be appropriate for oriented polymers or rigid-rod molecules, such as liquid crystals.
different from zero leads to orientations in the range of - to +. Packing occurs randomly within this range. Axial packing can be used in conjunction with another option that aligns the molecules along their principal axes; the long axis is reoriented along the z-direction.
Excluded atom constraints
For some systems it is important to consider chain connectivity and accessibility in packing. For example, in long-chain polymers, the ends of a polymer segment are normally connected to other segments, making some positions of the polymer segment inaccessible to other molecules. Blends allows you to exclude certain atoms from coming into contact with other atoms during packing. Configurations that place surrounding molecules in contact with an excluded atom are rejected.
Fitting the mixing energy and calculating Chi
Fitting Emix (T)
The interaction energy data contained in the interaction energy files can be plotted as a function of temperature, and fitted with a selected curve-fitting model. The energy of mixing is defined from the interaction energies and coordination numbers as follows:
Eq. 8
Boltzmann averaging the Eij values at each of the requested temperatures. For example, the values at each temperature T for Emix could be fitted using a model of the form A + BT + CT ln T, which, when divided by RT, gives the Kamide expression for the interaction parameter (T). Five curve-fitting models are provided.
The scatter plot and curve obtained from the fit are displayed in the graph window. Values for the fitted model parameters (A, B, and C) and the standard deviation for the curve are calculated and displayed in the text window. The model that gives the lowest standard deviation provides the best results.
(T)
The interaction parameter
can then be plotted as a function of temperature from the Emix (T) model data, using the fitted parameter values for A, B, and C:
Eq. 9
Calculating phase diagrams
The compatibility of binary mixtures can be illustrated by generating phase diagrams. These diagrams are obtained by calculating the free energy of mixing (G) as a function of composition at different temperatures. G (Gibbs free energy of mixing) versus 2 (volume fraction of the second component) has only one minimum. Two components are miscible for any composition. For temperatures below Tcr, the two points in contact with the straight line define two binodal points, A and D (below Tcr are two minima with ¹G/¹ = 0). Two inflection points with ¹2G/¹2 = 0, define spinodal points B and C. In the region of 2 < A and 2 > D, the two components are miscible. In the region between B and C, the system is unstable, separating into two phases with compositions equal to A and D. At the critical temperature, Tcr, A-D merge into a single point defined by ¹2G/¹2 = ¹3G/¹3 = 0.G
The Blends method for calculating G is based on a modified Flory-Huggins model (see "General methodology" on page 44).
The phase diagrams generated reflect the current choice of the
Emix (T) model, the degree of polymerization (X1, X2) of the two components, and the temperature points. These variables are specified on the Phase Diagram control panel (see the online help).
Bawendi, M. G.; Freed, K. F.; Mohanty, U. J.; Chem. Phys., 84, 7036 (1986)
References