Description of critical behaviour of model systems using non-Gaussian measures ( Yukhnovskii ’ s approach )

We present a brief review of the works devoted to a study of critical phenomena in three-dimensional model systems dwelling on the Yukhnovskii approach in detail. This approach which is based on the use of non-Gaussian measures allows one to obtain both universal and non-universal quantities. In order to illustrate the advantages of the approach proposed by I.R.Yukhnovskii we apply it to a study of non-universal quantities, namely: (1) the phase transition temperature of a 3D one-component lattice model, (2) the gas-liquid critical point properties of fluid systems.


Introduction
A description of phase transitions as well as critical phenomena connected with them remains a relevant problem.By now a number of theoretical approaches in this field have been proposed and the original schemes of numerical calculations of critical properties of model systems have been developed.A new stage in the development of the phase transition theory was related to the hypothesis of universality put by L.Kadanoff [1].The idea of an ǫ-expansion which appeared to be effective for the calculation of universal characteristics of statistical systems in the vicinity of their phase transition points was proposed in [2].This idea occupied a significant place in the theory of critical phenomena.Considerable progress towards the calculation of universal quantities was also made due to the technique of the resummation of asymptotic series proposed in [3][4][5].Among the approaches devoted to the study of critical phenomena one should distinguish works [6,7] which were not connected with the use of an ǫ-expansion.
One of the relevant approaches to the theoretical description of phase transitions was proposed by K.Wilson in [8] and was developed in his further paper, namely in [9].In these works the Kadanoff idea [1] was embodied in the concrete mathematic formulas.As the result of these investigations the calculation scheme of both critical exponents and critical amplitude ratios near the phase transition point was proposed.The theory was constructed by means of the moments of the Gaussian distribution.As is known these moments tend to the infinity when the system approaches the critical point.Consequently, the approach by K.Wilson did not allow one to obtain the explicit expressions for non-universal characteristics of the model systems under consideration.However, one can get reliable results for many universal quantities using the effective methods of the asymptotic series resummation (see [10]).
In this work we present some results of the latest calculations of non-universal characteristics for statistical models of the phase transition in the critical region using the approach proposed by I.Yukhnovskii [11].This approach similar to Wilson's work is based on the hypothesis of universality [1].The distinction between these approaches consists in using different basic distribution types.The Wilson method is based on the use of the Gaussian basic measures and therefore, results in the emerging of diverging diagrams near the phase transition point [12].In contrast to Wilson's works the Yukhnovskii approach is based on the use of the non-Gaussian measures.As a result, there arise no diverging diagrams in the calculating of the critical behaviour of three dimensional models.This approach allows one to calculate both the universal and non-universal quantities for a number of phase transition models.Among the systems investigated by the Yukhnovskii method is the Ising model [13][14], the n-component model [15][16][17][18], the hierarchical model [19], the cluster ferroelectric model [20], the fluids in the vicinity of their critical points [21][22][23] and others.The construction of the phase transition theory for each of these models has some special features.
In contrast to lattice systems, the description of phase transitions in continuous systems has a number of the important peculiarities.On the one hand, as one usually does in the liquid state theory, we should distinguish a reference system describing the behaviour at short distances.This will allow us to take into consideration the short-range and long-range interactions simultaneously.On the other hand, the grand canonical ensemble (GCE) should be used in order to describe the processes relating to the phase transitions in multi-component fluids in which the composition fluctuations play a crucial role (e.c., the gas-gas and liquid-liquid equilibria in binary fluid mixtures).The task of the development of the CV method for the case of the GCE is also caused by the problem of a selection of the CV phase space which includes the variable connected with the order parameter.Therefore, it is necessary to introduce the GCE in the CV method in the studies of the gas-liquid critical point in a one-component fluid.
We present herein the main aspects of the CV method with a reference system for a multi-component continuous system in the GCE as well as some results for simple and binary fluids obtained using the Yukhnovskii approach.
Yukhnovskii's approach is the synthesis of several basic constituents.First this is a choice of a phase space in which the system is described by means of a certain type of collective variables (CV) [11].For a magnetic system the CV are the variables connected with spin density fluctuation modes; for a ferroelectric they are connected with cluster state fluctuation modes; for a charged-particles system, with generalized charge fluctuation modes; for a binary alloy, with modes of one-particle distribution function; for a one-component fluid, with particle density fluctuation modes and etc.The description of phase transitions is connected with collective effects.Choosing the collective variables specific for a certain physical model we obtain a set of variables the averaged values of which are related to the order parameters.In this approach it is not necessary to introduce the quantities from the "outside" for a description of the ordering in the system.The following constituent of Yukhnovskii's approach is justifying and then using the non-Gaussian density measures as the basic ones [12].This feature distinguishes this approach from the widely known methods of the phase transition description based on the use of the Gaussian moments.It allows one to obtain the explicit analytical expressions for main thermodynamic functions near the phase transition point as the final result.The following and basic constituent of Yukhnovskii's approach is the way of calculating the partition function near the phase transition point.Although this original calculation method [24] as well as Wilson's approach exploit the renormalization group (RG) ideas, it is based on the use of the non-Gaussian measures.This allows one to obtain a qualitatively new form of the recurrence relations (RR) between the coefficients of the block Hamiltonians.
In the limiting case (corresponding to the Gaussian basic density measure) these RR reduce to the Wilson RR [9].As was shown in [25], while this limiting case does not allow one to perform the calculation of the expression for the free energy of the system under consideration, it provides reliable results for the critical exponents of thermodynamic characteristics.
Hereafter we apply the Yukhnovskii approach to the description of non-universal quantities, namely: (1) the phase transition temperature for a 3D one-component lattice model, (2) the gas-liquid critical point properties both of a one-component fluid and a binary mixture.

