Hydrothermal precipitation of arsenic compounds in the ferric–arsenic (III)–sulfate system: thermodynamic modeling

Minerals Engineering 16 (2003) 429–440 This article is also available online at: www.elsevier.com/locate/mineng

Hydrothermal precipitation of arsenic compounds in the ferric–arsenic (III)–sulfate system: thermodynamic modeling R.C.M. Mambote a, P. Krijgsman b, M.A. Reuter a

a,*

Department of Applied Earth Sciences, Section of Raw Materials Technology, Delft University of Technology, PO Box 5028, 2600 GA Delft, Netherlands b Ceramic Oxides International B.V., PO Box 8, AA Wapenveld, Netherlands Received 19 July 2002; accepted 7 October 2002

Abstract It is known that Fe (III) and As (III), through in situ oxidation, can precipitate to form a crystalline ferric arsenate under specific conditions of temperature and pH. A thermodynamic model is built to investigate this process. The development and the use of chemical speciation model have lead to simulate the distribution and the concentration of chemical species of the ferric–arsenate (III)–sulfate system at 293.15–553.15 K temperatures range. The modeling predictions are compared to experimental data measured by the authors. It was shown that a reasonable prediction could be made if different modeling approaches are used. The paper compares different approaches in thermodynamic modeling and discusses also the validity of different extrapolation models used at high temperatures.
2003 Published by Elsevier Science Ltd. Keywords: Modeling; Process control; Waste processing; Process optimization

1. Introduction Arsenic is a by-product of the mining of non-ferrous metals in particular, copper, zinc, lead, gold and silver. The major source of arsenic production and emissions is the copper processing industry (Loebenstein, 1994). Arsenic is viewed by the metal industry as a problem, not as a benefit due to its low market price, environmental issues associated with its production and the external costs which is a burden to the non-ferrous metal processing industry. During the past decade a combination of both lower grade ores and environmental concerns has resulted in accelerated technological development to meet this arsenic challenge. In addition to this existing problem, in the near future, there will be an increase of arsenic wastes due to the processing of refractory gold-ores with high arsenic contents.

* Corresponding author. Address: Department of Applied Earth Sciences, Delft University of Technology, Mijnbouwstrat 120, Delft 2628 RX, Netherlands. Tel.: +31-15-278-1636; fax: +31-15-278-2836. E-mail addresses: [email protected] (R.C.M. Mambote), [email protected] (M.A. Reuter).

These reasons can justify why literature about arsenic abounds (Riveros et al., 2001); however, its chemistry is still not well understood. The major advances in this field over the last years have been the hydrothermal precipitation of arsenic with iron, to synthesize a ferric arsenate crystalline phase. From Fe(III)–As(III)–SO4 – H2 O system, it was noticed that the immobilization of arsenic in crystalline form, by in situ oxidation of As3þ to As5þ under hydrothermal conditions, is technically feasible. Although FeAsO4 can be precipitated at a temperature less then 473.15 K with initial solution of As(V) and Fe(III) (Monhemius and Swash, 1999), it was found that all precipitates obtained at temperatures less then 483.15 K were amorphous. The results revealed that crystalline ferric arsenate can be obtained by hydrothermal precipitation using the Fe–AsO2 –SO4 system at temperatures above 483.15 K, a pH less then 2.5 and hydrogen peroxide as an oxidant (Mambote et al., 2001). Process design requires a solid theoretical framework, which will allow prediction of the chemical reactions. As suggested by Sandler (1994), thermodynamic modeling and phase equilibria are the heart of chemical process design. Existing theoretical approaches to modeling the

0892-6875/03/$ - see front matter
2003 Published by Elsevier Science Ltd. doi:10.1016/S0892-6875(02)00253-4

430

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

Nomenclature a Ac A/ b b Bc Bc B_ bM=X c C c CM=X f ðIÞ

fi fc Fi DHT00 I K Kc

constant equal to 0.001875 Debye–H€ uckel parameter Debye–H€ uckel parameter is related to Ac (Eqs. (41) and (42)). constant equal to )12.741 parameter assigned to a constant value of 1.2 in the Eq. (41). Bromley B ion–ion interaction term Debye–H€ uckel parameter B-dot parameter ionic binary interactions term constant equal to 7:84  104 Pitzer ion–ion interaction term is the ternary interaction term Debye–H€ uckel function in the Eq. (40) (function of temperature, solvent properties and ionic strength, I) fugacity coefficients fugacity coefficient amount of inflowing material i (mol) standard reaction enthalpy ionic strength of the solution equilibrium constant, generally a function of temperature and pressure equilibrium constant at concentrated solution, designated by Fillippou et al. (1995) as a mass stability constant (bT ).

