Calculation of Component Activities for Se-Based Binary Systems by NRTL Equation

The interaction energy parameters of the Nonrandom two-liquid (NRTL) equation for the Bi–Se, Pb–Se, Sb–Se, and Tl–Se binary systems were determined. The difference in the goodness of fit between the systems was discussed for each system. The activities of components at corresponding temperatures for the four binary systems were calculated and discussed using the NRTL equation. The mixing enthalpies and Gibbs free energies of these systems were obtained. The calculations for the four binary systems are in agreement with the experimental data from the literature. It shows that the NRTL equation is reasonable and accurate based on the concept of local composition and two-fluid theory, especially for binary systems containing liquid–liquid equilibria. The results we obtained using the NRTL equation can be very useful in performing vacuum experiments for crude selenium.

To produce high-grade selenium in a more economical and environmentally friendly way, Zha applied a vacuum distillation method to purify crude selenium [8]. Shiryaev removed the Barium impurities from Selenium by vacuum distillation, and studied the behavior of barium impurities during evaporation process [9]. Vacuum distillation is a successful and rapidly expanding method to separate various elements from alloys, with the advantage of simplified process, good environmental protection, low operation costs, and simple equipment. Thermodynamic properties of Se-based binary systems should guide the optimum condition of the vacuum distillation process.
It is necessary to study the thermodynamic properties of melts both in metallurgy and designing new materials. Due to the complexity of high-temperature experiments and the limitation of precision, experimental measurement is difficult, and the thermodynamic data obtained from the experiment often have deviations. Therefore, theoretical studying is a significant and effective approach to obtain thermodynamic properties of alloys in metallurgical processes. General expressions for partial Gibbs energies are derived from the Redlich-Kister polynomial for substitutional and interstitial solutions in CALPHAD [10]. Since most systems are non-ideal, a modified concentration (activity) should be applied to replace the actual concentration to analyze the thermodynamic behavior of solutions accurately. Activity of component has become one significant research topic in the field of thermodynamic properties [11][12][13].
Selenium is a nonmetallic element, a binary system containing selenium and one metal is unique comparing with other binary systems of metals. As shown in Fig. 1, there are stable compounds and unstable compounds in the solid-liquid phase diagrams of Selenium-based binary systems [14]. In addition, Se-based systems usually contain liquid-liquid miscibility phase regions. In some systems, the miscibility phase region even accounts for one third of the concentration range of the system, making it difficult to accurately measure the component activities experimentally for the Se-based systems.
Since the complexity of the systems and the difficulty of high-temperature experiment, therefore, it is necessary to correlate the activity in the whole concentration range through thermodynamic model.
Wilson proposed an expression based on the concept of local composition for activity coefficients, which has been used successfully to describe the thermodynamic properties of multicomponent mixtures with only binary data [15]. For most binary alloy systems, the Wilson equation can accurately calculate the activities of components. In our previous work, we precisely calculated activities of components for the Pb-Sn and Pb-Ag binary systems and correlated the activities of components for the Pb-based and Snbased binary and ternary systems using the Wilson equation, which agrees with experimental data [16][17][18]. However, the Wilson equation cannot represent liquid-liquid equilibria. In 1968, based on the concept of local composition and two-fluid theory, Renon and Prausnitz proposed a new equation with similar properties to the Wilson equation, the NRTL (Nonrandom two-liquid) equation [19]. Due to the differences in molecular interactions, the mixing is not entirely random, so the local composition differs from the overall composition. The NRTL equation is widely used to predict the activities and phase equilibrium in binary and multicomponent systems including liquid-liquid equilibria [20,21]. Zhang et al. applied the NRTL equation to calculate the activity coefficient and other thermodynamic properties of the Zn-Pb system [22] and got satisfactory results.
This work aims to show the application of the NRTL equation to calculate component activities of the Bi-Se, Pb-Se, Sb-Se, and Tl-Se binary systems. The accuracy of the NRTL equation in calculating the activities of compounds for binary systems containing liquid-liquid immiscibility was discussed. The obtained activities are important thermodynamic data for separation and purification of crude selenium in vacuum distillation, especially for the process of ordinary distillation, stage batch distillation, or even pressure swing distillation in the vacuum furnace.

Redlich-Kister (R-K) Polynomial
The thermodynamic properties of binary systems in the liquid state, such as Gibbs free energy, enthalpy, and entropy, could be described by the Redlich-Kister polynomial commonly [10].  [13]. a Bi-Se system, b Pb-Se system, c Sb-Se system, d Se-Tl system Based on the Redlich-Kister polynomial, the molar excess Gibbs free energy G E can be expressed: The partial molar excess Gibbs free energy functions G E 1 and G E 2 are as follows: The relation between the activity coefficient i and the partial molar excess Gibbs free energy is defined as follows: where T is the temperature of the binary system, R represents the ideal gas constant, respectively.

