PHREEQC release notes

Version 3.6.0: December 2, 2019

Bug fix in pointers to the input structure that could produce an error in the .._MIX keywords.

Bug fix in transport of colloids from mobile to stagnant cells and vice versa (special basic function CHANGE_SURF). It may have produced incorrect results in versions 3.5.1 and 3.5.2.

Version 3.5.2: August 10, 2019

Implicit calculation of multicomponent diffusion was extended with 1 stagnant layer in the column, and options were added:
-- When the mixing factor (defined with keyword MIX) is omitted for the boundary solution and the adjacent stagnant cell, the implicit model calculates a mixing factor for the boundary stagnant cell, multiplying the mixing factor for the mobile cell (this is cell 1 or cell 'cells') with the ratio of the kgs of water in the stagnant cell and the mobile cell, see e.g. Fig. 2 of Appelo et al., 2010, GCA 74, p. 1205. Diffusion into the stagnant cell can be omitted, defining a mixing factor of 0. For example for a column with 'cells' = 4, stagnant cell = 6, boundary solution 0: MIX 0; 6 0.
-- Thermal diffusion with the stagnant cells will be calculated when temperatures differ by more than 0.1 oC. Multicomponent diffusion coefficients decrease with the viscosity of the solution, markedly affecting the results. File ex12b.phr in c:\phreeqc\exmpls compares traditional and multicomponent diffusive transport of heat and solutes with temperatures changing from 0 to 25 oC.

The activity-coefficient effect for multicomponent diffusion f = 1 + Δlog10(γ) / Δlog10(m) is now applied only if 0.2 < f < 2.25. The factor is not applied in implicit calculations.

Added an option to define the minimal log10(molality) of a species for including it in multicomponent diffusion. The definition is done with option -implicit false/true 1 -30 under keyword TRANSPORT, where 1 and -30 are the default max_mixf and default minimal log10(molality) (= min_dif_LM), respectively. The recommended value for min_dif_LM is -12.

Removed the limit, set in Version 3.5.1, that transported moles had to be > 1e-15 in multicomponent diffusion.

Bug fix in Version 3.5.1 that could bypass the calculation of diffusion through the Donnan layer.

Added an option to use the chemical model structure of the previous calculation for specific cells with option -same_model under keyword TRANSPORT. For example, -same_model 2-5 8-11 will use for cells 2-5 the model defined for cell 1, and for cells 8-11 the model defined for cell 7. For cells 0 and 1, and 6 and 7, the model is set-up with the default 'check_same_model' in prep.cpp.

Corrected memory leaks in implicit calculations of diffusive transport.

Version 3.5.1: May 29, 2019

Added an implicit algorithm that allows for large timesteps when diffusion is calculated. It will be used when -implicit true 1 is set under keyword TRANSPORT, where 1 = 'max_mixf' = D * Δt / Δx2. The default value of max_mixf = 1, but it can be increased often to > 4, depending on the required accuracy and the chemistry of the system, and the calculations will be about an order faster.
With implicit true, electro-migration is calculated more stable and usually without (possibly) disturbing warnings of 'negative moles', see Electro-migration and electro-remediation. The implicit calculations have been implemented for a regular column, not for stagnant cells, and not for interlayer diffusion.

The activity-coefficient effect for multicomponent diffusion was modified, previous versions could give an erroneous large correction.

The column boundary cells 0 and count_cells + 1 are now print-/punched by default when the boundary condions are not 'closed'.

Version 3.5.0: March 26, 2019

New IPhreeqc and PhreeqcRM modules may facilitate interfacing with other computer programs via libs and dlls. Download from the USGS site.