high-temperature ionic solution chemistry may be divided into two major categories: • Ion association known as chemical speciation model, it assumes that ions form stables complexes. The free energy of the complex formation reaction may then be evaluated at elevated temperatures by approximating the dependence of ionic heat capacities on temperature from regression or by considering forces involved in the binding of the ions in complexes. • Ion interaction, known as the Pitzer approach, assumes no association, but interactions between simple ions. All physico-chemical processes, including solubility or deposition of different salts, are known to be determined by the ionic composition. The theoretical prediction of thermodynamic properties of multiple electrolyte solutions remains one of the outstanding problems of physical solution chemistry. One of the most important tasks in aqueous solution chemistry is speciation, i.e. the exact description of a system with different metal ions and ligands in terms of the concentrations of species present at given conditions (tem-

KT mi mi Nc NR Ns pi R DST00 r T V X

c h m þ , m x yi zi

represents the equilibrium constant at infinite dilution ðI ¼ 0Þ. concentration of the species molality of the specie i number of components represents the number of independent chemical reactions of the solution number of species stoichiometric coefficient of product i constant equal to 8.314 J K1 mol1 standard reaction entropy change stoichiometric coefficient of the reactant i reaction temperature (K) total vapor (mol) concentration of other species present in the 3 system (SO2 4 , H2 O2 and AsO4 ) in the Eq. (44) activity coefficient constant equal to 219 value of the stoichiometric coefficient of the cation and anion (Eq. (39)) constant equal to 1.00322 vapor-phase mole fraction for i charge of the ith ion in the solution

perature, ionic strength and inert salt) and known total concentrations of the components. The composition of coexisting phases may in suitable cases be used to obtain important information on the P , T and chemical conditions of a chemical reaction. Two other approaches to speciation computation exist as reported by Stumm and Morgan (1996). An alternative approach based on mass action laws that consider the Gibbs free energy of the chemical system. The speciation based on the equilibrium constant approach is formally equivalent to minimizing the free energy of the same chemical system. The constrained optimization problem requires a more complicated algorithm for its solution (Holstad, 1999). Another approach is based on the algorithmic combination of the overall reaction kinetics with a procedure that minimizes the Gibbs energy of the system. The time course of the reactant consumption is assumed to follow the experimental Arrhenius rate law. The reaction path computation may give valuable insight into the behavior of the system, but the strategy requires more computation and its accuracy is more difficult to control (Koukkari et al., 1997; Holstad, 1999). Another thermodynamics approach to model a system concern the atomistic analysis (statistical mechan-

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

ics) by studying the canonical partition function that carries all the thermodynamic information of system. Once the partition function has been determined, the thermodynamic properties of the system may be readily obtained. The use of this approach has been the subject of criticism (Lupis, 1983). The limitations involved in applying this approach are: • The models are built with many assumptions about the energy contribution of a system of atoms. • The models do not distinguish between liquids and disordered solids. The most widely known modeling method for this approach is the ‘‘quasi-chemical approximation’’. This approach can help to predict the excess functions of the system. Unfortunately, it fails to give more insight through the system under investigation. The present work aims to understand and to optimize the ferric arsenate reaction at elevated temperatures. Several works have studied the thermodynamics of the ferric arsenate reaction (Robins, 1981; Welham et al., 2000a,b), most of these works were more focused on the stability aspect of ferric arsenate compound. This work presents some novel features that were not discussed before, these are: • The modeling of arsenic and iron (III) system in temperatures range above 373.15 K data measured; • The study of arsenic (III) oxidation at high temperatures; • The comparison of the semi-empirical method versus the quasi-Newton method for numerical solution.

specify the various species, components present and reactions involved in solution. Casas et al. (2000) defined the term species as any chemical entity present in solution such, ion, complex or molecule. The term components as a minimum number of species that allows complete description of the system. The number of components for a given system is constant and expressed by the following equation: Nc ¼ Ns  NR

ð1Þ

The selected major reactions involved in this system are summarized below. Reaction equations in the Fe(III)–As(III)–SO4 system Fe2 ðSO4 Þ3 þ As2 O3 þ 2H2 O2 þ H2 O () 2FeAsO4 þ 6Hþ þ 3SO2 4 þ

ð2Þ

H2 O () H þ OH



H2 OðgÞ () H2 OðlÞ

þ

K1 ¼ ½H ½OH=½H2 O

K2 ¼ ½H2 O=ðpH2 OÞ

H2 O2 () 2H2 O þ O2

K3 ¼ ½H2 OðpO2 Þ=½H2 O2 

3 K4 ¼ ½Fe3þ 2 ½SO2 4  =½Fe2 ðSO4 Þ3 

FeOH

2þ

þ

þ H () Fe

3þ

2.1. Basic mathematical model for the three phase Fe(III)–As(III)–SO4 The conceptualization of this basic model is based on the methodology proposed by Rafal et al. (1994). In order to describe an ionic system, the first step is to

ð4Þ ð5Þ ð6Þ

þ H2 O

K5 ¼ ½Fe ½H2 O=½FeOH2þ ½Hþ  3þ

þ FeðOHÞ2

þ

þ H () FeOH

2þ

ð7Þ

þ H2 O þ

K6 ¼ ½FeOH ½H2 O=½FeðOHÞ2 ½Hþ  2þ

4þ Fe2 ðOHÞ2

ð8Þ

þ 2Hþ () 2Fe3þ þ 2H2 O 2

2

4þ

2

K7 ¼ ½Fe3þ  ½H2 O =½Fe2 ðOHÞ2 ½Hþ  þ

3þ

K8 ¼ ½Fe3þ ½H2 O =½FeðOHÞ3 ½Hþ 3 þ

FeOOH þ 3H () Fe

3þ

ð10Þ

þ 2H2 O

K9 ¼ ½Fe ½H2 O =½FeOOH½Hþ 3 3þ

ð9Þ

þ 3H2 O

3

2

ð11Þ

Fe2 O3 þ H2 O () 2FeOðOHÞ 2

K10 ¼ ½FeOðOHÞ =½Fe2 O3 ½H2 O

ð12Þ

HSO 4

þ () SO2 4 þH  þ K11 ¼ ½SO2 4 ½H =½HSO4 

2AsO 2

ð13Þ

þ

þ 2H () As2 O3 þ H2 O 2

þ K12 ¼ ½As2 O3 ½H2 O=½AsO 2  ½H 

2. Conceptualization of the thermodynamic model

ð3Þ

Fe2 ðSO4 Þ3 () 2Fe3þ þ 3SO2 4

FeðOHÞ3 þ 3H () Fe Nevertheless, applying speciation model in aqueous system at elevated temperature can provide a wealth of information that is difficult to obtain in experimental way. Indeed, the determination of the speciation in aqueous system presents considerable difficulties, due to lack of analytical techniques for in situ determination of some ions and metallic complexes. At the risk of oversimplifying the system, a combination of experimental measurements and non-ideal thermodynamic model must be used in order to establish the solution speciation in the real system. In this study, emphasis is carefully given to this aspect.

431

2

ð14Þ

  AsO3 4 () AsO2 þ O2 þ 2e 3 K13 ¼ ½AsO 2 ðpO2 Þ=½AsO4 

AsO3 4

þ

þ 2H ()

K14 ¼

AsO 2

ð15Þ

þ H2 O2

3 þ 2 ½AsO 2 ½H2 O2 =½AsO4 ½H 

ð16Þ

3 þ HAsO2 4 () AsO4 þ H 2 þ K15 ¼ ½AsO3 4 ½H =½HAsO4 

ð17Þ

432

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

2 þ H2 AsO 4 () HAsO4 þ H 2 þ K16 ¼ ½HAsO2 4 ½H =½H2 AsO4 

ð18Þ