A partition function in the CV representation for a one-component lattice model
Let us consider a simple one-component system of spins on a three dimensional crystal lattice with a period c.The Hamiltonian of such a system can be written as The interaction potential Φ(r ij ), where where b is the effective interaction radius of the potential.The Fourier transform of (2) has the form [27] Φ In the case of a simple cubical lattice with a period c we have the Brillouin zone 1 where i = x, y, z, N is the total number of particles (N = N x N y N z ).
In the region of small values of the wave vector k we shall use the parabolic approximation for ( 3) which takes place for k B ′ .The quantity B ′ is found from the condition Φ(k) = 0 and is equal to Let us introduce the parameter S 0 determining the size of interval (B ′ , B] Corresponding to (7) we have S 0 1 and hence b b min , where For most models of statistical physics the quantity b is considerably larger than b min .Let us write a partition function of the model (1) in the CV representation [28] 1 The calculation of (3) is performed using spherical coordinates.For a lattice system the following condition should be satisfied 1 In the spherical coordinates this condition can be rewritten in the form where the module of the wave vector k is changed within the interval (−B, B]2 , J(ρ) is a Jacobian of the transition.Variables ρ k and ρ l are related by the relations We assume that as it is in the case of Ising-like systems.
We shall calculate (9) according to the method proposed in [27].This method is based on the Wilson approach [8] which consists in the layer-by-layer calculation of a partition function by means of the successive exclusion of short-wave-length fluctuations from consideration.This is the realization of the Kadanoff idea of the construction of effective block lattices described in detail in [9].The essence of the calculation of (9) reduces to a certain approximation for the potential Φ(k) by its replacement with a set of constant values Φ n , different for each interval of k ∈ (B n+1 , B n ].Such a procedure allows one to calculate the free energy of the system in the non-Gaussian measure approximation and to describe the critical behaviour without using the traditional methods, particularly, an ǫ-expansion.As a result, within the framework of a unified scheme one can perform the calculation of the full expressions for thermodynamic functions and thus obtain both universal and non-universal characteristics of the model under consideration.