NRTL Equation
Considering the non-randomness of mixing in the NRTL equation, the relation between the local mole fractions x 21 and x 11 is given for a 1-2 binary mixture [19]: where α 12 is a constant that represents the non-randomness of the mixture (α 12 = α 21 ). T represents the temperature of the binary system, R the ideal gas constant, respectively. g 21 and g 11 are interaction energies between a 1-2 and 1-1 pair of molecules (g 12 = g 21 ), respectively.
Interchanging subscripts 1 and 2, Eq. 5 becomes The local mole fractions have the following relations: Rearrangement of Eqs. 5-8, the local mole fractions may be solved as follows: Based on the two-liquid theory, Renon assumed that there are two types of cells in a binary mixture-One centered on molecule 1 and one centered on molecule 2 (as shown in Fig. 2) [19]. Each center molecule is surrounded by assortments of the same molecules, with each surrounding molecules surrounded similarly, and so on.
The whole system can be equivalent to a virtual mixture formed by mixing these two kinds of cells. For the first type, the excess Gibbs free energy for the cell is given by For the second type, the excess Gibbs free energy for the cell is given by The molar excess Gibbs energy is dominated by the excess Gibbs free energy of transferring x 1 molecules from a cell of pure liquid 1 into cell 1 of the solution and that of transferring x 2 molecules from a cell of pure liquid 2 into cell 2 of the solution. Therefore, The residual Gibbs free energy of the binary system is given by Substituting Eqs. 9-13 into Eq. 14 produces the expression of molar excess Gibbs free energy of a binary mixture, which is easily generalized to solutions containing multiple components: The nonrandom factor α is assumed to be independent of temperature and composition, a structural factor of solution. The value of α should be carefully adjusted for different kinds of solutions. At present, α parameter should be determined according to the general rules proposed in the literature [23], or fitted to the experimental data.
.  The activity coefficients of components for a binary system are obtained by differentiation as follows: With the reasonable selection of the adjustable parameter α, the calculation by the NRTL equation is in agreement with the experimental data for highly non-ideal mixtures, especially for partially miscible systems. A rigorous evaluation of the derivation for the NRTL equation shows that the equation is more suitable for excess mixing enthalpy than for excess Gibbs free energy.
The relation between the excess Gibbs free energy and the excess mixing enthalpy is as follows: By associating NRTL parameters with excess mixing enthalpy, we have For comparing with experimental data, the average relative error is expressed by and the average standard error is given by where a i,exp is the experimental activity and a i,cal is the calculated result using the NRTL equation, respectively. m is the number of experimental data.

Results and Discussion
The parameters 12 and 21 of binary systems can be determined by fitting the experimental data [24,25] through Eqs. 18 and 19. An objective function was chosen to calculating the optimal 12 and 21 . The objective function is (18) The α parameter is associated to the non-randomness of the mixture. Considering the inherent requirements of the α parameter, to get an idea of quality of the correlations and to choose the binary interaction parameters that lead to a better agreement between experimental and calculated data through Eq. 23. This work, α parameter was set to different values between 0.01 and 0.5 during parameter optimization so that the best results are presented in Table 1. Some researchers suggest that α is related to the coordination number of the system [26,27]. For the structure of a melt, we can take a unit volume of the melt to analyze its approximate coordination number.
Substituting the corresponding parameters α, τ 12, and τ 21 into Eqs. 18 and 19, the activities of components for the Bi-Se, Pb-Se, Sb-Se, and Tl-Se binary systems are obtained. From Table 2, the average relative deviations of obtained activities for four Se-based systems by the NRTL model are reasonable. For Bi-Se, Sb-Se, and Tl-Se systems, the average standard deviations are less than 10%. Obviously, the NRTL equation can calculate the activities of these three systems more accurately. For the Pb-Se system, the average standard deviation of activity by the NRTL equations is very large, since at the concentration range of x Se = 0-0.25, the activity values of component selenium are less than 10 -3 , a small change in the calculated value will cause a large deviation.

Bi-Se System
As shown in Fig. 3, the calculated activities of component selenium are satisfactory using the NRTL equation, especially precisely describing the transformation from negative deviation to positive deviation at x Se = 0.7-1. Compared with the ideal solution, the activities of the binary Bi-Se system present an asymmetrical distribution.
Assuming that the pair-potential energy interaction parameters (g 12 -g 22 ) and (g 21 -g 11 ) are independent of temperature, the values of τ 12 and τ 21 at corresponding temperatures can be obtained from Eq. 16 if the values of τ 12 and τ 21 at a given temperature are known. The parameters τ 12 and τ 21 at 882 and 998 K are shown in Table 3. Substituting the corresponding parameters α, τ 12, τ 21, τ 12′ , and τ 21′ into Eq. 21, the enthalpies of mixing for the Bi-Se system at 882 and 998 K are obtained. From Fig. 4, the calculated values are consistent with the experimental values from the literature [28] at 998 K and show the same trend as the experimental result at 882 K, which illustrates the accuracy of the NRTL equation in activity calculation. There is a V-shaped geometry and a minimum value, which indicates the formation of a stable melt structure around x Se = 0.6 at 998 K. Therefore, the interaction between bismuth and selenium atoms will slow down the separation of alloying elements during vacuum distillation. This feature is in agreement with the fact that the existence of Bi 2 Se 3 .