þ H3 AsO4 () H2 AsO 4 þH

K17 ¼ Fe

3þ

þ

yH2 O þ yO2 ¼ 1

þ ½H2 AsO 4 ½H =½H3 AsO4 

AsO3 4

ð19Þ

() FeAsO4

K18 ¼ ½FeAsO4 =½Fe3þ ½AsO3 4 

ð20Þ

3þ FeHAsOþ þ HAsO2 4 () Fe 4

K19 ¼ ½Fe

3þ

FeH2 AsO2þ 4

() Fe

ð21Þ

2þ K20 ¼ ½Fe3þ ½H2 AsO3 4 =½FeH2 AsO4 

ð22Þ

The overall Eq. (2) will not be considered for the model, the independent aqueous intraphase reactions should be taken in account. By checking the formation of DG for species and considering that most As(III) is oxidized to As(V), some equilibrium reactions were eliminated. Subsequently, the equilibrium constant equations for the remaining 20 reaction equations can be written by utilizing the Eq. (23). K¼

Now, the material balance for the total amounts of arsenic, iron, sulfate and hydrogen can be written by the following balance equations: Iron:

þ mFeH2 AsO2þ þ mFeOOHðaqÞ þ 2mFe2 ðOHÞ4þ ÞWk 4

þ H2 AsO4

p1

p2

pP

r1

r2

rR

ðcP 1 mP 1 Þ ðcP 2 mP 2 Þ ðcpP mpP Þ

ð25Þ

ðmFe3þ þ mFeOH2þ þ mFeðOHÞþ2 þ mFeHAsOþ4

þ ½HAsO2 4 =½FeHAsO4  3þ

It known that the mole fractions of all the vapor species must add up to unity, the vapor balance is depicted by the Eq. (25).

ðcR1 mR1 Þ ðcR2 mR2 Þ ðcRR mRR Þ

ð23Þ

Based on these equations computed with the Eq. (23), 20 non-linear equations with 26 unknowns have to be solved. The unknowns are given in the Table 1. 2.2. Charge and mass balance equations

2

¼ FFeAsO4 þ 2FFe2 O3 þ FFeðOHÞ3

In principle the mass balance of arsenic (III) is given by Eq. (27), according to the speciation of arsenic (III) in aqueous solution as reported by Pettine and Millero (2000). Arsenic (III): ðmAsO2 þ mAsO3 þ mH2 AsO3 þ mHAsO2 3

3

þ mH3 AsO3 ÞWk ¼ FAs2 O3

ð27Þ

Arsenic (V): ðmAsO3 þ mH2 AsO4 þ mHAsO2 þ mFeH2 AsO2þ 4

4

AsðT Þ ¼ x½mAsðIIIÞ ðT Þ þ ð1  xÞ½mAsðVÞ ðT Þ

½ðmHþ Þ þ 3ðmFe3þ Þ þ 2ðmFeOH2þ Þ þ ðmFeðOHÞþ2 Þ

ðmSO2 þ mHSO4 ÞWk ¼ 3FFe2 ðSO4 Þ3

4

þ 3ðmAsO3 Þ þ 2ðmHAsO2 Þ þ ðmH2 AsO4 Þ 4

ð24Þ

ðmHþ þ mOH þ 2mFeðOHÞþ2 þ mFeOH2þ þ mFeOOHðaqÞ þ mHSO4 þ mHAsO2 þ 2mH2 AsO4 þ 2mFeH2 ASO2þ 4

2

Unknowns

Neglected.

4

þ mFeHAsOþ4 þ mFe2 ðOHÞ4þ ÞWk þ V ð2yH2 O Þ

Table 1 List of variables used in the model

a

ð30Þ

Hydrogen:

2

¼ ½2ðmSO2 Þ þ ðmOH Þ þ ðmHSO4 Þ þ ðmAsO2 Þ

nFe2 ðSO4 Þ3 nAs2 O3 nH2 O2 nFeAsO4 nH2 O yH2 O mHþ mSO2 4 yO2 mOH mFe3þ mFeOH2þ mFeðOHÞþ2

ð29Þ

4

þ ðmFeHAsOþ4 Þ þ 2ðmFeH2 AsO2þ Þ þ 4ðmFe2 ðOHÞ4þ Þ

1 2 3 4 5 6 7 8 9 10 11 12 13

ð28Þ

The system involves a redox equilibrium, therefore, Eqs. (27) and (28) will form one equation of the total arsenic mass balance that will be considered in the system equations. Sulfate:

4

4

þ mFeHAsOþ4 ÞWk ¼ FH3 AsO4 þ FFeAsO4

The charge and the mass balance equations are given by Eqs. (24)–(31).

4

ð26Þ

Unknowns 14 15 16 17 18 19 20 21a 22a 23 24 25 26

nFe2 O3 mAsO2 mHSO4 mAsO3 4 mHAsO2 4 mH2 AsO4 mH3 AsO4 mFeHAsOþ4 mFeH2 AsO2þ 4 V mFeOOHðaqÞ mFe2 ðOHÞ4þ 2 mFeðOHÞ3 ðaqÞ

¼ 2FH2 O þ 2FH2 O2 þ 3FH3 AsO4

ð31Þ

where Wk is the constant representing the ratio of the number of moles of water in this system to the number of moles of water in 1 kg. Thus Wk is approximately equal to H2 O/55.508, which converts all molalities to moles. In this material balance, the oxygen balance is eliminated because the writing of the electroneutrality balance makes the writing of both hydrogen balance and oxygen balance redundant. Therefore, a system of 20 equilibrium equations, one charge balance, four material balances, one vapor balance and 26 unknowns, describe the system to model the chemistry mathematically. Thus,

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

the model developed involve a square (m m) solution matrix, i.e. the number of unknowns equals the number of independent equations.