Calculation of the partition function
The long wave-length fluctuations (which correspond to small values of the wave vector k) play an important part in determining critical properties of the model.For this reason the Fourier transform of the potential Φ N (k) can be replaced by its parabolic approximation (5).However, values of the phase transition temperature T c depend on how accurately an effect of Φ(k) is taken into account (especially for large values of the wave vector).Here we cannot restrict our consideration to small values of k.Therefore we use the following approximation for the interaction potential Now we integrate in (9) over variables ρ k for k ∈ B\B 0 .For this we present Z in the form We pass from variables ρ k (in the first term of ( 14)) to variables η k by means of the transformation exp and use the integral representation for delta function The integration in ( 15) and ( 16) is performed over the N ′ variables, where The Brillouin zone corresponding to the crystal lattice with the period c ′ (c ′ = cS 0 ) and with N ′ cites is determined by means of the relations In further calculations we shall use the spherical coordinates system.We obtain from (14), taking into account ( 15) and ( 16) Here the following notation is introduced It should be noted that the transformation (15) was first used in [28].In the following the modification of this transformation proposed in [19] will be used.Following this work we introduce the variable We have the relations for ω k Taking into consideration (20) we can write Let us integrate over variables ρ l .We introduce the notation Taking into account (11) we find It is known that cos(2π ωl ) is the Fourier image of the probability measure δ(x 2 − 1).
According to Martsynkevich's theorem [29], the function exp[P (ω)] (where P (ω) is a polynomial within powers of ω) can be the exact Fourier transformation of a probability measure if P (ω) is the Gaussian type polynomial, namely P (ω) = P 2 (ω) = −a 2 ω 2 .Therefore, function exp(P n (ω)), where P n (ω) is a higher order polynomial, is not an exact Fourier transformation of the probability measure.On the other hand, it is proved that Gaussian type distributions do not allow one to describe correctly the behaviour of three dimensional systems in the critical region.The higher order correlation functions play a crucial role in the vicinity of T c and one should describe the physical processes using non-Gaussian distributions of fluctuations.
For small ω l the Fourier image of the probability measure δ(x 2 − 1) can be presented in the following form: We find for cumulants , and etc.The relation ( 24) is an exact one for small ω l -values.For small ω l -values the Fourier image of the probability measure δ(x 2 − 1) can be written in the form exp(P 2n (ω l )), where is the non-Gaussian polynomial of an infinite order.Let us perform an analytical continuation of P 2n .The coefficients standing at all the powers of ωl in the exponent (24) are negative.Such a function is analytical, finite and it decays rapidly if ω l increases.We choose the function as a Fourier image of some measure.This measure will be considered to correspond to some three-dimensional lattice system with a one-component order parameter.It is easy to verify that such a measure corresponds to a model with an unbounded spin and the Fourier transform of this model is given by (25).Such an approximation took place implicitly in a series of the works, namely in [13][14][15][16].The work [30] is also devoted to similar problems.
We use (25) as a Fourier transform of the measure and write the partition function (18) in the form Here C = exp( 1 2 βΦ 0 N), ωl is found from ( 20)-( 22) and certain type relations for ηk hold.Let us perform some transformation is a Kronecker symbol on the lattice with Brillouin zone B. We assume that .
It allows us to pass to a new (block) lattice with Brillouin zone B 0 (17).This conventional approximation [8,11,15] is the realization of the Kadanoff's idea [1] about the block structure generation in the vicinity of the phase transition point.As a result, where Let us calculate the integral over ω.For this we pass to the ω l -representation We have (26) Integrals of this type were calculated in [27].Expression (26) can be presented in the form3 where and we find for the coefficients a 2n where Since the coefficients a ′ 2n depend on the values of the renormalized cumulants M ′ 2n they are functions of b/c.We have for b/c large enough and a ′ 2n are all zero at n 2. For large values of interaction radius b we can use the function instead of the function f (ω) from (27), where Here g 0 = 1/3, f 0 = 8/45, k 0 = 34/315, the quantity S 0 is presented in (7).We find for the coefficients a ′ 2n where Formulas (29) are equivalent to (28).The existence of the small parameter S −d 0 in (29) (for large values of b/c) makes them more advantageous ones than (28).
We perform the approximate calculation of the coefficients a ′ 2n from (29) taking into consideration that S −d 0 is the small quantity.We have We find for the small S −3 0 -values (it corresponds to the large b/c-values, b/c > 0.5) ) .

Introducing the notation
we have (with an accuracy to O(S −4d 0 )) According to (28) we find The values of the coefficients a 2s can be calculated for each value S 0 and therefore for each b/c (S 0 = π √ 2( 6 π ) 1/3 b/c).In table 1 the values of the coefficients a ′ 2s from (30) as well as their values obtained from (29) are presented.As is seen, the approximate formulas (30) allow us to obtain these coefficients with reasonable accuracy even with the comparatively small values of b/c.Therefore, for all b/c > 0.5 we can use the approximate formulas (30) for the calculation of the coefficients a 2n .
The partition function can be written as where coefficients a ′ 2n are determined by (28) or by (29).For b/c > 0.5 the coefficients a ′ 2n can be calculated by means of the approximate formulas (30).