CVODE`s progress variables (overall time, time-step, calls) are displayed in the status window when kinetic rates are calculated.

Bug fix of GAS_PHASE_MIX with -fixed_pressure members.

Revised numerical method for solid solutions to respond better to KNOBS; -step and -pe_step convergence parameters.

Version 3.4.8: January 18, 2019

Bug fix of incorrect mass added in solid solutions with a component containing an element not present initially, but added with INCREMENTAL_REACTIONS.

Updated the PHREEQC help file with descriptions of the .._MIX keywords (SOLUTION_MIX, EQUILIBRIUM_PHASES_MIX, etc.).
The .._MIX keywords combine fractions of previously defined keyword blocks and facilitate the modification of keyword blocks. For an example, see keyword GAS_PHASE in the help file.

Enabled active_fraction factors for EXCHANGE_SPECIES in Pitzer activity models, and generally improved the convergence of exchange calculations with active_fraction factors.
Active_fraction factors are useful for modeling non-ideal ion exchange, see Active fraction model for varying exchange selectivity.
The aquia example now runs 3 times faster than before.

The number of mixruns in a diffusion/dispersion model is limited to the largest 32 bit integer (= 2^31) to avoid an integer fault. When exceeded, an error message is printed with the suggestion to decrease the time_step, and/or increase the cell-lengths.

Version 3.4.7: October 2, 2018

Bug fix in multicomponent diffusion with changing porosities. Without interlayer diffusion, an error resulted when the porosity in a cell was below the limiting porosity.

The activity calculation of exchange species can be omitted with option '-exchange_gammas false' in keyword EXCHANGE (default is 'true'). Previous versions ignored this option in non-Pitzer activity models.

Bug-fix for the activity calculation of Gapon-type exchange species in non-Pitzer activity models. Ca0.5X now uses (γCa2+)0.5, and log gamma (Ca0.5X) = 0.5 * log gamma (Ca2+) as should be. Previous versions calculated the gamma with the charge of the exchanger (z = -1 of X in Ca0.5X).

Version 3.4.6: June 20, 2018

The solid solution calculation in the Pitzer-model was copied from the numerical derivatives used in the ion-association model.

Bug-fix for porosities, which could be set incorrectly when porosities were not defined in TRANSPORT simulations with more cells than before.

Modified the diffusion properties for (possible) boundary cells in stagnant calculations with enhanced transport through a Donnan layer. In version 3.4.2, the harmonic mean of all the diffusional properties was introduced for all the stagnant cells.
However, many models use cell 3 and cell [2 + stagnant cells] as outer, well-mixed reservoir solutions, and diffusion is then determined by the properties of the boundary cells in the column (similar to the regular column, where cells 0 and [cells + 1] are well-mixed). If cells 3 and [2 + stagnant cells] are without a surface with a Donnan layer, they are taken as well-mixed now.
  cell 3 without donnan layer, cell 4 with a donnan layer: diffusion is determined by the properties of cell 4.
  cell 3 without donnan layer, cell 4 without a donnan layer: diffusion is determined by the harmonic mean of the properties of cells 3 & 4.
PHREEQC's choice can be manipulated, adding a surface with a very small number of sites and a tiny Donnan layer in cell 3 and/or the cells it diffuses into, and similar for cell [2 + stagnant cells]. See for example opa_col3.phr.

A Donnan layer on SURFACEs can now be calculated with Pitzer activity models.

Bug-fixes for option -numerical_derivatives in KNOBS, setting it true may improve convergence of donnan layer calculations (or not).
Overall iterations (this is, all the calls to CL1 in a reaction step) are printed in the output file if different from the iterations in the final set. It can be USER_PRINTed/PUNCHed with the variable 'iterations'.

Streamlined the output of -debug_prep in KNOBS. Added options -debug_mass_balance (total moles of elements), and -debug_mass_action (species data: log_k, dz, reaction stoichiometry).

The numbers in the matrix for CL1 can be inspected in my_array[row * (count_unknowns + 1) + col]

Version 3.4.5: May 7, 2018

Added the Soreide and Whitson interaction coefficients of H2O(g) and Mtg(g) in the code.

Bug fix of option '-new_line false/true' in SELECTED_OUTPUT, new_line was not initialized to 'true' in the phreeqc.exe.

Version 3.4.4: April 30, 2018

Removed warnings about coupling of surface and equilibrium_phase. For example, montmorillonite can be made a cation exchanger:
     Summ_Montmorillonite; Al2.33Si3.67O10(OH)2-0.33 + 12 H2O = 2.33 Al(OH)4- + 3.67 H4SiO4 + 2 H+; -log_k -44.4
     Summ Summ-
     Summ- = Summ-
  SOLUTION 1; Na 10; Cl 10; pH 7 charge; C(4) 1 CO2(g) -2
     Ca-Montmorillonite 0 1e-3
     Summ_Montmorillonite 0 0
     Summ Summ_Montmorillonite 0.33 3.11e5
     -donnan; -equil 1

Added option '-new_line false/true' in SELECTED_OUTPUT. If true, starts a punch declaration on a new line (default). If false, the new line must be written explicitly, for example: 10 punch eol$ + str_e$(tot("Na"), 18, 13)

Bug fix in the Donnan calculation of a positive surface: the surface charge could go to unsolvable numbers, but is now limited to +5e3 equivalents, similar to the limit of -5e3 equivalents for the negative surface charge.

Bug fix in the Donnan calculation of a surface coupled to a kinetic reactant: the amounts of water could be distributed incorrectly over free and DL water.

CVODE and RK now calculate the combination of an irreversible reaction and kinetics in the same way. In previous versions, CVODE added the reaction after the kinetic step was calculated, RK adds the reaction at the start.

Version 3.4.3: March 9, 2018

Bug fix in TRANSPORT with stagnant zones: solutions were not updated when mixing was not defined for the highest number of the stagnant cells.

Bug fix in multicomponent diffusion through a Donnan layer: a SURFACE without any charged-surface component could give a segmentation fault.

Version 3.4.2: February 27, 2018

Input files with experimental data obtained at high T, P and salinity are installed by the phreeqc-installer in c:\phreeqc\high_P_T.

Multicomponent diffusion in stagnant zones now uses the harmonic mean of the cell properties, where earlier versions took the average. If the properties are equal, the results will be identical; if the properties vary, the cells with slow diffusion are now given stronger weight.
In Example 21 (radial diffusion of tracers through Opalinus Clay), the clay has other properties than the filters that hold the sample. For matching the experimental results, the fit parameters were adapted (CEC, EDL enrichment of 22Na+ and Cs+, tortuosity of 36Cl-, and the geometrical factor for interlayer diffusion). For neutral species like HTO, the calculations are not affected.

The multicomponent diffusion function was rewritten to include diffusion through the EDL's of all the differently charged surfaces. (Earlier versions calculated diffusion through one EDL only, of the first surface in alphabetical order.)

Bug fix with pitzer database, when running with McInnes scaling, but without Cl in the solution.

Version 3.4.1: February 15, 2018

For compiling with Visual C++, replaced conditional statements 'if (variable == NAN)' with 'if (isnan(variable))', and similar for !=.

Bug fix for calculating density of the initial solution. For example,
    SOLUTION; -units moles/L; Na 1.5; Cl 1.5; -density 1.0 calculate
will now give:
   Na, Cl 1.549 mol/kgw, Density (g/cm³) = 1.05613

Total Carbon and CO2, when zero, are not printed in "Description of solution".

Print and punch of cells in transport calculations with stagnant zones now follows the order of the cell numbers.

Enabled multicomponent diffusion among boundary and stagnant cells, and of cells within a stagnant layer.
For example, with 2 cells in the column, 1 stagnant layer: MIX 0; 4 0.5; MIX 4; 5 0.1

In agreement with IAPWS data 0 - 250°C, changed the temperature damping coefficient of Dw for H+ to 1000 in phreeqc.dat and pitzer.dat:
    -dw 9.31e-9 1000 0.46 1e-10

Bug fix for electro-migration through a column (membrane) that contains pure water.

Updated the HELP file.

Version 3.4.0: November 9, 2017

Notepad++ installs the phreeqc-adapted files (userDefineLang.xml, shortcuts.xml, session.xml) for all users.

The effect of an electrical field on multicomponent diffusion (electro-migration) is calculated with the Nernst-Planck equation if one of the 2 column end-cell solutions has a potential different from 0 (SOLUTION 0, SOLUTION cells + 1). For examples of electro- migration, see

A number of new features have been implemented for electro-migration.

(1) The electric potential can be defined in keyword SOLUTION:

    -potential 3 # Volt

(2) The electric current in a column experiment can be defined in keyword TRANSPORT:

    -current 1e-3 # Ampere

(3) When using -multi_d and -interlayer_d in TRANSPORT, the porosity for calculating interlayer diffusion is added to the porosity defined with -porosities or -multi_d. For example:

-porosities 0.3 0.29 0.28 ...
    -multi_d true 1e-9 0.3 0.05 1.0 false
    -interlayer_d true 0.07 0.01 250

The total porosities will be 0.37, 0.36, 0.35, ... in cells 1, 2, 3, ...
In previous versions, the interlayer porosity was subtracted from the porosity defined with -multi_d.

(4) If the properties in a regular column (surface areas, lengths, double layers, tortuosities, diffusion coefficients) are different for the two cells, the harmonic mean is used for calculating the mass transfer (Appelo et al., 2010, Geochim. Cosmochim. Acta 74, 1201-1219).

(5) The following Basic functions have been added:

CURRENT_A      # the current through the column, amps.
POT_V              # potential in a cell, volts.
T_SC("Cl-")    # The transport- or transference-number of the ion, equal to the fraction of the specific conductance contributed by the species(-).
VISCOS         # At present, returns the viscosity of a pure water solution, same as VISCOS_0, milliPascal second.
VISCOS_0       # Returns the viscosity of water at current conditions, milliPascal second. Ref.: Huber et al., 2009, J. Phys. Chem. Ref. Data, Vol. 38, 101-125.

(6) A temperature damping factor and two parameters for ionic strength dependence were added for the calculation of the diffusion coefficient of a solute species:

H+ = H+
    -dw 9.31e-9 763 0.46 1e-10

where the first number is the diffusion coefficient at 25 C, and the second number is a damping factor for the temperature correction, as proposed by Smolyakov, according to Anderko and Lencka, 1997, Ind. Chem. Eng. Res. 36, 1932-1943:

Dw(TK) = 9.31e-9 * exp(763 / TK - 763 / 298.15) * TK * 0.89 / (298.15 * viscos).

The ionic strength correction is:
Dw(I) = Dw(TK) * exp(-0.46 * DH_A * |z_H+| * I^0.5 / (1 + DH_B * I^0.5 * 1e-10 / (1 + I^0.75)))

The correction is applied when the specific conductance is calculated, and can be invoked in electro-migration models, adding TRUE to -multi_d:

-multi_d true 1e-9 0.3 0.05 1.0 TRUE

ColdChem.dat database has been expanded to include sulfate and sulfate salts. The database is described in a new journal article by Toner and Catling: Journal of Chemical & Engineering Data, 2017,

Error in dist_x function and in counting of cells for some transport problems.

Added new database, core10.dat, contributed by Marc Neveu. The database is derived from phreeqc.dat and llnl.dat, with careful data checking, use of SUPCRT for temperature dependence, and addition of molar volumes for pressure dependence. Downloaded July 18, 2017 from the following URL

Please see the following report for details: Neveu, Marc, Desch, S.J., and Castillo-Rogez, J.C., 2017, Aqueous geochemistry in icy world interiors: Equilibrium fluid, rock, and gas compositions, and fate of antifreezes and radionuclides, Geochimica et Cosmochimica Acta, v. 212, p. 324-371,

Changed the first parameter in the analytical expression for ethane C2H6. The negative sign was apparently omitted.

GAMMA and LG functions were not correct for exchange species.

Bug fix for when -no_check followed mineral name in PHASES. Now catches error if equation does not the next line following the mineral name.

Bug fix for -donnan, -debye, and surface related to mineral that has zero moles.

Added a new option for a density calculation to SOLUTION and SOLUTION_SPREAD. The speciation calculation for a solution is iterated until the density used in the conversion of units is the same as the density calculated for the solution. It affects solutions with concentrations defined as per liter (mg/L, for example).
Here are three examples for specifying the calculation, where \t indicates a tab.
-density 1.1 calculate

-density 1.2 calculate

pH\t density\t
charge\t calculate\t
7.0\t 1.03\t
If the calculation is performed, the density in the description of solution part of the output will be identified with "(Iterated)".

Added new Basic function TITLE, which returns the last definition by a TITLE keyword data block.

Version 3.3.11: March 10, 2017

Added the ColdChem.dat database developed by Toner and Catling, Journal of Chemical & Engineering Data, 2017, DOI: 10.1021/ acs.jced.6b00812. It is a comprehensive Pitzer model in the Na-K-Ca-Mg- Cl system valid for low temperature, from 298.15 to < 200 K. PHREEQC was modified to allow definition of the temperature dependence of A(phi) in the Pitzer formulation through a new option, -APHI, in the PITZER keyword. The coefficients of the APHI expression are the same as the temperature dependence of other Pitzer parameters. A new Basic function, APHI, has been added to give the value of
the parameter in the current calculation.

Added logic to limit the range of the fugacity coefficient.

Change to avoid potential NULL pointer in -multi_D calculations.

Error check that the same element name was not used for both a SURFACE and EXCHANGE, which caused a NULL pointer. Some variables in the unknown structure are now initialized.

PhreeqcRM: The Fortran version of RM_SpeciesConcentrations2Module did not provide the correct return value. Generally, -3 was returned (IRM_INVALIDARG), even if the method worked properly and should have returned 0 (IRM_OK).

Version 3.3.10: January 12, 2017

Updated Basic function SYS("equi"...). The function now gives the correct amounts of equilibrium phases at the end of a reaction. Previously, the function erroneously gave the amounts of equilibrium phases before reactions.

Version 3.3.9: November 15, 2016

Updated sit.dat to version 9b0 from Changed documentation of some RATES definitions in Amm.dat to have the correct diameter assumed for mineral grains.

Bug when surface was related to kinetics and the proportion of sites to kinetic reactant was zero (not common). The bug caused a segmentation-fault crash.

PhreeqcRM: Previously, the guess for the unknown for mass of water was the value from the last time step. For widely varying saturations, the old value caused non-convergence for the cell. Now the mass of water is estimated from the updated value of the moles of oxygen in the cell.

Version 3.3.8: September 13, 2016

PhreeqcRM: Added comment to output that user-grid cell numbers are 0 based (C-like).

PhreeqcRM: Removed tabs from RM_interface.F90.

Fixed a potential bug where memcpy target and source locations overlapped.

Implemented some changes for redox calculations with the Pitzer model. Revising master variable values when large concentrations occur, and basis switching between redox environments. Algorithm for initial estimates of variables was also changed, as well as removing activity coefficient equations when not needed.

The file dw.cpp was removed from the source code. It contained obsolete code inherited from PHRQPITZ to calculate the density of water.

PhreeqcRM: Order of header files was changed to put mpi.h first to avoid potential header conflicts with some MPI implementations.

The Peng-Robinson gas calculation failed with very small gas partial pressures. Code was modified to set a minimum of 1e-10 atm for the calculation.

PhreeqcI: File/registry virtualization was disabled. Windows allowed users to write into system directories (Program Files for instance) even though they did not have sufficient permissions. The files were actually stored in another directory in the user's profile. The files are difficult to find, and, under some conditions, PhreeqcI would fail. Now, it is not allowed to write into directories for which the user does not have write permission.

Added new Basic function PHASE_VM("Calcite"). The function returns the molar volume for a mineral. The molar volume is defined for the mineral in PHASES with the -vm option. Use the Basic function GAS_VM for gas components.

Added two new Basic functions related to KINETICS.

n$ = KINETICS_FORMULA$("Albite", count, elt$, coef)

This function searches for a kinetic reaction definition (Albite, in the example). If found, n$ is set equal to the first argument (Albite), otherwise an empty string is returned. The function returns these values: count is the number of items in the arrays elt$ and coef; elt$ is a list of element names in the formula for the kinetic reaction; and coef is a numeric array containing the number of atoms of each element in the kinetic-reaction formula, in the order defined by elt$, which is alphabetical by element.

m = SYS("kin", count , name$ , type$ , moles)

This function identifies all of the kinetic reactants in the current KINETICS definition and returns the sum of moles of all kinetic reactants. Count is number of kinetic reactants. Name$ contains the kinetic reactant names. Type$ is "kin". Moles contains the moles of each kinetic reactant. The chemical formula used in the kinetic reaction can be determined by using a reaction name from Name$ as the first argument of the KINETICS_FORMULA$ Basic function.

Updated the online manual ( to have revisions noted in italic text and better graphics for figures.

Version 3.3.7 April 21, 2016

Cerussite misspelled in 5 databases (Amm, minteq, minteq.v4, phreeqc, and wateq4f. Modified Calcite RATES to handle M0=0. Added rate for quartz to wateq4f.dat.

PhreeqcI: The selected-output file was missing headers under some conditions, for example : SOLUTION;SELECTED_OUTPUT;END.

PhreeqcI: Fixed bug that overwrote selected output files when (re-)opening a .pqi file

Version 3.3.6 April 18, 2016

PhreeqcRM: InitialPhreeqcCell2Module method had incorrect indexing when using a mapping of model cells to chemistry cells that was not 1-to-1. The result was incorrect values for porosity, rv, and saturation.

Additional error checking. In SOLUTION_MASTER_SPECIES the master species must contain the element name.

If an aqueous species lacked a molar volume parameter (Vm), then the change in molar volume for reactions involving that species were not calculated correctly. Now summations of molar volumes, including zero for the species, are calculated.

Added Vm values for Gibbsite, Kaolinite, Albite, Anorthite, K-feldspar, Ca-Montmorillonite, Talc, Illite, Chrysotile, Sepiolite, Hematite, Goethite, Pyrite, and Mackinawite for the phreeqc.dat and Amm.dat databases.

PHREEQC and PhreeqcI: PhreeqcI has updated documentation through the Help menu. All entries in green font are modifications, corrections, or new features.

Corrected description of "EQUI" option of the SYS function.

-porosities is now an option in TRANSPORT for defining a heterogeneous distribution of porosity for -multi_d calculations.

RHO_0 is a new Basic function that gives the density of pure water at the current temperature.

EDL_SPECIES(surf$, count, name$, moles, area, thickness)

Returns the total number of moles of species in the diffuse layer. The arguments to the function are as follows: surf$ is the name of a surface, such as "Hfo", excluding the site type (such as "_s"); count is the number of species in the diffuse layer; name$ is an array of size count that contains the names of aqueous species in the diffuse layer of surface surf$; moles is an array of size count that contains the number of moles of each aqueous species in the diffuse layer of surface surf$; area is the area of the surface in m^2; thickness is the thickness of the diffuse layer in m. The function applies when -donnan or -diffuse_layer is defined in SURFACE calculations.

Corrected description of "EQUI" option of the SYS function.

Corrected description of "EQUI" option of the SYS function. PhreeqcI has updated documentation through the help. All entries in green font are modifications, corrections, or new features.

PhreeqcI: Changed the way the current working directory is set. Now, the current working directory is not changed when browsing for and changing the database.

PhreeqcI: SOLUTION_SPREAD did not import correctly when using the P4W tab of PhreeqcI to generate the an example problem.

In converting units of mass to moles per kilogram of water, it is possible to calculate a negative mass of water. This problem is now trapped, and an error message is issued.

Version 3.3.5 February 3, 2016

Created an option to enter cell porosities in keyword TRANSPORT, which was motivated by the need to define porosities for a variable-porosity multicomponent diffusion problem. Previously the Basic function CHANGE_POR had to be used to allow variable porosity.

-porosities 0.3 0.29 0.28 # porosities of cells 1, 2, 3-n, used for calculating tortuosity in multicomponent diffusion calculations. The numbers entered here take precedence over the value given with -multi_D. If one stagnant layer is defined together with an exchange-factor > 0 ('-stagnant 1 4e-6 0.3 0.1'), the mobile (here = 0.3) and immobile (here = 0.1) porosities defined with -stagnant are used.

IPhreeqc: Added optional argument sl to GetSelectedOutputValue in Fortran. The revised method is as follows, where sl is a new optional argument.

status = GetSelectedOutputValue(id, i, j, vt, dv, sv, sl)

In the example, if the character string generated in IPhreeqc and stored in sv is 40 characters or less, the value of sl will be 0. If the character string generated in IPhreeqc was greater than 40 characters, the truncated string is stored in sv, but the value of sl will contain the length of the string before truncation.

All programs: Now using Jenkins for continuous integration.

Erroneous input could cause a segmentation fault for surfaces related to minerals or kinetics.

Minor format change for specific conductance.

PhreeqcRM: Advection example in C and C++ had wrong MPI type for integer transfers. Example code in documentation also had wrong MPI type for some integer transfers. Integer transfers in Fortran should use MPI_INTEGER, while integer transfers in C and C++ should use MPI_INT.

PhreeqcRM: Revised transfer of data for MPI when re-balancing. Serializes lists of characters, integers, and doubles, instead of writing full "dump" version of each reaction. Is much faster for cell transfers when re-balancing among processes.

PhreeqcRM: Fixed bug when transferring data for MPI. Surface related to kinetics lost concentrations of sorbed elements.

New set of convergence parameters that delay the removal of an unstable phase by 1 iteration. Also new KNOBS identifier, -equi_delay n, where n is the number of iterations to retain an unstable phase.

Comments for K-spar and Albite RATES were incorrect for phreeqc.dat and wateq4f.dat

PhreeqcRM: Provided better error message for CreateMapping errors.

Version 3.3.3 October 23, 2015

PhreeqcRM: C function clock() does not provide sufficient resolution. Efficiency calculation generated NaN, and re-balancing was poor. Revised to use MPI and OpenMP timers when available.

PhreeqcRM: Efficiency calculation generated NaN for OpenMP calculations.

Fixed formula for pressure dependence of B1, B2, F1, F2 in the Pitzer formulation. Previously had a limit of -10C; temperatures less than -10 produced a floating point exception and convergence failure.

IPhreeqc and PhreeqcRM: Debug version asserted when trying to write an error message to std::cerr.

PhreeqcRM: Added methods and documentation for RM_GetEndCell and RM_GetStartCell for C and Fortran. C++ methods already existed. The methods give the range of chemistry cell numbers that are assigned by each worker (OpenMP) or process (MPI). If re-balancing is enabled, these numbers may change as cells are shifted among workers or processes.

Version 3.3.2 October 2, 2015

PhreeqcRM: Bug in configure scripts required developer zlib to be installed. Now zlib is only needed if the --with-zlib option is given in order to write compressed dump files (USE_GZ).

Fixed bug with SOLUTION_SPREAD. Element names starting with "[" were not handled correctly.

Handled case where GAS_PHASE components were not defined in PHASES.

New R version submitted.

PhreeqcRM: Added new method SetScreenOn (RM_SetScreenOn) to control messages about re-balancing and any messages written with ScreenMessage (RM_ScreenMessage).