3. Solution model The speciation computation involves the calculation of the concentration of each individual chemical species in a multicomponent-multiphase chemical system. The numerical problem is to solve the system as described above of coupled linear and non-linear equations subject to the constraint that all unknowns are positive. Literature in chemical speciation model abound, surprisingly most workers use the NewtonÕs method when dealing with this computation. However, there are two major disadvantages with this method (Penny and Lindfield, 2000); • The method may not converge unless the initial approximations are good. • The method requires the user to provide the derivatives of each function with respect to each variable. The user must therefore provide n2 derivatives and any computer implementation must evaluate the n functions and the n2 derivatives at each iteration. To deal with this problem a number of techniques have been proposed but the group of methods which appears most successful is the class of methods known as the quasi-Newton methods (Holstad, 1999; Penny and Lindfield, 2000; Brassard and Bodurtha, 2000). The quasi-Newton methods avoid the calculation of the partial derivatives by obtaining approximations to them involving only the function values. The set derivatives of the functions evaluated at any point X r may be written in the form of the Jacobian matrix; Jr ¼ ½ofi ðX r Þ=oxj  for i ¼ 1; 2; . . . n and j ¼ 1; 2; . . . n ð32Þ The quasi-Newton method provides an updating formula which gives successive approximations to the Jacobian for each iteration. In this paper, the BroydenÕs method was chosen to solve this system of non-linear equations. Details about this numerical technique can be found in the numerical method textbooks (Dennis and Schnabel, 1996). 3.1. Thermodynamic data In building a deterministic mathematical model of aqueous systems K plays an important role. Mathematically, K are coefficients of the polynomial equations, the solution will depends how accurate the KÕs are. The equilibrium constant K is a function of both temperature and pressure but not of composition (concen-

433

tration). The thermodynamic equilibrium at different temperatures can be extrapolated by several methods reported in literature (Helgeson, 1967; Zemaitis et al., 1986; Martell and Smith, 1982; Rafal et al., 1994). The extrapolation of equilibrium constants beyond the rage of experimental data has an uncertainty associated with the failure of the chosen model to represent accurately physical reality, but also which is from random errors originating from the parameters in the model. The thermodynamic data used in this work are derived from different sources (Martell and Smith, 1982; Zemaitis et al., 1986; Rafal et al., 1994). The computer database of the Outokumpu HSC chemistry was also used, this software uses the Criss–Cobble method for the Cp extrapolation. Several methods exist to extrapolate K at elevated temperature: • The HelgesonÕs method  
DST00 h log KT ¼ 1  exp expðb þ aT Þ T0  x 2:303RT  DHT00 T  T0 cþ ð33Þ  h 2:303RT • If the Cp is known K can be computed by the method reported by Zemaitis et al. (1986);     Z T DG0R DHR0 1 1 ln KT ¼ DCpR0 dT    T T0 RT0 R T0 ð34Þ • The vanÕt Hoff equation   DHR0 1 1 ln KT ¼ ln K0   T T0 R

ð35Þ

The reliability of the solution model depends of the extrapolation method used and the thermodynamic sources. There are some simple calculations that can be made to check the continuity of formation once the species of have been established. It follows that the stepwise application of equilibrium constants are related as: Kn > Kn1 > Kn2 There are examples in the databases where this requirement is not met. For example the Fe(III) hydrolysis constants given by Stumm and Morgan (1996) or other sources. Another problem is related to the uncertainty due to the discrepancy between sources as depicted in Table 2. To adopt an extrapolation method to calculate the K constant for aqueous complexes at temperatures and pressures greater than 298.15 K and 0.1 MPa, several extrapolation methods were reviewed, the discussion of these methods have been reported in literature (Helgeson, 1967; Rafal et al., 1994 and Sverjensky et al., 1997). Predictions of KðT ; P Þ rely on the estimated or measured

434

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

Table 2 Discrepancy between thermodynamic sources Sources

K values for the Eq. (7) at 298.15 K Fe3þ þ H2 O () FeOH2þ þ Hþ

HSC database Butler, 1964; Stumm and Morgan, 1996 Martell and Smith, 1982

6:479  103 8:9  104 6:456  103

standard partial molal enthalpies or entropies of dissociation, together with assumptions about the temperature dependence of the standard partial molal heat capacities of the reaction which are generally unreliable at temperatures above 423.15–523.15 K, depending the assumption. After evaluation of these extrapolation methods, the Helgeson equation was adopted to compute K in this work. The activity coefficient of inorganic solutes can be calculated by using different models reported by VasilÕev (1962), Zemaitis et al. (1986) and Pitzer (1991). According to the Debye–H€ uckel expression, the activity coefficient for an individual ion is given by pffiffi Ac z2i I pffiffi  log ci ¼ ð36Þ 1 þ a0i Bc I uckel parameters, Here Ac and Bc are the Debye–H€ which are characteristic of the solvent at a specified temperature and pressure. The quantity a0i has a value dependent on the effective diameter of ion. The ionic strength is defined by 1X mi z2i I¼ ð37Þ 2 Reasonable estimates of activity coefficients can often be made with the Eq. (36) at ionic strengths extended up to approximately 0.1. Unfortunately, in many industrial processes higher ionic strengths are encountered. A method for estimating ionic coefficients in concentrated aqueous solutions has been developed by Helgeson. This procedure is based on the modification of Debye– H€ uckel model, the B-dot equation, expressed as follows: pffiffi Ac z2i I pffiffi þ B_ I log ci ¼  ð38Þ 1 þ a0i Bc I Presently, there is no single model that can describe the solution speciation over a wide range of concentrations and temperatures, however, PitzerÕs activity coefficient model is widely accepted (Rafal et al., 1994). This equation is: " #   1:5 2mþ m c c 2 2ðmþ m Þ ln c ¼ jZþ Z jf þ m B þm C c m m

ln c ¼ ZM ZX f ðIÞ þ þ

X

h

X  i X h mX bM=X þ mZ CM=X X

mM bM=X þ

X

 i mZ CM=X

ð40Þ

M