Calculation of the phase transition temperature
We use the results for the coefficients a ′ 2n obtained above for calculating the phase transition temperature.Using the results of the work [31] we can write the partition function in the ρ 4 -model approximation where The coefficients a ′ 2n obtained in [31] are functions of temperature (see formulas ( 14)-( 16) in [31]).This is due to the calculation method proposed in [31].The constant C ′ also depends on temperature.
Using the calculation method for the coefficients a ′ 2n proposed in this work we obtain the only temperature depending coefficient, namely a 2 .All the rest coefficients with n = 0 and n 2 do not depend on temperature.This result is achieved due to the extension of the collective variable phase space (equation ( 20)) proposed first in [19].
Let us use the values of the initial coefficients a 2 , a ′ 4 from (31) for the calculation of the phase transition temperature.The quantities entering into (31) have the form where a ′ 2 and a ′ 4 are determined from (29) or from (30).Consider the coefficient a 2 in detail.According to ( 12)-( 13) we have For the calculation of the phase transition temperature we use the equation written in [31]: Here r * and u * are the coordinates of the fixed point of the recursion relations between the coefficients a (l) 2n of the two successive block structures where the quantities f 0 and ϕ 0 are calculated in [31].For α ′ = 1.00 we find f 0 = 0.5426, ϕ 0 = 0.6894.
In [31] it is also found that and R (0) = −0.5333.Taking into consideration (32) and (33) we obtain the equation for where Using the above expressions for the coefficients from (34) we obtain the value of temperature β c Φ(0) as a function of b/c.The results of this calculation are presented in figure 1.

Functional representation of a grand partition function of a multi-component continuous system
Let us consider a classical multi-component continuous system of interacting particles consisting of N a 1 particles of species a 1 , N a 2 particles of species a 2 , . . .and N am particles of species a m .The system is in volume V at temperature T .
Let us assume that an interaction in the system has a pairwise additive character.The interaction potential between particle γ at r i and particle δ at r j may be presented as a sum of two terms: where ψ γδ (r) is a potential of a short-range repulsion that can be chosen as an interaction between the two hard spheres σ γγ and σ δδ .φ γδ (r) is an attractive part of the potential which dominates at large distances.
Let us start with a grand partition function where (dΓ) = γ dΓ Nγ , dΓ Nγ = d r γ 1 d r γ 2 . . .d r γ Nγ is an element of the configurational space of the γth species; z γ is the fugacity of the γth species: where N γ is the average number of the γth species.Further consideration of the problem is done in the extended phase space: in the phase space of the Cartesian coordinates of the particles and in the CV phase space.An interaction connected with the repulsion (potential ψ γδ (r)) is considered in the space of the Cartesian coordinates of the particles.We call this multi-component hard-spheres system a reference system (RS).The thermodynamic and structural properties of the RS are assumed to be known.Although it is known that mixtures with only repulsive interactions might undergo a phase transition [32], we assume that in the region of temperatures, concentrations and densities we are interested in, the thermodynamic functions of the RS remain analytic.The interaction connected with an attraction (potential φ γδ (r) ) is considered in the CV space.
Let us introduce the grand partition function of the RS where µ γ 0 is the chemical potential of the γth species in the RS.Then the grand partition function (35) can be written as [22,33]: where Ξ 0 is given in (36).The part of the grand partition function which is defined in the CV phase space has the form of a functional integral: Here, 1) φγδ (k) is a Fourier transform of the attractive potential φ γδ (r).
k y , k z , and n i = 0, ±1, ±2, . . ., the thermodynamic limit L → ∞ is assumed.The function φγδ (k) satisfies the following requirements: φγδ (k) is negative for the small values of k and lim k→∞ φγδ (k) = 0.The behaviour of φ γδ (r) in the region of the core r < σ γδ should be determined from the conditions of optimal separation of the interaction.

A one-component fluid
For the case of a one-component continuous system of the N classical particles equations ( 36)-(41) reduce to the forms [21,23,36]: where µ 0 is the chemical potential of the particle in the RS; Here and is the n-particle structure factor of the RS.One can obtain the following expressions for S n (0, . .., 0) (n 4) [37]: Here the CV ρ k is connected with the density fluctuation modes.