Pb-Se System
Due to the lack of sufficient experimental data for the Bi-Se, Sb-Se, and Tl-Se binary systems, only the activities of the Pb-Se binary system were calculated by the R-K polynomial. As can be seen from Fig. 5, the calculations of the activities of component Se by the two models are in agreement with the experiment [25]. For the activities of component Pb, at the Fig. 3 Comparison of the calculated activities by the NRTL equation (lines) with the experimental data [24] (symbols) for the Bi-Se system at 994 K Compared with the ideal solution, the negative deviation of the binary Pb-Se system is pretty strong and presents an asymmetrical distribution. This is mainly due to the significant difference in the properties of lead and selenium atoms, such as electron density, electronegativity, and atomic volume.
The calculations and experimental data of Gibbs free energies for the Pb-Se system at 1360 K are given in Fig. 6a. The calculated Gibbs free energies and partial Gibbs free energies of selenium show no difference from the experimental ones; in contrast, the calculated partial Gibbs free energies of lead show a significant deviation from experimental data at x Se = 0.5-1. This phenomenon may be attributed to the prediction error of the activities of lead and the presence of intermediate compound PbSe . From Fig. 6b, the mixing enthalpies calculated in the present work are close to those reported by Lin [25]. It can be seen clearly that based on the good results of the activities of components, the Gibbs free energies and the mixing enthalpies should be reasonable and reliable.

Sb-Se System
For the Sb-Se binary system, the activities of antimony show a strong negative deviation in the whole concentration range from Fig. 7, while the activities of selenium present a strong negative deviation in the range of x Se = 0-0.8, and a weak positive deviation at x Se from 0.8 to 1. The calculated values are in agreement with the experimental data [24].
The calculated enthalpies of mixing of the Sb-Se system at 935 K are given in Fig. 8. Compared with the experimental values from the literature [28], a V-shaped geometry and a minimum value in Fig. 8 are evident, which indicates a strong interaction between antimony and selenium atoms and the formation of a stable melt structure around x Se = 0.6 at 935 K. This feature is in agreement with the existence of the compound Sb 2 Se 3 . However, there is an obvious deviation between the calculated enthalpies of mixing and the experimental values. The reason is that the lack of activities in the antimony-rich region, which causes deviation of the enthalpies for the Sb-Se system. Fig. 6 Comparison of the calculated partial Gibbs free energy using the NRTL equation (a) and mixing enthalpy (b) with experimental data [25] for the Pb-Se system at 1360 K  Fig. 7 Comparison of the calculated activities by the NRTL equation (lines) with the experimental data [24] (symbols) for the Sb-Se system at 994 K   Fig. 8 Comparison of the calculated excess mixing enthalpy (lines) using the NRTL equation with experimental data [27] (symbols) for the Sb-Se system at 935 K 1 3

Tl-Se System
In Fig. 9, the NRTL equation presents a precise calculation of the activities for the Tl-Se binary system, especially accurately describing the transformation from negative deviation to positive deviation at x Se = 0.7-1.
The calculated and experimental excess mixing enthalpies of the Tl-Se system at 735 K are given in Fig. 10. The calculations are lower than the values reported in [28]. This Fig. 9 Comparison of the calculated activities by the NRTL equation (lines) with the experimental data [24] (symbols) for the Tl-Se system at 860 K Fig. 10 Comparison of the calculated excess mixing enthalpy (lines) using the NRTL equation with experimental data [27] (symbols) for the Tl-Se system at 753 K phenomenon is attributed to two reasons. One is the lack of experimental activities of thallium, and the other is the existence of TlSe and Tl 2 Se, which will affect the calculated excess enthalpy of mixing of the binary system.

Conclusions
The Se-based binary systems show a strong negative deviation, which indicates a strong interaction between the selenium atom and other atoms. It is consistent with the fact that there are compounds in these binary systems. The component activity and excess enthalpy data for the Pb-Se, Bi-Se, Sb-Se, and Tl-Se binary systems are quite successful in simultaneously correlating by the NRTL equation. The NRTL equation can accurately describe the transformation from negative deviation to positive deviation for the system, especially in the region containing two liquid phases. To obtain reliable results, we should correctly apply this equation to particular classes of systems. This work will provide reliable thermodynamic data for the separation and purification of crude selenium in vacuum distillation.