Moreover, the PitzerÕs equations are based on virial expansion theory and require a fairly large set of interaction coefficients to model systems of great complexity. The Debye–H€ uckel model used in PitzerÕs equations is not the usual Debye–H€ uckel charging formulation exemplified in Eqs. (36)–(38), but a different one derived by Pitzer and called the Debye–H€ uckel-osmotic model. The function f ðIÞ is given by   pffiffi 4A/ I f ðIÞ ¼  lnð1 þ b I Þ ð41Þ b 2:303Ac ð42Þ 3 The advantage of PitzerÕs model is that it can predict ionic speciation for concentrated systems, i.e. ionic strength up to 10 M. Its disadvantage however consists of the large number of parameters required to quantify the physico-chemical interactions among all the components present in solution. Additionally, several parameters are temperature dependent. Due to the reasons as explained above, this model will not be used in this study. According to Fillippou et al. (1995), aqueous equilibria in relatively concentrated solutions that cannot be described adequately with the use of thermodynamic equilibrium constants. From both the experimental and numerical points of view, it is convenient to include the non-ideality of the solution in the equilibrium constants of the reactions. A semi-empirical relationship proposed by Davies in 1962 reported by Casas et al. (2000), derived from Debye–H€ uckel theory, can be used to correct the equilibrium constant:  pffiffi  I 2 pffiffi  0:2I logðKc Þ ¼ logðKT Þ þ AðT ÞDZ ð43Þ 1þ I

A/ ¼

3.2. Solving the mathematical model

ð39Þ

The following points describe the solution method used to solve the chemical speciation model.

The Eq. (39) can also be written as follows for the computation of an electrolyte MX :

• Define the system in term of key components; • assuming that the concentration is equal to activity;

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

• estimate both the concentration and I; • correction of the K values by using Davies equation; • following the recommendation of Garrels and Christ (1965), the activity coefficient for neutral aqueous species of a polar nature are set to unity, i.e. the equation is log c ¼ 0; • activity coefficient calculation by using the extended Debye–H€ uckel equation (computation with Mathcad); • the non-linearity of the equations must be solved by using the logarithmic form (quasi-linear system); • to ensure convergence and to avoid negative roots, the Broyden numerical algorithm should be used; • the programming and numerical solutions by using Matlab; In this stage the system of 26 equations can be reduced in a non-linear equations system of five equations and five unknowns. First approach, let us consider Eq. (43) for all equilibrium constants of Eqs. (2)–(22) and the concentration of all species in Eqs. (26)–(31) can be defined in terms of the concentration of key components in respect to the Eq. (1), namely mHþ , mSO2 , mAsO2 , mAsO3 , mFe3þ . The entire 4 4 problem is reduced to find the roots of this 5  5 system. • Reducing the 5  5 system to 1  1 system (pH method) The pH is obtained by reducing the system equations by one non-linear equation with one unknown, the activity of hydrogen is calculated from the Eq. (31), by substituting the experimental values of the other components of the system (mSO2 , mAsO2 , mAsO3 , mFe3þ ). The 4 4 mathematical transformation of Eq. (31) leads to the expression below, by using KðT Þ values; (Hþ ) can be computed numerically by NewtonÕs method to obtain the correct root. This pH method has three advantages: it helps to find species that really exist in the system, it leads to find mathematical equations having solutions that will further be used in the system equations. It helps also to find initial guesses to be used subsequently when finding the roots of the system equations. 3

2

ðHþ Þ  X ðHþ Þ  K5 ðFe3þ ÞðHþ Þ  K6 ðFe3þ Þ  2½K7 ðFe3þ Þ2  ¼ 0

ð44Þ

435

The Eq. (44) converges by neglecting the fifth term. After finding the pH of the system at different temperatures (Table 3), the distribution activities of the other species can be calculated as functions of ðaHþ Þ, as shown by the following equations:  Iron: ðHþ Þ1 ½ðHþ Þ þ K6 ðHþ Þ1 þ K5 ðFe3þ Þ  ðFe3þ ÞT ¼ 0 ð45Þ  Arsenic: If a3 ¼

½H3 AsO4  ½AsO3 4  ; . . . ; a0 ¼ C C

ð46Þ

Then "

K17 K17 K16 K17 K16 K15 a3 ¼ 1 þ þ þ þ 3 ðH Þ ðHþ Þ2 ðHþ Þ a2 ¼ a3

#1

K17 K17 K16 K15 ; . . . ; a0 ¼ a3 3 ðHþ Þ ðHþ Þ

293.15 373.15 423.15 483.15 523.15 553.15 a b

Activity of Hþa 3

2:40  10 9:55  103 2:50  102 2:50  102 2:50  102 2:50  102

ð48Þ

From the Eqs. (46)–(48) all arsenic species can be derived in respect with pH results  Sulfate: 2 2 þ ðSO2 4 Þ þ K11 ðSO4 ÞðH Þ  ðSO4 ÞT ¼ 0

ð49Þ

• Uncertainty The uncertainty is computed by using the derivative approach proposed by Cabaniss (1999), due to the fact this method is faster then the Monte Carlo simulation. The derivative approach estimates by assuming that:  the input constraints have a Gaussian uncertainty distribution with a known mean (l) and a standard deviation (r);  the equations relating to the calculated concentration of input are sufficiently linear that the second and higher order derivatives can be disregarded. If the input uncertainties are mutually independent, so that cross-terms are negligible, the standard deviation of an output concentration (Sc ) can be calculated from the first partial derivatives of the concentration C with respect to the constraints Zi .

Table 3 Computed pH values of the system Temperature (K)

ð47Þ

pH

Experimental values of pHb

2.62 2.02 1.60 1.60 1.60 1.60

2.5 – – 2.0 1.9 1.82

This calculation is made with the assumption that all arsenic ions react with ferric ions at 423.15. Experimental data are obtained after cooling the hot solution.

436

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

Table 4 Comparison of results obtained by solving the system equation and the pH method 5  5 system equations BroydenÕs method

1  1 system equation Newton method (pH method for initial guesses)

Difference

Temperature (298.15 K) pH a: (Fe3þ )/(Fe3þ )T a: (H3 AsO4 )/(AsO3 4 ) 2 a: (SO2 4 )/(SO4 )T

2.62 0.7 0.308 0.80

2.51 0.73 0.359 0.73

0.11 0.03 0.051 0.07

Temperature (373.15 K) pH a: (Fe3þ )/(Fe3þ )T a: (H3 AsO4 )/(AsO3 4 ) 2 a: (SO2 4 )/(SO4 )T

2.02 0.933 0.707 0.075

2.01 0.866 0.707 0.076

0.01 0.067 0 0.001

Temperature (423.15 K) pH a: (Fe3þ )/(Fe3þ )T a: (H3 AsO4 )/(AsO3 4 ) 2 a: (SO2 4 )/(SO4 )T

1.60 0.94 0.81 0.0037