A binary fluid mixture
Now let us consider a two-component system consisting of N a particles of species a and N b particles of species b (γ 1 , γ 2 , . . ., γ n = a, b in (37)-( 41)).Having passed in (38) to CV ρ k and c k by means of the orthogonal linear transformation we obtain for Ξ 1 : Here the following notations are introduced: ρ k and c k are CV connected with the total density fluctuation modes and the relative density (or concentration) fluctuation modes respectively.
Functions µ + 1 and µ − 1 have the forms: and are determined from the equations Functions Ṽ (k), W (k) and Ũ(k) are combinations of Fourier transforms of the initial interaction potentials φγδ (k): where Index i n indicates the number of variables γ k in the cumulant expansion (51).Cumulants M (in) n are expressed as linear combinations of the partial cumulants M γ 1 ...γn (see (41)) and are presented for n 4 in [22].

Gas-liquid critical point: a one-component fluid and a binary fluid mixture
Gas-liquid critical points of both a one-component fluid and a binary symmetrical mixture were studied using the approach proposed for the 3D Ising model.Hereafter we present some results of this investigation.

A one-component fluid
The gas-liquid critical point of a one-component continuous system was studied within the framework of the CV method with a RS in [21,23,[36][37][38].Based on the relations (42)-(44) the expression for the grand partition function in the vicinity of the gas-liquid critical point was obtained [21,37].In this case the form of this expression is similar to (31).But the main difference is the presence of odd powers of CV ρ k in the exponent.In the vicinity of the gas-liquid critical point ρ k is connected with the density fluctuation modes.
Using the layer-by-layer integration method [27] for the calculation of the partition function of a one-component fluid in the vicinity of the gas-liquid critical point the following results concerning the non-universal quantities were obtained [36,38]: • the equation of state for T T c and T T c ; • the equation for the parameters of the critical point; • the chemical potential isotherm and the equation for the density jump below T c .

The microscopic Ginsburg-Landau-Wilson (GLW) Hamiltonian for a binary fluid mixture
In order to derive the partition function of a binary fluid mixture in the vicinity of their phase transition points we should first find the CV connected with the order parameter.To this end we restrict ourselves to the consideration of the Gaussian approximation of the functional integral (45)-(51) by setting n = 2 in (50)-(51).While this approximation yields the classical critical behaviour, it provides the correct qualitative picture of the phenomenon under consideration.In order to determine the phase space of the CV connected with the order parameter we introduce independent collective excitations by diagonalizing the square form by means of the orthogonal transformation: where A(k), B(k), C(k) and D(k) are certain functions of the microscopic parameters, temperature, density and concentration of the system (see appendix B in [35]).
As a result, we obtained the following expression for the square form: The explicit expressions for ε 1 (k) and ε 2 (k) are given in [35].The detailed analysis of the coefficients ε 1 (k) and ε 2 (k) allows us to draw the following conclusions: 1.Only one of the two quantities, ε 1 (k) and ε 2 (k), is critical, no matter whether the system approaches the gas-liquid or mixing-demixing phase transition point.This is branch ε 1 (k).
2. Because ε 1 (k) has the minimum at | k| = 0, the CV connected with the order parameter is the variable η 0 in the case of the gas-liquid critical point, as well as in the case of the mixing-demixing phase transition.The particular form of η 0 for each of these phenomena can be determined by means of the relations between the microscopic parameters, temperature, density and concentration of the system or by means of the thermodynamic relations.
3. In the plane (ρ 0 , c 0 ) we distinguished two directions: the direction of strong fluctuations η 0 and the direction of weak fluctuations ξ 0 .The direction of strong fluctuations (the order parameter) is determined from the equation where θ is the rotation angle of axes η 0 and ξ 0 in the (ρ 0 , c 0 ) plane, δ = ρ(v a − v b ), v i = ( ∂V ∂N i ) T,P,N j =i is a partial volume.Based on the Gaussian distribution we have determined the critical branch and, correspondingly, CV η 0 connected with the order parameter.Now we can accomplish the second important task, namely, we can derive the basic density measure with respect to η k .To this end we shall follow the program: (i) having passed from CV ρ k and c k to CVη k and ξ k in (45), we shall integrate over irrelevant variables ξ k (which do not include the variable connected with the order parameter) with the Gaussian density measure; (ii) then we shall construct the effective GLW Hamiltonian (we shall restrict our consideration to the η 4 -model).
The result of the first step of the proposed program is the functional with respect to the variables η k (and χ k conjugate to η k ): where Ξ G ξ and ∆ ξ are the results of integrating over ξ k which do not include χ k , J(η) has the form: Here M n (0, . ..) = M(n) n (0, . ..) + △M n (0, . ..). △M n (0, . ..) are corrections obtained as a result of integrating over ξ k .
It can be shown that the behaviour of function P (k) in the neighbourhood of the point of phase transition is similar to the behaviour of the initial potential Φγδ (k): P (k) takes both negative (at small | k|) and positive (at large | k|) values.In the region | k| > B we can integrate over χ k and η k with the Gaussian measure density as the basic one.As a result, we obtain the similar expression for Ξ in which the summation is performed over | k| B and coefficients M 2 (0), M 2 (0), and M 4 (0) are replaced by new coefficients M2 (0), M2 (0), and M4 (0) [35].Similar to the onecomponent system we consider the quantity B as the size of the first Brillouin zone of the block lattice.
After the integration over χ k we can present (53) in the form (within the framework of the η 4 -model): where the effective GLW Hamiltonian E 4 (η) has the form: Here the following notations are introduced: Mn / N .U(a, z) is the parabolic cylinder function.The radius B in (54) is found from the condition P (k = B) = 0 and it is considered as the size of the first Brillouin zone of a certain block lattice.
E 4 (η) has the form analogous to the basic density measure of the 3D Ising model in an external field.But the main difference is the dependence of coefficients a 0 , a 2 and a 4 on the microscopic parameters of the system.

Gas-liquid critical temperature of the symmetrical square-well mixture
We consider a binary symmetrical fluid mixture, i.e., a system of equal-size particles interacting via the same attractive potentials between "like" particles (φ aa (r) = φ bb (r) = φ(r)) and via different attractive potentials between "unlike" particles (φ(r) = φ ab (r)).The concentration x = 0.5 is a critical one for this model mixture.Notwithstanding its simplicity, the symmetrical mixture exhibits all the three types of two-phase equilibrium which are observed in real binary fluids: gasliquid, liquid-liquid and gas-gas.
We consider the binary symmetrical mixture in the vicinity of the gas-liquid critical point.In this case the basic density measure has the same form as (54) but the main quantities from (54) are reduced to the following ones [26]: and the cumulants are reduced to the cumulants of a one-component system.The interaction in the system was described by the potential r < σ, −ǫ γδ , σ r < λσ, 0, r λσ.
where σ is a hard sphere diameter, λ is a range of the potential, and ǫ γδ is a welldepth of the interaction between the particles of types γ and δ.
For a symmetrical mixture the following relations occur: ǫ aa = ǫ bb = ǫ = ǫ ab .We introduce parameter r measuring the relative strength of the interaction between the "like" and "unlike" particles: r = ǫ ab /ǫ.The case r = 1 corresponds to a onecomponent fluid.Gaussian approx.--2.0 1.0 2.753 2.684(51) [40] 0.129 0.123(43) [40] Using the expression for the GLW Hamiltonian and the method of layer-by-layer integration proposed for the three dimensional Ising model [27,31] we calculated the critical parameters (critical temperature and critical density) of the symmetrical mixture square-well mixture [26].The results are presented in figure 2

Figure 1 .
Figure 1.The values for the phase transition temperature for the exact (curve 1) and approximate (curve 2) values of the partition function coefficients.
• •) is the Dirac delta function.The prime means that the product over k is performed in the upper semi-space.ρ k,γ = ρ c k,γ − iρ s k,γ are collective variables of the γth species, where the indices c and s denote the real part and the coefficient of the imaginary part of ρ k,γ .ρ k,γ describes the value of the k-th fluctuation mode of the number of γ-th species particles.Each ρ c k,γ and ρ s k,γ takes all the real values from −∞ to +∞. (dρ) is a volume element of the CV phase space:

Table 1 .
Coefficients of the partition function functional (a 2g , a 4g , a 6g ) and their approximate values (a 2s , a 4s , a 6s )

Table 2 .
Critical properties of a symmetrical binary mixture, calculated by the CV method and MC simulations (r = 1.0 corresponds to a one-component system).