Geochimica et Cosmochimica Acta, 1999
Condensation in Dust-enriched Systems, by D.S. Ebel and L. Grossman,
Geochimica et Cosmochimica Acta, 1999
The condensation code described here, "VAPORS", is described more completely by Ebel et al. (1999). All calculations are normalized to a total of one mole of atoms in the system.
A typical condensation run at fixed Ptot and bulk composition is begun with only the vapor phase present at 2400K. Most solutions
are obtained at 10K intervals, using the result at the previous temperature as a first approximation. At each fixed pressure, temperature and bulk composition of the system, the partial
pressures of the pure monatomic gaseous elements (the basis components of the gas phase) are obtained by calculating the distribution of the elements among 374 species in the gas phase,
using standard techniques (Lattimer et al., 1978; Smith and Missen, 1982). The stability of each potential, stoichiometrically pure, single component condensate phase is then evaluated
from the partial pressures of the elements and the Gibbs energy for that phase, by considering the energy balance of the formation reaction of the condensate from the monatomic gaseous elements.
In the case of a liquid or solid solution phase, the "best" composition is determined by finding that composition at which the activities of the components describing the solution phase most
closely match equivalent activities in the gas, using the algorithms of Ghiorso (1994). This composition is then tested for stability in much the same way as a stoichiometrically pure condensate,
but also accounting for the thermodynamic mixing properties of the components in the solution phase. In some cases where silicate liquid is present, this algorithm failed to find the "best"
pyroxene solid solution composition, and the program proceeded with a pure diopside end-member composition instead. In such cases, the program was restarted with a "seed", Ti, Al-bearing
diopsidic pyroxene substituted for the pure diopside at and above the temperature step at which pure diopside had been found to be stable. In all such cases, a complex, Ti-, Al-bearing diopsidic
pyroxene was found to be stable at the temperature where pure diopside had been found, or, at most, 20K higher. The Gibbs energy of the system was always on the order of 0.5 J lower per mole of
elements in the system for the assemblage with the pyroxene solid solution than for the one with pure diopside. This problem occurs only with the pyroxene solid solution model, probably because
of the difficulty in determining both the composition and ordering state of the near end-member pyroxene in equilibrium with a gas phase highly depleted in some of the pyroxene-forming elements.
Once a phase is determined to be stable, it is added to the stable assemblage in a seed amount (10-7 moles), which is subtracted from the gas.
The next step is to distribute mass between the phases to minimize the total free energy of this new system.
In this work, the second order technique of Ghiorso (1985), following Betts (1980), was adapted to the problem of distributing mass between phases to minimize directly the total Gibbs free energy of
a system consisting of gas and multiple pure and solution phases, both solid and liquid. The Gibbs energy of the entire system can be imagined as a surface in m dimensions, where m is the
total number of components independently variable in each of the phases present. The components of the gas are the monatomic elements, while those of solution phases are the end-members of these phases.
Each distribution of elements between these components at fixed temperature and pressure defines a state of the system, and corresponds to a point on the Gibbs surface. In successive iterations, information
about the local slope and curvature of the Gibbs surface at the current state of the system is used to determine the direction toward a minimum on this surface, along which the next iterative solution must lie.
Then atoms are redistributed accordingly among the gas and condensates, that is among the m components, so that this minimum is approached as closely as possible. From the perspective of this new state of
the system, the Gibbs surface "looks" different, so a new minimum must be sought in a further iteration. Convergence is declared when the vector norm of all the changes in composition in the m directions
does not change by > 10-12 between iterations. The VAPORS program usually converges in less than ten iterations in this part of the algorithm.
Upon convergence to a free energy minimum, the stabilities of non-condensed phases are assessed as described above, and if additional phases are found to be stable relative to the gas, they are added as
described above and the minimization algorithm is repeated. Even trace phases such as perovskite are typically present at levels >10-6 moles per mole
of elements in the complete system. If the amount of a phase drops below a minimum value, set at 10-10 moles, that phase is removed from the condensate
assemblage, and the minimization algorithm is repeated. If no phase must be added or removed after the minimization, the system is considered solved for that temperature, pressure and bulk composition, and a new
temperature step is initiated.
Convergence of each solution is assessed independently by calculation of the difference in the chemical potential of each condensed component between the gas and condensates. For temperatures >1400K, these
differences for each component are always < 10-7 of the chemical potential in the gas, and usually very much better,
e.g. ~10-12. At lower temperatures, particularly in dust-enriched systems, these differences in some cases increase for components containing the
elements Ca, Al, and Ti, and no results are reported here for any temperature step in which the difference exceeds 10-4 for any condensate component.
Even in an example where these differences are ~3x10-4, they would record uncertainties corresponding to a shift of only
~10-10 of the total Ca in the system between the gas and the condensate assemblage. These reaction imbalances occur because the algorithms call for
numerical approximation of the first and second derivatives of the Gibbs energy of the gas with respect to the concentration of each of the condensing elements in it, and this approximation becomes increasingly
sensitive to machine numerical precision at very low concentrations of elements in the gas (e.g., 10-20 moles per mole of elements in the system).
Mass balance is preserved to within <10-27 of the moles of atoms present throughout all calculations.