1.56 0.97 0.86 0.0038

0.04 0.03 0.05 0.0001

Temperature (483.15 K) pH a: (Fe3þ )/(Fe3þ )T a: (H3 AsO4 )/(AsO3 4 ) 2 a: (SO2 4 )/(SO4 )T

1.60 0.94 0.81 0.0037

1.56 0.97 0.86 0.0038

0.04 0.03 0.05 0.0001

Temperature (533.15 K) pH a: (Fe3þ )/(Fe3þ )T a: (H3 AsO4 )/(AsO3 4 ) 2 a: (SO2 4 )/(SO4 )T

1.60 0.94 0.81 0.0037

1.56 0.97 0.86 0.0038

0.04 0.03 0.05 0.0001

Temperature (553.15 K) pH a: (Fe3þ )/(Fe3þ )T a: (H3 AsO4 )/(AsO3 4 ) 2 a: (SO2 4 )/(SO4 )T

1.60 0.94 0.81 0.0037

1.56 0.97 0.86 0.0038

0.04 0.03 0.05 0.0001

s ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2  2  2 oC oC oC 2 2 Sc ¼ rz1 þ rz2 þ r2z3 þ oZ1 oZ2 oZ3 ð50Þ

4. Results and discussion 4.1. Solving the system of non-linear equations (5  5 system equations) When all these equations are established with pH method, a more interesting methods is to set up the system equation made with the Eqs. (44), (45), (47) and 49 by using the numerical method discussed above and the program developed with Matlab. From the approximate solution of pH method (1  1 system) developed above, we can obtain the initial guesses to be applied for solving the 5  5 system equation. However, from this stage even the Newton method may converge. The results obtained are compared with pH method and

are showed by Table 4. Consequently, we may state the following: • Even the BrodenÕs method gives the correct solution of the model, it can be noticed that the results obtained with the pH method are of about the same magnitude. • Above 423.15 K no change is noticed in the aqueous phase of the system, the reaction is occurred and the pH of the system is constant. 4.2. pH of the system Fig. 1 shows the computed and the experimental values of pH of the system at various temperatures. Results indicate that the pH of the system decreases with the increase of temperatures. It is noticeable that no significant change in pH when temperatures of the system increase above 473.15 K. The comparison of model prediction and experimental values has shown that there is a good agreement of the model at low temperature.

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

437

4.4. Arsenic species in the system

3 2.5

4.3. Ferric species in the system

4.5. Sulfate species in the system

In agreement with the Eq. (44), the revised iron mass balance will include only species that really exist in the 2þ system. However, the FeHAsOþ are 4 and FeH2 AsO4 neglected due to lack of thermodynamic data, the only data for these two species ever published are for the DG at 298.15 K by Robins (1990). The lack of DH and DS makes the extrapolation impossible to calculate K values at elevated temperatures. By using the experimental data of other components and the pH computed above, the speciation of iron can be established by using the Eq. (45). Fig. 2 shows the ferric ions speciation at different temperatures. The predominance of Fe3þ over the pH range of reaction can be noticed. When the temperature increases and the activity of ferric ions decreases due to the reaction, the ferric ions (Fe3þ ) increases due to the decrease of pH. This calculation is in full agreement with the previous study conducted by Welham et al. (2000b).

In agreement with Eqs. (44) and (30), sulfate species in solution can be derived with the Eq. (49). Fig. 4 shows the calculated species of sulfate mass balance at different temperatures. At 298.15 K, (SO2 4 ) ions are predominant in the system. The equivalence point where the sulfate

pH

1.5 Series1

1

Series2 0.5 0 273.15

323.15

373.15

423.15 473.15

523.15

573.15

Temperature (K)

Fig. 1. Calculated and experimental pH at different temperatures (series 1: pH computed, series 2: pH measured).

distribution coefficient

This observation is in full agreement with the uncertainty calculation (Eq. (50)). The outputs obtained at 298.15 K have uncertainty of 1–3% and at high temperatures the outputs become extremely sensitive up to 10%.

Many researchers have investigated the oxidation reaction of arsenic (III) to arsenic (V). Pettine and Millero (2000) give the theoretical ratio [As(V)/As(III)]. From this ratio it can be noticed a predominance of the oxidized species and negligible concentrations of the reduced form in oxic conditions. Therefore, it is assumed that all As(III) is oxidized to As(V); subsequently only the speciation of arsenic (V) will be accounted in the modeling work. From the Eqs. (46)–(48) all arsenic species can be derived in respect with pH results. Fig. 3 shows the chemical speciation of arsenic at different temperatures. From Fig. 3, it can be noticed that at temperature 293.15 K two species predominate H3 AsO4 and H2 AsOþ 4 in the system. When the temperatures increase and the reaction advances, the pH decreases, therefore, the activity of H3 AsO4 decreases and H2 AsO 4 becomes the major species in the system.

2

1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 273.15

alpha3 alpha2 alpha1 alpha0

323.15

373.15

423.15

473.15

523.15

573.15

Temperature (K)

1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 273.15 323.15 373.15 423.15 473.15 523.15 573.15

Fe3+ Fe(OH)2+ Fe(OH)2(+)

Temperature (K)

Fig. 2. Calculated distribution of ferric species at different temperatures (series 1: (Fe3þ ), series 2: (FeOH2þ ) and series 3: (Fe(OH)þ 2 )).

distribution coefficient

distribution coefficient

Fig. 3. Calculated distribution of arsenic species at different temper2 3 atures (a3 : H3 AsO4 , a2 : H2 AsO 4 , a1 : HAsO4 and a0 : AsO4 ).

1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 273.15

alpha alpha1

323.15

373.15

423.15

473.15

523.15

573.15

Temperature (K)

Fig. 4. Calculated distribution of sulfate species at different tempera2  2 tures (a ¼ ðSO2 4 Þ=ðSO4 ÞT , a1 ¼ ðHSO4 Þ=ðSO4 ÞT ).

438

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

and bisulfate concentrations are equal is obtained at temperature of 323.15 K. At temperatures above 373.15 K, (HSO 4 ) become the major specie in the system and subsequently the pH of the system is less than 1.75.

temperatures have the same effect to both solutions and the precipitation of FeAsO4 occurs at room temperature.

4.6. Temperature versus pH diagram

Let us now look the phase and stability diagrams of the ferric–arsenic system. The stability fields of phases can be quickly checked in term of Eh–pH diagram. From the Eh–pH diagram of As–Fe–H2 O system, it can be noticed that H3 AsO4 is more stable phase in an oxidizing region, at low temperatures and acidic pH. Oppositely, ferric arsenate becomes the more stable phase in the same region at temperatures above 423.15 K. The phases analysis for this system as experimentally determined (Mambote et al., 2001) as shown that in low pH solutions precipitating above 483.15 K, the solid phase identified was a well crystallized dark green, unknown ferric arsenate phase. At lower temperature below 483.15 K and at low pH as well, the solid phase was a pale weak yellow amorphous ferric arsenate. The experimental results are in disagreement with the predicted phase diagram and the T –pH graph due to the fact both graphs do not account the in-situ oxidation of As (III) to As (V).

The pH versus temperature graph was computed by using the speciation model established above for the solution of 3  102 mol/l of Fe2 (SO4 )3 and 2  102 mol/l of AsO3 4 both alone. The results of these calculations are shown by Figs. 5 and 6. It can be noticed for ferric sulfate solution, the shift of solution equilibrium to lower pH with an increase in temperature. The pH of the solution decreases with the increase of temperatures. The optimal reaction conditions are reached when temperature is above 423.15 K and pH is less then 2.5. The intersection point can be considered as an optimal reaction temperature. This finding is in full agreement with the experimental results obtained by Monhemius and Swash (1999). This can be understood due to the fact this graph consider only the reaction between arsenic (V) and iron (III). There is a disagreement by using such diagram to predict FeAsO4 precipitation. Fig. 5 does not give information for the temperature at which precipitation occurs in the system rather than the optimal reaction temperature. There are two reasons why the T –pH graph does not give information about the thermal precipitation prediction:

8.5 7.5 6.5 pH

5.5 pH -As

4.5

pH-Fe

3.5 2.5 1.5 0.5 273.15

323.15

373.15

423.15

473.15

523.15

573.15

Temperatures (K)

4.7. Phase diagrams

4.8. Can speciation model help us in understanding crystallization in the system? From the experimental work conducted by Mambote et al. (2001), it was shown that a stable crystalline ferric arsenate phase is not synthesized at temperature below 483.15 K from liquors containing Fe(III) and As(III), through in situ oxidation. A question arises if a thermodynamic modeling can help to predict the temperature at which transformation of amorphous to crystalline occurs in the system. The development model was used above to simulate the speciation in aqueous solutions which contain arsenic (III) and ferric sulfate (III) in the 298.15–553.15 K temperatures range, 0–20 g/l As(III) and 0–20 g/l Fe(III) concentration ranges. From the results discussed above it can noticed that:

Fig. 5. Temperature versus pH diagram.

8.5

pH

7.5 6.5

pH -As

5.5

pH-Fe pH-As -

4.5

pH-Fe-

3.5

pH-As+

2.5

pH-Fe+

1.5 0.5 273.15

323.15

373.15

423.15

473.15

523.15

573.15

Temperature (K)

Fig. 6. Temperature versus pH diagram (extrapolation including inaccuracies, þ: plus uncertainty, ): minus uncertainty).

• The ferric arsenate phase is formed at pH less then 2.5, this result is in good agreement with experimental results. • The distribution and activities of dissolved species of the system as a function of temperatures is obtained by the outputs of the model. • The model shows the temperature effect on the composition and the pH of the system. The activity of hydrogen ions increases with the increase of temperature. • From the model developed, the supersaturation index can be calculated in order to determine the precipitation of ferric arsenate in the system. In our case, the precipitation of As(V) and Fe(III) occurs at tempera-

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

ture above 298.15 K, therefore, this index does not give additional information leading to understand the system. • The speciation model helps also to calculate the fugacity of various gases present in the system. The fugacities calculated in this way are those found in a gas phase which is in equilibrium with the system. • The model gives also information to compute the Eh due to the fact that the redox equilibrium is involved in the system. The model described above does not predict the transformation process of precipitation in the current system. The transformation is governed by the rate of kinetics. When a solid phase forms, a number of crystalline modifications such as polymorphs, hydrates or other solvates, a mixture of forms may be produced and, over a period of time, the thermodynamically less stable forms will transform to the stable form, known as OstwaldÕs law of stages.

5. Conclusions In this work it has been evaluated the thermodynamic modeling of the ferric–arsenic(III)–sulfate system in aqueous solution under elevated temperature and pressure. The formulations are based on the chemical speciation model and the solution model is found by numerical method. In aqueous ferric–arsenic–sulfate system, arsenic is 3 distributed as soluble species: AsO 2 , AsO4 , H3 AsO4 ,  2 H2 AsO4 , HAsO4 , etc. The predominant species depend strongly on pH and weakly dependent on concentration. It is important to note that even theoretically, it is obvious that ferric ions reacts with arsenate to form a series of þ complex ions (FeH2 AsO2þ 4 , FeHAsO4 , FeAsO4 (aq.) and 3 Fe(AsO4 ) ), but little information is given in the literature and their existent were not been confirmed by spectrometric measurements. Despite, the model predict the existent of FeH2 AsO2þ and FeHAsOþ 4 4 as plausible complexes to be formed at specific pH conditions of the system. These complexes were not accounted in the modeling process due to lack of thermodynamic data. The solution of the model revealed the following features: • The optimal precipitation conditions in order to precipitate the ferric arsenate are in pH below 2.5. • The predominant species depend on the pH and also the type and amount of ligands present in the solution. At low pH and high temperature, the metal (Fe3þ ) becoming decreasingly hydroxylated and the 2 anions (AsO3 4 , SO4 ) protonated. • The temperature has an effect on pH, when the temperature increases, when the pH decreases subsequently the precipitation of ferric arsenate is favored.

439

Therefore, the hydrothermal precipitation of the ferric–arsenic–sulfate system can be understood in thermodynamic terms as a conjugation of two separate effects: • The shift of the solution–solid equilibrium interface to lower pH with the increase in temperature. • The change in pH of the solution as the temperature is increased. The results of the computational performed in this work were confirmed by pH measurements, phases analysis and by an extensive review of published work. It was found a good agreement between predicted and measured results. However, the optimal reaction temperature determined by experiment is higher (>483.15 K) than the optimal reaction temperature predicted by the thermodynamic model (423.15 K). It is due to the fact the T –pH graph does account the oxidation equilibrium reaction of As(III) to As(V) in the system. Therefore, the full picture of the system can be obtained only if both thermodynamic and kinetic models are used. Although the chemical speciation model is the preferred thermodynamic modeling approach used in the literature, it is also the most difficult modeling techniques to perform, especially at elevated temperatures. The thermodynamic database must be well selected, since inappropriate database selection will significantly alter the result. The accuracy of the solution model will depend on the accuracy of K and the extrapolation of K at high temperatures. An attempt to evaluate the extrapolation models to predict K values at elevated temperatures has been considered, it was noticed that the HelgesonÕs method has given more reliable values. This can be understood due to the fact this method covers more properties over a wider range of temperatures and pressures.

Acknowledgements The authors wish to acknowledge the contribution of Professor Dr. J. Monhemuis, Dr. P. Swash, Dr N. Welham, Dr. R. Van der Weijden and F. Henry.

References Brassard, P., Bodurtha, P., 2000. A feasible set for chemical speciation problems. Computers and Geosciences 26 (3), 277–291. Butler, J.N., 1964. Ionic Equilibrium, a Mathematical Approach, Addison-Wesley, London. Cabaniss, E.S., 1999. Uncertainty propagation in geochemical calculations: non-linearity in solubility equilibria. Applied Geochemistry 14, 255–262.

440

R.C.M. Mambote et al. / Minerals Engineering 16 (2003) 429–440

Casas, J.M., Alvarez, F., Cifuentes, L., 2000. Aqueous speciation of sulfuric acid-cupric sulfate solutions. Chemical Engineering Sciences 55, 6223–6234. Dennis Jr, J.E., Schnabel, R.B., 1996. Numerical Methods for Unconstrained Optimization and Non-linear Equations. PrenticeHall, Englewood Clifs, Philadelphia. Fillippou, D., Demoplous, G.P., Papangelakis, V.G., 1995. Hydrogen ion activities in mixed metal sulfate aqueous systems. AICHE Journal 41 (1), 171–184. Garrels, R.M., Christ, C.L., 1965. Solutions, Minerals and Equilibria. Freeman, Cooper and Company, San Francisco. Helgeson, H.C., 1967. Thermodynamics of complex dissociation in aqueous solutions at elevated temperatures. Journal of Physical Chemistry 71, 3121–3136. Holstad, A., 1999. Numerical solution of nonlinear equations in chemical speciation calculations. Computational geosciences 3, 229–257. Koukkari, P., Laukkanen, I., Liukkonen, S., 1997. Combination of overall reaction rate with Gibbs energy minimization. Fluid Phase Equilibria 136 (1–2), 345–362. Loebenstein, J.R., 1994. The materials flow of arsenic in the United States. Bureau of Mines Information Circular, US department of the interior, pp. 1–9. Lupis, C.H.P., 1983. Chemical Thermodynamics of Materials. NorthHolland, New York. Mambote, R.C.M., Reuter, M.A., van Sandwijk, A., Krijgsman, P., 2001. Immobilization of arsenic in crystalline form from aqueous solution by hydrothermal processing above 483.15 K. Minerals Engineering 14 (4), 391–403. Martell, A.E., Smith, R.M., 1982. Critical Stability Constants, vol. 5. Plenum Press, New York. Monhemius, A.J., Swash, P.M., 1999. Removing and stabilizing As from copper refining circuits by hydrothermal processing. Journal of Minerals Metals and Materials Society 51 (9), 30–33. Penny, J., Lindfield, G., 2000. Numerical Methods Using Matlab, second ed. Ellis Horwood, New York.

Pettine, M., Millero, F.J., 2000. Effect of metals on the oxidation of As(III) with H2 O2 . Marine Chemistry 70 (1–3), 223–234. Pitzer, K.S., 1991. Ion interaction approach: theory and data correlation. In: Pitzer, K.S. (Ed.), Activity coefficients in electrolyte solutions, second ed.. CRC press, Boca Raton, pp. 75–153. Rafal, M., Berthold, J.W., Scrivner, N.C., Grise, S.L., 1994. Models for electrolyte solutions. In: Sandler, S.I. (Ed.), Models for Thermodynamic and Phase Equilibria calculations. Marcel Dekker, New York, pp. 601–670. Riveros, P.A., Dutrizac, J.E., Spencer, P., 2001. Arsenic disposal practices in metallurgical industry. Canadian Metallurgical Quartely 40 (4), 395–420. Robins, R.G., 1981. The solubility of metals arsenates. Metallurgical Transactions B 12B, 103–109. Robins, R.G., 1990. The stability and solubility of ferric arsenate: an update. In: Gaskell, D.R. (Ed.), EPDC Congress Õ90, TMS. pp. 93– 104. Sandler, S.I., 1994. Models for Thermodynamic and Phase Equilibria calculations. Marcel Dekker, New York. Stumm, W., Morgan, J.J., 1996. Aquatic Chemistry, third ed. WileyInterscience, New York. Sverjensky, D.A., Shock, E.L., Helgeson, H.C., 1997. Prediction of the thermodynamic properties of aqueous metal complexes to 1000 C and 5 kb. Geochmica et Cosmochimica Acta 61 (7), 1359–1412. VasilÕev, V.P., 1962. Influence of ionic strength on the stability constants of complexes. Russian Journal of Inorganic Chemistry 7, 924–927. Welham, N.J., Malatt, K.A., Vukcevic, S., 2000a. The stability of iron phases presently used for disposal from metallurgical systems––a review. Minerals Engineering 13 (8–9), 911–931. Welham, N.J., Malatt, K.A., Vukcevic, S., 2000b. The effect of solution speciation on iron–sulphur–arsenic–chloride systems at 298 K. Hydrometallurgy 57, 209–223. Zemaitis, J.F., Clark, D.M., Rafal, M., Scrivner, N.C., 1986. Handbook of Aqueous Electrolyte Thermodynamics. American Instute of Chemical Engineers, New York.