Bernarda Lovrinčević

Adrien Bella

Isham Le Tenoux-Rachidi

Martina Požar

Franjo Sokolić

Aurélien Perera

Methanol-ethanol "ideal" mixtures as a test ground for the computation of Kirkwood-Buff integrals

Mixtures of 1-alkanols are a textbook example of the concept of ideal mixtures. Yet, such mixtures have a very strong local order due to the hydrogen bonding interactions, with a strong tendency for chain formation. Despite this apparent non-ideality, the Kirkwood-Buff integrals of such system exhibit near ideal behaviour. This dual property can be used to test the calculations of the Kirkwood-Buff integrals in a controlled mixing situation, and clarify many points, in particular the statistical problems that can be encountered. By studying the methanol-ethanol mixtures, we uncover an interesting physical asymmetry between low methanol and low ethanol concentrations, which can produce statistical artifacts in the calculation of Kirkwood-Buff integrals, illustrating and exemplifying some of the difficulties encountered in such calculations. Finally, liquid state integral equations results for these mixtures are reported. They help demonstrate that thermodynamic ideality hides complex correlations and microscopic non-ideality.

Introduction

The relationship between microscopic interactions and statistical correlation functions is essential to understand the physical and structural properties of liquids from a theoretical point of view [1,[START_REF] Chaikin | Principles of Condensed Matter Physics[END_REF]. The Kirkwoof-Buff integrals (KBI) [START_REF] Kirkwood | [END_REF] are a perfect tool to analyse this relationship, since they correspond to the integral of the pair correlation functions, and they can be related to various thermodynamic properties, such as the isothermal compressibility, the molar and partial molar volumes and the chemical potentials [START_REF] Kirkwood | [END_REF]4,5]. In particular, they allow to measure the deviation from ideality of mixtures in terms of the excess chemical potentials [START_REF] Kirkwood | [END_REF]5]. However, the ideality of the KBI for this type of mixtures poses several problems. A binary mixture can be considered ideal when both components interact quite similarly. Mixtures of alcohols, particularly mono-ols, are considered as ideal [START_REF] Koga | Solution Thermodynamics and its Application to Aqueous Solutions[END_REF][START_REF] Rowlinson | Liquids and Liquid Mixtures[END_REF] since the interaction between the hydroxyl groups dominates entirely all the other atomic interactions, and therefore those between the methyl groups can be considered as a small perturbation. Indeed, in a previous study [START_REF] Pozar | [END_REF] we have found that the KBI of the methanol-ethanol mixtures are nearly ideal, although there was considerable local structural ordering, due to the hydrogen bonding interactions.

We revisit here these mixtures, with the particular aim of sorting out how the differences in alkyl tails affects the pseudo-ideality of these mixtures, as well as the statistics. We report a surprising difference in statistics between low methanol and low ethanol concentrations, to which we provide a physical explanation. For the same level of statistical sampling, the low methanol content mixture seems to show methanol segregation, contradicting the expected ideality. This behaviour is absent from the mixtures with low ethanol content. Instantaneous snapshots seem to confirm this difference in homogeneity. However, when longer statistics are performed from the low methanol content side, the expected ideality is restored. This finding indicates that ideality is a global concept and could be violated at short times. Since this is not observed for low ethanol content, it introduces a kinetic parameter in the statistics. In addition, this difference is supported by a marked difference between the oxygen-oxygen and oxygen-hydrogen pair correlations, which are the markers of the role played by the hydroxyl groups in the hydrogen bonding clustering and chaining phenomena. We will discuss this situation in the Results section below.

Another way to test the ideality of such mixtures is to solve approximate integral equation theories (IET), such as the reference interaction site method (RISM) theory [START_REF] Hansen | Theory of Simple Liquids[END_REF]. Since they are approximations which neglect high order correlations, they are very sensitive to local order and global organisation.

Hence, they are indirect indicators of the complexity of a system. Indeed, the resulting equations cannot be solved for disordered systems which exhibit strong local order, such as aqueous mixtures, for example. In order to solve these IET, it is necessary to make approximations so drastic, that the resulting correlations are nearly featureless, when compared to those from computer simulations [START_REF] Perera | [END_REF]11]. In addition, these theories do not have the size dependence problems which computer simulations have. Anticipating the results of Section 4, we can state that methanol-ethanol mixture exhibit very strong correlations and local order, which these theories cannot capture. This poses the question of what thermodynamic ideality really means.

The present paper differs from the previous one [START_REF] Pozar | [END_REF] in that it presents a detailed account of the numerical determination of the KBI, a subject that has been of interest in recent investigations by many authors, as well as the above mentioned IET study.

The remainder of this paper is structured the following way. In the next section, we summarise the methods used in this work, and in particular the technical problems concerning the evaluation of the KBI. In the Results section, we examine in detail the calculation of the KBI by computer simulation and IET methods. Our conclusions are gathered in the last section.

Methodology

Computer simulations

The simulation protocol used herein is the same used in our previous study of the methanol-ethanol mixtures in Ref. [START_REF] Pozar | [END_REF]. We focus on the TraPPe models of both mono-ols [12], since the previous study has indicated that the OPLS model [13] results are not so much different, as far as ideality properties of the mixtures are concerned. We use the GROMACS 4.5.5 package [14]. N=2048 total number of molecules is considered, with a simulation time step of 2fs. Ambient conditions of T=300K and 1bar (0.98 atm) are imposed through the Nose-Hoover [15,16] and Parrinello-Rahman [17,18] algorithms, with relaxation times of 0.1ps and 1ps, respectively. The simulation box, with an initial volume estimated from the interpolation of experimental molar volumes of respective neat alcohols, is filled with both alcohols with the use of the PACKMOL program [19]. The system is then energy-minimised and 200ps NVT run is performed for an initial equilibration, followed by a 1ns NPT step, in order to ensure ambient conditions. Typically 2ns runs are performed for statistics, with 1000 independent configurations saved for calculation of the structural properties. In the cases of low methanol content 10ns runs were required, as discussed in the next section. These statistics are good enough to provide smooth looking correlations g iaj b (r) between atom i a belonging to molecular species a and atom j b of molecular species b. The atom-atom structure factors defined as :

S iaj b (k) = δ ab + ρ √ x a x b ˆd r [g iaj b (r) -1] exp i r. k (1)
are calculated by numerical fast Fourier transform. This calculation implies that the g iaj b (r) have been properly treated to have the correct asymptote of 1. Neglecting to make this correction will produce serious artifacts at very small k values. The Xray and neutron scattered intensities have been equally calculated using the expression derived previously [20] from the Debye expression [21,[START_REF] Debye | The collected papers of Peter J.W. Debye[END_REF]:

I(k) = ρ ab √ x a x b iaj b f ia f j b S (T ) iaj b (k) (2)
where the summations runs over all the species a and b, and atoms i a and j b within each species, ρ = N/V is the total density, f m is the form factor of atom m and depends on the nature of the scattered radiation (taken from crystallographic tables [START_REF] Prince | International Tables for Crystallography[END_REF]) and S (T ) iaj b (k) is the total structure factor, which accounts for the intra-molecular correlations and selects atoms which belong to the same molecular species, and is defined similarly to Eq.( 1)

S (T ) iaj b (k) = W iaj b (k)δ ab + ρ √ x a x b ˆd r [g iaj b (r) -1] exp i r. k (3)
where W mn (k) = j 0 (kd mn ) is the intra-molecular correlation function, which is related to the spherical Bessel function j 0 (x) = sin(x)/x in case of rigid molecule (an approximation we shall use here), and where d mn = | r m -r n | is the distance between atomic sites m and n belonging to the same molecule.

It is noteworthy that the scattered intensities reported in our earlier work [START_REF] Pozar | [END_REF] were based on the Pings-Waser expression [START_REF] Pings | [END_REF], which is misleading since it is missing the intra-molecular part W ab (k), which is crucial to match the experimental data. This deficiency is often "cured" by plotting kI(k), which allows to minimise the small-k problems produced by this incomplete expression. Hence, the present paper offers a more appropriate account of scattering properties than the previous one, showing the importance of including intra-molecular correlations.

In this work, we express the Xray scattering intensities in cm -1 . In this case, the formula for the scattering in Eq(2) needs to be multiplied by r 2 0 , where r e ≈ 2.8210 -13 cm is the electron radius. The appropriate factor to convert I(k) expressed as in Eq.(2) to cm -1 is the multiplicative factor κ = r 2 e N/L 3 /100, where N is the number of molecules in the simulations (here N = 2048) and L is the box size for a given mixture. The neutron scattering intensities are plotted in units the squared form factors, as given by the cristallographic tables [START_REF] Prince | International Tables for Crystallography[END_REF].

Integral equations

We briefly recall the RISM IET methodology [START_REF] Hansen | Theory of Simple Liquids[END_REF]. It consists in solving the site-site Ornstein-Zernike (SSOZ) equation , together with the hypernetted chain (HNC) equation, in order to find the atom-atom correlation functions g ab (r), as well as the corresponding direct correlation functions c ab (r). The SSOZ equation can be written in the Fourier space in a matrix form. For this purpose, we introduce the matrix element notation 3) and the matrix M = W -1 -C, the SSOZ equation can be written in a compact form:

F ij (k) = ρ √ x i x j fij (k), where f ij (r) is either h ij (r) = g ij (r) -1 or c ij ,
ij (k) = W ij (k) + H ij (k) defined in Eq.(
S (T) M = I ( 4
)
where I is the identity matrix.

The closure relation can be written as

g ij (r) = exp [-βv ij (r) + h ij (r) -c ij (r) + b ij (r)] (5)
where b ij (r) is the site-site bridge function which contain higher rank direct correlation functions [25]. Neglecting this function leads to the HNC closure

g ij (r) = exp [-βv ij (r) + h ij (r) -c ij (r)] (6)
In a previous work [26], we have extracted the bridge function for several neat liquids, and in particular methanol and ethanol. In the present work, we test the application of equations ( 6) and ( 5). This second equation is solved by weighting the neat alcohol bridge functions with the concentration as follows:

b i M j M (r) = a M M (1 -x)b (neat) i M j M (r) (7)
b

i E j E (r) = a EE (1 -x)b (neat) i E j E (r) (8)
for methanol and ethanol, respectively, and the cross bridge term is approximated as a mix in terms of the neat liquid bridge functions as:

b i M j E (r) = a M E (1 -x)b (neat) i M j M (r) + a EM xb (neat) i E j E (r) (9)
where the index M and E stand for atoms of the methanol and ethanol molecules, respectively, x is the ethanol mole fraction, and the 2 neat liquid bridge functions in Eq.( 9) are chosen such that the same pairs of atoms match. In the case of the methanol methyl and ethanol methylene united atoms, we have used the neat methanol methyl-methyl bridge function instead. The 4 coefficients a M M , a EE , a M E and a EM should ideally be 1. However, in practice, in many cases, it is not possible to bring some or all of these coefficients to 1 without loss of convergence. The reason for this problem is found in the anomalous raise of the S ij (k = 0), indicating that the mixture is subject to large concentration fluctuations from the HNC theory point of view. This anomaly is typical of the HNC closure, and the KH closure form is supposed to cure it. This closure is similar to the HNC closure, except that, in the exponential Eq.( 6), whenever the argument is positive, it is replaced by the MSA expression for g ij (r). This way, the tendency of HNC to overestimate positive correlations, is tamed down by a lower level closure. This procedure is known to provide numerical solutions whenever HNC fails to do so. However, the resulting correlations are weak and lack the typical structure of strongly correlated liquids. There is no systematic documentation of these problems, since it is not understood why these IET generally fail or are insufficient for complex liquids [27].

On the calculation of the KBI

In recent years, there has been a increased interest in the computation of Kirkwood-Buff integrals (KBI) from computer simulations. These quantities are defined as:

G ab = ˆd r [g iaj b (r) -1] (10)
where g iaj b (r) is the pair distribution function between atoms i a and j b , belonging to molecular species a and b, respectively, and separated by the distance r. The KBI G ab do not depend on specific atomic sites since the integral must be invariant to the choice of any particular centre in the molecules. The atom-atom distribution functions g iaj b (r) can be computed from the standard histogram method [START_REF] Allen | Computer Simulation of Liquids[END_REF], and the functions are usually sampled in the interval r ∈ [0, L/2], where L is the system cubic box size. Therefore, in view of Eq.( 1), it does not seem particularly difficult to numerically compute the desired integrals, which come down to evaluating:

G ab = 4π ˆL/2 0 dr r 2 [g iaj b (r) -1] (11)
Yet, this simple task is made difficult by the fact, discovered decades ago by Lebowitz and Percus [START_REF] Lebowitz | [END_REF], that g ij (r) in a finite size system would not converge to the expected limit, which is 1, and that the distortion would induce an asymptotic shift, which in the canonical ensemble would be:

g ij (r → +∞) → 1 - ij N ( 12
)
where ij is related to the macroscopic concentration fluctuation < N i N j > -< N i >< N j > through the relation, expressed in the Grand canonical ensemble as [START_REF] Kirkwood | [END_REF]30]

ij = V < N i N j > -< N i >< N j > < N i >< N j > = 1 ρx i x j ∂ρ i ∂βµ j T V µ k ( 13
)
where V is the volume of the system, ρ =< N > /V is the total number density (N = i N i where N i is the number of particles of species i), x i =< N i > / < N > is the mole fraction of species i, ρ i = x i ρ and β = 1/k B T is the Boltzmann factor (T is the temperature and k B the Boltzmann constant).

The relation for the isobaric ensemble is more involved [START_REF] Allen | Computer Simulation of Liquids[END_REF], and although the asymptote shift remains similar to Eq.( 12), the ij term is more complicated.

If the spurious shift in the asymptote is not properly corrected, then the so-called running KBI (RKBI) defined as

G ab (r) = 4π ˆr 0 ds s 2 [g iaj b (s) -1] (14)
instead of tending to the limit defined by Eqs. [START_REF] Perera | [END_REF]11) would have a r-dependent cubic curvature ∼ ij r 3 set by the spurious limit of Eq.( 12).

It is important to underline that it is the problem mentioned above which is responsible for the difficulties met with the obtention of the KBI through simulations, and not the finite size problems which are inherent to these methods. These problems can be suppressed by considering large enough systems. This problematic was previously analysed by Sarracuze and coworkers [31,32], who clearly state that the fluctuation problems play the most important role. Recent works on KBI [33,34,35,36] go back to putting emphasis on finite size problems, and considering them in par with the fluctuation problems. We believe it is more important to correct for the asymptote problem of Eq.12.

In many papers [37,38,[START_REF] Perera | chapter Micro-Heterogeneity in Complex Liquids[END_REF], we have proposed to correct for the shifted asymptotes by simply shifting the atom-atom pair correlation functions g ab (r).

There are some rules that should be obeyed, however. First, it seems reasonable to assume that the short range part of the g ab (r) are ensemble independent and should not be affected by fluctuations. This could cover the first and the second neighbour distances 0 < r < 2σ ab . Then, there should be a smooth continuation to the distance from which the asymptote may be defined. In earlier publications, we propose to use a multiplicative hyperbolic tangent function to provide the correct asymptote, with the shifted function given by g (shifted)

ab (r) = g ab (r) × 1 2α ab [1 + tanh( r -λ ab κ ab )] (15)
with α ab = lim r→L g ab (r) is the unshifted asymptote, λ ab ≈ σ ab and κ ab chosen such that the rate of change of the hyperbolic tangent covers 1 range of σ ab . These 3 parameters are very empirical, but testing the method shows that these parameters should be chosen in a very narrow range of values if the final result should have a physical meaning. The shift is monitored by trial and error choice of the shift parameter α ab by plotting the running KBI (RKBI) defined as

G ab (r) = 4π ˆr 0 ds s 2 g (shifted) ab (s) -1 (16)
a quantity which should reach an asymptote defined by Eq.( 14) which is the KBI for the pair function g ab (r).

In practice, it turns out that the shifting operation through the simpler recipe

g (shifted) ab (r) = 1 α ab g ab (r) (17)
is equally efficient, since the shift is so close to 1 that small error due to shifting the g ab (r) in the first neighbour region does not lead to noticeable changes, and is always much smaller than the error bars of estimating the proper KBI. In practice, and for all the cases we have tested so far, one has α ab = 1 + ab with | ab | 0.01. In the next section we will show a detailed analysis of the application of this method to the methanol ethanol mixtures, as well the various problems encountered.

In the cases when one of the species is in low content, the proper determination of the long range correlation may be difficult to obtain because long statistics may be required to ensure the sampling of the macroscopic homogeneity. There seems to be no precise rule for this. In such cases, the asymptote will be affected by the noisy sampling of the real homogeneity, and it may be possible to determine it properly. In such cases, one can use the following trick, which we have used in many of our previous papers. It consists in noting that, when one species is in small concentration, the other is in majority, and consequently its KBI is well defined through the asymptote shifting method. On such case, one can profitably use the analytical expression of the KBI, which can be approximated as follows

G 11 = G 12 + 1 1 -x V 2 D -V (18)
G 12 = V 1 V 2 V D (19)
G 22 = G 12 + 1 x V 2 D -V (20)
where x = x 2 is the mole fraction of species 2, V 1 and V 2 are the molar volumes of the neat species, V = (1-x)V 1 +xV 2 is the volume of the mixture, and D = ρ a ∂βµ a /∂ρ a (a = 1, 2) is the term which describes the deviation from ideality of the mixture. When D = 1, the chemical potentials µ a = k B T ln x a ρ (a = 1, 2) are ideal and the mixture is ideal. In the expressions above, we have neglected several terms, such as the compressibility of the mixture, which is usually very small since dense liquids are generally incompressible, and we replaced the partial molar volumes by the neat liquid volumes. In practice, we found that the simpler equations above are quite accurate against experimental data [START_REF] Perera | [END_REF]30]. These expressions allow to extract the term D from the well behaved KBI (the one corresponding to majority species) and replace it in the analytical expression for the problematic KBI. We found empirically that this method gives excellent results in place of the problematic asymptote extraction.

Results

Many structural properties for this particular mixture have been previously reported in Ref. [START_REF] Pozar | [END_REF]. Therefore, in this section, we will focus on the KBI problem and the integral equation approach to these mixtures, in relation to their ideality.

KBI from simulations

In this section we present the information related to the KBI of the methanolethanol mixtures, by using the asymptote shifting method Eq. (17).

First of all, it is important to verify that the asymptotes of the various atom-atom correlations are affected by fluctuations the same way for each species pairs. This is an internal check of the consistency of the calculations. We illustrate this consistency in Fig. 1(a-b) for the 60% ethanol mixture.

Fig. 1a shows all of the atom-atom structure factors (Eq.( 1)), with different colors for each species pairs, blue for all the methanol pairs, magenta for the ethanol pairs, and green for the cross correlations. The atomic details do not matter for this demonstration, which is why these are not specified. We see that the various atom-atom structure factor S ab (k) are different until nearly k ≈ 0, where they converge to 3 distinct values related to the 3 KBI. It is noteworthy that, to the numerical estimate, the convergence point seems to be k = 0, and there is no superposition starting at some non-zero k value. This is the case when micro-segregation occurs [30,41], which is not the case here: the methanol-ethanol mixtures are therefore microscopically homogeneous.

Fig. 1b shows the same information, but from the RKBI point of view, as defined in Eq.( 16). All the RKBI corresponding to all the appropriately shifted atom-atom pair correlation functions, are shown with the same color convention as Fig. 2a. It is seen that there are 3 distinct asymptotes, corresponding to the 3 KBI, namely G mm (blue), G me (green) and G ee (magenta). Although the individual atom-atom pair correlations differ at short and medium distance, they all converge to the same fluctuation influenced KBI value. Both the shifted and unshifted values are shown. These correspond to 10ns statistics.

Fig. 1c illustrates the same RKBI, but for 2ns statistics. Since this mixture corresponds to rather low methanol content (40%), we see how the spurious clustering kinetics influences the statistics and the asymptotes, which cannot be made horizontal, even by shifting the shifting procedure. There is a residual curvature, which witnesses the inhomogeneity in the outer reservoir provided by large r-values, as explained in Section 3.1. The statistic problem is so serious that we even observe an nonphysical crossover of the methanol and cross correlation RKBI curves.

Fig. 2 shows the RKBI G OmOm for the methanol oxygen atom O m from the correlation functions g OmOm (r), for various ethanol mole fractions x ranging from x = 0(pure methanol) to x = 0.95. These functions are shown to illustrate the asymptote shifting method (the G OmOm corresponding to the non-shifted g OmOm (r) are shown in gray lines, and are always tilting down to various degrees, depending on the magnitude of the shift). It is seen that the shifted RKBI reach a reasonably horizontal asymptotes for the small x values corresponding to the rich methanol side, where we expect less concentration fluctuations of this species. However, above the equimolar mixture x > 0.5, it is seen that the asymptotes deteriorate and start to show spurious noisy behaviour, which correspond to large heterogeneity in the distribution of the methanol molecules and small methanol clusters. 10ns runs are not sufficient to clarify the sampling problems. In such cases, the determination of the asymptotes is impossible. However, we can use the trick explained at the end of Sub-section 3.1. The ideal KBI are equally shown in Fig. 3 as green horizontal lines for all concentrations ranging from x = 0 to x = 1. It is seen that for x<0.6 the asymptotes nearly superpose to these ideal KBI, indicating that the mixture is indeed ideal. This cannot be asserted for x>0.7, but we will see below that these mixtures are equally ideal.

Fig. 3 shows the RKBI G OmOe for the cross methanol-ethanol oxygen atom O m and O e obtained from the correlations g OmOe (r). In this case, we see that all the RKBI converge to their ideal values, although some noise is still present in some asymptotes, particularly at low methanol content when x > 0.7. These results again assert that the mixtures are nearly ideal.

Fig. 4 shows RKBI G OeOe for the ethanol oxygen atom O e , obtained from the correlations g OeOe (r), for various ethanol mole fractions x ranging from x = 0.05 to x = 1(pure ethanol). This time around, it is the asymptotes for small ethanol content which are noisy and not well defined, since this is now the minority species, and dispenses problems to sample the macroscopic homogeneity within a small system and microscopic sampling times. However, we note that the problem appear less serious than for the case of low methanol content.

This asymmetry between the low ethanol and low methanol content, witness the near ideal character of the mixture, since each molecules differ by one methyl atom. Nevertheless, Fig. 5a shows, that when trick for correcting the deficient asymptotes is used, then the KBI are near perfectly ideal, as they superimpose to the ideal curves obtained from Eqs. (18)(19)(20) by setting D = 1. If one wants to get the proper KBI from a computer simulation, for these particular cases where it seems difficult to obtain them without the trick, then one need to simulate very large systems, and for even longer times, in order to sample the proper homogeneity of these mixtures, which tend to have long lived heterogeneity. This is a real operational problem, which can be circumvented in the present case, but which are very much problematic for aqueous mixtures [42,43,44,45], for example. The present study help understand the nature of the practical problem which appear under other context.

Fig. 5b illlustrates the LP shifts obtained from the simulations as ij x i versus ethanol concentration, where ij is given by Eqs. (12,13). In practice, we numerically shift by the quantity θ ij = ij /N which is very small because of N (here N = 2048). We can always express the chemical potential of species i in a binary mixture as βµ i = ln(ρ i ) + ρ i B ii + ρ j B ij which represent the ideal term and the contribution from excess term [1,[START_REF] Koga | Solution Thermodynamics and its Application to Aqueous Solutions[END_REF][START_REF] Rowlinson | Liquids and Liquid Mixtures[END_REF], and where B ii and B ij contains all the high order density terms for a binary mixture. Then, using Eq.( 13), we see that, for like species (i = j) we have ii ≈ (1 -ρ i B ii )/x i to first order in density. Fig. 5b shows that, for methanol (in blue) and in ethanol (in magenta), ii x i ≈ 1 -x i B ii despite the noisy data. If the mixture were ideal, we would expect this relation to hold exactly, and where B ii would be the second virial coefficient. As for the cross term, the cross derivatives in Eq.( 13) can be evaluated through the Gibbs-Duhem equation and, for ideal mixtures, one indeed expects that x i ij = -1. Therefore, Fig. 5b confirms directly the ideality of the mixture, despite the noisy data coming from the empirical procedure of shifting the asymptote visually to 1. It serves to illustrate both the procedure, its shortcomings and its validity in this particular case. We wish to stress that Eq.( 13) does not hold for the NPT ensemble, and yet we could get very convincing results.

Integral equation results

Fig. 6, shows the 3 oxygen-oxygen correlation functions for the 20% ethanol mixture, as obtained from computer simulations (in blue), HNC (in green), HNC+bridge (in red). It is seen that the HNC approximation alone cannot reproduce the strong chain correlations which lead to the sharp first peak indicated by computer simulations, and the second peak is entirely missed. However, adding the approximated bridge (as explained in Section 2.2) allows to catch up to the first peak to some extent while fitting the second peak perfectly. This means that the neat liquid bridge alone are not enough to reproduce the strong correlations which lead to the first peak in the mixtures. As for the structure factor, the HNC+bridge approximation is relatively good, except for a notable rapid increase very close to k=0, which is not observed for the HNC approximation. This means that adding the new bridge, in addition to bringing a better agreement, also increases the concentration fluctuations. We find this coupling between the good description of local order with the spurious increase of concentration fluctuations, quite noteworthy. It certainly reveals how local heterogeneity and disorder are coupled. Fig. 7 shows the same information as Fig. 6, but for the equimolar mixture. In addition, we show the results of the KH+bridge closure (dashed purple lines), which is seen to be overall weaker than the HNC+bridge closure, but is seems better for the small k part of the structure factors. However, it suffers also of the spurious k = 0 raise, hence showing that the trick underlying this approximation works only in particular cases. Since HNC+bridge is clearly superior, this analysis shows that the KH closure should be used, since it has no theoretical justifications. The values obtained for the coefficients in Eqs. [START_REF] Rowlinson | Liquids and Liquid Mixtures[END_REF][START_REF] Hansen | Theory of Simple Liquids[END_REF][START_REF] Pozar | [END_REF] are a M M = 1, a M E = a EM = 0.90 and a EE = 0.85 for the HNC+bridge solution. It was not possible to raise them further to 1, and trying to do so led to dramatic increase of the k=0 part of the structure factors, suggesting anomalous concentration fluctuation. However, this could be achieved for the KH closure, with the modest k=0 increase seen in Fig. 7.

Fig. 8 shows IET results for the 80% ethanol content. Very similar observation to Fig. 6 and Fig. 7 can be made. It is interesting that the ethanol correlations for all 3 concentrations are much better reproduced by the bridge correction than the methanol oxygen correlations. This is much against intuition, since one would expect the less complex molecule to be better described than the more complex one. The bare HNC results are very similar between the 3 cases. These 2 results combined, suggest that the methanol bridge is more complex than that provided by our approximation in Eq.( 7) and to some extend Eq.( 9).

The difficulties of solving IET for these mixtures without the need of the bridge function, clearly suggest the crucial role of the irreducible higher order rank correlations. For simply disordered liquids -typically Lennard-Jones mixtures, neglecting the bridge function does not alter dramatically the structural properties. However, for complex disorder liquids -typically hydrogen bonded liquids and mixtures, it is often impossible to solve IET. In the present case, it is still possible to solve HNC, but with very weak structural features. This is somewhat in contradiction with the fact such mixtures should be thermodynamically ideal, which usually suggests that the interactions are not strong. So, these mixtures are both ideal from thermodynamics, and highly non-ideal as far as structural features are concerned. The thermodynamic ideality arises from a possible cancellation of very strong but opposing non-ideal contributions. The present IET analysis brings an interesting complementary view on this problem, which we hinted to in a previous study [START_REF] Pozar | [END_REF].

Xray and Neutron scattering

Fig. 9a-b shows the Xray and neutron intensities, as computed from the structure factors derived from the atom-atom correlation functions obtained from computer simulations. First of all, we observe that the scattering pre-peak, witnessing chain clustering of the hydroxyl groups, is clearly observed for the Xray intensities, and seen to continuously evolve from a shoulder at k ≈ 1 Å -1 for pure methanol, to a well defined peak at k ≈ 0.7 Å -1 for pure ethanol. It is seen in both cases that the shape of the main peak, in the range k ≈ 1.5 -1.8 Å -1 , has an interesting cross-over when going from pure methanol to pure ethanol. Upon increasing ethanol mole fraction, this main peak first decreases, then shifts to smaller k-values, and then increases monotonously. This cross-over corresponds to a change in the prepeak, which goes from a shoulder to a pronounced pre-peak. So, it should correspond to a structural change when going from methanol type to the ethanol type clustering. This should be detected in spectroscopy studies, but the currently available information [46] does not provide any conclusive evidence. Interestingly, the Xray diffraction experiments shown in Fig.9a [47,48], clearly show the crossover since the main peaks agree with the calculated data. We observe that the 2 experimental data for ethanol differ considerably, indicating that unambiguous experimental Xray diffraction is not achieved. The discrepancy with the calculated data can be justified by modeling inaccuracies. Nevertheless, the shapes of the pre-peak behaviour for methanol and ethanol are quite well reproduced, which is not obvious since these correspond to emergent (cluster) structures, and should be more sensitive to model dependence. The TraPPE model [12] seems to capture these features quite well.

Conclusion

In this work, we have revisited the methanol-ethanol mixtures, in order to demonstrate how strong directional hydrogen bonding interactions affect the evaluation of the correlation functions, through the statistical problem posed by the existence of long lived hydrogen bonded structure. These structures create an artificial local inhomogeneity which will not melt on usual simulation time scales, and requiring longer runs and statistics. This problem is illustrated by the low concentration asymmetry observed for small statistics, more biased towards the small methanol content side. The existence of such structures, mostly chain-like oxygen clusters, homogeneously mixing both types of oxygen atoms, has been demonstrated by a cluster study in our previous work [START_REF] Pozar | [END_REF]. These structures create a strongly non-ideal system, from the point of view of correlations, but their integrals, namely the KBI, behave as if it was an ideal mixture. The apparently contradictory behaviour indicates that the thermodynamic signature is incomplete, and can be biased by canceling contributions. We have also illustrated the methodological problems in evaluating the KBI from computer simulations, together with the problems related to the misconception of the importance of the LP shift. We have shown how the evaluation of the KBI are affected both by the LP shift in Eq.( 12) as well by the statistical problems inherent to this mixture. This study, through the problems posed by the statistics and approximate IET approaches, indicates that substantial non-ideality is hidden behind the apparent thermodynamic ideality. Blue lines for data from simulations (from bottom to top x = 0.05, 0.1 to 0.9 by steps of 0.1, and 0.95); green lines for ideal KBI (from bottom to top x = 0, 0.5, 0.1 to 0.9 by steps of 0.1, 0.95 and 1), and gray lines for RKBI from uncorrected asymptotes (see text). The scaled LP shift ij x i (see Eqs. (12,13). Symbol for values from simulations, lines for ideal behaviour. In both panels, blue is for methanol, magenta for ethanol and green for cross contributions. for reference, the green curve is HNC Eq.( 6), the red curve is the HNC+bridge described by the Eqs. (5,[START_REF] Rowlinson | Liquids and Liquid Mixtures[END_REF][START_REF] Pozar | [END_REF][START_REF] Hansen | Theory of Simple Liquids[END_REF] The insets show the corresponding 3 structure factor S O X O Y (k) with same color conventions.

Figure captions

Fig.7

Same as Fig. 6, but for the equimolar mixture. The dashed purple curve is for the KH approximation with the same bridges as the HNC+bridge approximation. (12,13). Symbol for values from simulations, lines for ideal behaviour. In both panels, blue is for methanol, magenta for ethanol and green for cross contributions. .

Fig. 1 1 Fig. 2 RKBI

 112 Fig.1 Illustration of the self-consistent convergence of all the atom-atom RKBI of a mixture to the 3 common asymptotes which represent the 3 KBI values (here for ethanol mole fraction x = 0.6 of the methanol-ethanol mixture.). Blue for methanol, green for ethanol and magenta for methanol-ethanol cross correlations. (a) Convergence of the atom-atom structure factors at k=0 into the 3 independant KBI related value. (b) Proper RKBI asymptotes for 10ns statistics. (c) Incorrect RKBI asymptotes for 2ns statistics. See text in Section 3.1 Fig.2 RKBI for the methanol-methanol oxygen atom correlations, for all the ethanol concentrations (curves shifted by +10 for clarity).Blue lines for data from simulations (from bottom to top x = 0.05, 0.1 to 0.9 by steps of 0.1, and 0.95); green lines for ideal KBI (from bottom to top x = 0, 0.5, 0.1 to 0.9 by steps of 0.1, 0.95 and 1), and gray lines for RKBI from uncorrected asymptotes (see text).

Fig. 3

 3 Fig.3Same as Fig.2, but for methanol-ethanol oxygen atom correlations. Blue and Grey lines for x = 0.05, 0.1 to 0.9 by steps of 0.1, and 0.95. Green lines as in Fig.3.

Fig. 4

 4 Fig.4 Same as Fig.2, but for ethanol-ethanol oxygen atom correlations. Blue and Grey lines for x = 0.05,, 0.1 to 0.9 by steps of 0.1, 0.95 and 1. Green lines as in Fig.3. Fig.5 (a) KBI as obtained from computer simulations as full symbols, and lines for ideal values. Open symbols represent the anomalous values (from direct extraction from Fig.2-Fig.4 -see text). (b)The scaled LP shift ij x i (see Eqs.(12,13). Symbol for values from simulations, lines for ideal behaviour. In both panels, blue is for methanol, magenta for ethanol and green for cross contributions.

Fig. 6

 6 Fig.6 Integral equation results for the 20% ethanol mixture, for the 3 oxygen-oxygen correlation functionsg O M O M (r) (a), g O M O E (r) (b) and g O E O E (r) (c).Blue curves are computer simulation results for reference, the green curve is HNC Eq.(6), the red curve is the HNC+bridge described by the Eqs.(5,[START_REF] Rowlinson | Liquids and Liquid Mixtures[END_REF][START_REF] Pozar | [END_REF][START_REF] Hansen | Theory of Simple Liquids[END_REF] The insets show the corresponding 3 structure factor S O X O Y (k) with same color conventions.

Fig. 8

 8 Fig.8 Same as Fig.6, but for the 80% ethanol mixture.

Fig. 9 Fig. 5 -

 95 Fig.9Xray (a) and Neutron (b) scattered intensities from computer simulation structure factors. Neat methanol data (x = 0) is in dark brown and neat ethanol data (x = 1) in black, x = 0.05 in red, 0.1 in orange, 0.2 in gold, 0.3 in lime, 0.4 in green, 0.5 in jade, 0.6 in azur, 0.7 in blue, 0.8 in cobalt blue, 0.9 in magenta and 0.95 in purple. Symbols are experimental Xray diffraction data from Ref.[47] for methanol (orange triangles) and ethanol (blue triangles), from Ref.[48] for ethanol (red dots).

..Fig. 7 -Fig. 8 -

 78 Same as Fig.6, but for the equimolar mixture. The dashed purple curve is for the KH approximation with the same bridges as the HNC+bridge approximation. . . . Same as Fig.6, but for the 80% ethanol mixture. .

Acknowledgments A. B. and I. L. T.-R. thank the LPTMC for providing computational facility for their internship. This work has been supported in part by the Croatian Science Foundation under the project UIP-2017-05-1863 ''Dynamics in microsegregated systems''.

. . Fig. 2 -RKBI for the methanol-methanol oxygen atom correlations, for all the ethanol concentrations (curves shifted by +10 for clarity). Blue lines for data from simulations (from bottom to top x = 0, 0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.9 and 0.95); green lines for ideal KBI (from bottom to top x = 0, 0.05, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8 and Grey lines for RKBI from uncorrected asymptotes (see text). . [47] for methanol (orange triangles) and ethanol (blue triangles), from Ref. [48] for ethanol (red dots).


Angela Felisia

Chiarello

An overview of non-local traffic flow models

An overview of non-local traffic flow models

Felisia Angela

Introduction

There are three different important approaches to model traffic flow phenomena: microscopic models, mesoscopic models and macroscopic models. Microscopic models are based on the individual vehicles and their interactions; particular micromodels are cellular automata that prescribe stochastic rules how vehicles behave into cells in which the road is divided; mesoscopic models apply gas-kinetic equations to describe the traffic flow; macroscopic models consider quantities that characterize the collective behaviour of vehicles. In this paper we will consider the framework of macroscopic models. The first macroscopic traffic flow model, based on fluiddynamics equations, has been introduced in transportation literature since the midfifties of last century. It is the Lighthill, Whitham and Richards (LWR) model [START_REF] Lighthill | On kinematic waves. II. A theory of traffic flow on long crowded roads[END_REF][START_REF] Richards | Shock waves on the highway[END_REF]. The LWR model consists in one scalar equation that expresses the conservation of the number of cars:

∂ t ρ + ∂ x (ρ v(ρ)) = 0,
with ρ = ρ(t, x) representing the mean traffic density, i.e. the number of vehicles per unit length and v denoting the mean velocity. Due to the difficulty to simulate the Felisia Angela Chiarello Inria Sophia Antipolis Mediterranée, Université Côte d'Azur, 2004, route des Lucioles-BP 93, 06902 Sophia Antipolis Cedex, France, e-mail: felisia.chiarello@inria.fr 1 formation of traffic jams, several approaches have been developed during the years, addressing the need for more sophisticated models. In particular, the classical LWR model does not match the experimental data because it is based on the assumption that the mean traffic velocity is a function of the traffic density, which is not realistic in congested regimes. For these reasons, high order models were developed. In particular, in [START_REF] Payne | Models of Freeway Traffic and Control[END_REF][START_REF] Whitham | Linear and nonlinear waves[END_REF] and [START_REF] Aw | Resurrection of "second order" models of traffic flow[END_REF][START_REF] Zhang | A non-equilibrium traffic model devoid of gas-like behavior[END_REF], the authors considered second order models consisting in a mass conservation equation for the density and an acceleration balance law for the speed. Recently, "non-local" versions of the LWR model have been proposed in [START_REF] Blandin | Well-posedness of a conservation law with non-local flux arising in traffic flow modeling[END_REF][START_REF] Sopasakis | Stochastic modeling and simulation of traffic flow: asymmetric single exclusion process with Arrhenius look-ahead dynamics[END_REF]. In this type of models, the speed depends on a weighted mean of the downstream traffic density. As a consequence, the speed is a Lipschitz function with respect to space and time variables, ensuring bounded acceleration and overcoming the limitation of classical macroscopic models that allows for speed discontinuities. Non-local traffic models are intended to describe the behaviour of drivers that adapt their velocity with respect to what happens in front of them. For this reason, the flux function depends on a "downstream" convolution term between the density of vehicles and a kernel function supported on the negative axis ] -∞, 0] . As in classical (local) models, the speed is a monotone non-increasing function, because the higher is the density of cars on a road, the lower is their speed. There are general existence and uniqueness results for non-local equations in [START_REF] Amorim | On the numerical integration of scalar nonlocal conservation laws[END_REF][START_REF] Keimer | Existence, uniqueness and regularity results on nonlocal balance laws[END_REF] for the scalar case in one-space dimension, in [START_REF] Colombo | Control of the continuity equation with a non local flow[END_REF][START_REF] Keimer | Existence, uniqueness and regularity of multi-dimensional nonlocal balance laws with damping[END_REF] for multi-dimensional scalar equations, and in [START_REF] Aggarwal | Nonlocal systems of conservation laws in several space dimensions[END_REF] for the multi-dimensional system case. To prove the existence of solutions for these non-local models, there are mainly two different approaches. One is providing suitable compactness estimates on a sequence of approximate solutions constructed through finite volume schemes, as in [START_REF] Aggarwal | Nonlocal systems of conservation laws in several space dimensions[END_REF][START_REF] Amorim | On the numerical integration of scalar nonlocal conservation laws[END_REF]. Another approach relies on characteristics and fixed-point theorems, as proposed in [START_REF] Keimer | Existence, uniqueness and regularity results on nonlocal balance laws[END_REF][START_REF] Keimer | Existence, uniqueness and regularity of multi-dimensional nonlocal balance laws with damping[END_REF]. The first approach [START_REF] Amorim | On the numerical integration of scalar nonlocal conservation laws[END_REF] requires Kružkov-type entropy conditions to prove the L 1 -stability with respect to the initial data through the doubling of variable technique, while in [START_REF] Keimer | Existence, uniqueness and regularity results on nonlocal balance laws[END_REF][START_REF] Keimer | Existence, uniqueness and regularity of multi-dimensional nonlocal balance laws with damping[END_REF] the uniqueness of weak solutions is obtained directly from the fixed point theorem. An interesting point is the numerical approximation of the solutions of non-local models. It results challenging due to the high non-linearity of the system and the dependence of the flux function on integral terms, which highly increase the computational cost. There are different ways to numerically integrate non-local conservation laws. A first order Lax-Friedrichs-type numerical scheme is used to approximate the problem and to prove the existence of solutions in [START_REF] Aggarwal | Nonlocal systems of conservation laws in several space dimensions[END_REF][START_REF] Amorim | On the numerical integration of scalar nonlocal conservation laws[END_REF][START_REF] Betancourt | On nonlocal conservation laws modelling sedimentation[END_REF][START_REF] Blandin | Well-posedness of a conservation law with non-local flux arising in traffic flow modeling[END_REF][START_REF] Goatin | Well-posedness and finite volume approximations of the LWR traffic flow model with non-local velocity[END_REF] . Another first order numerical scheme is proposed in [START_REF] Friedrich | A Godunov type scheme for a class of LWR traffic flow models with non-local flux[END_REF], where an upwind numerical scheme for a class of non-local flux problems is presented. Alternative first order schemes based on the splitting of the non-local equation in two equations, i.e. the Lagrangian and the remap steps are proposed in [START_REF] Chiarello | Lagrangian-antidiffusive remap schemes for non-local multi-class traffic flow models[END_REF]. In comparison with Lax-Friedrichs scheme and upwind scheme, Lagrangian-Antidiffusive Remap schemes are much less diffusive. Concerning high-order numerical schemes, it is worth citing the papers [START_REF] Chalons | High-order numerical schemes for one-dimensional nonlocal conservation laws[END_REF][START_REF] Chiarello | High-order Finite Volume WENO schemes for non-local multi-class traffic flow models[END_REF][START_REF] Friedrich | Maximum principle satisfying CWENO schemes for nonlocal conservation laws[END_REF]. In [START_REF] Chalons | High-order numerical schemes for one-dimensional nonlocal conservation laws[END_REF], the authors propose discontinuous Galerkin and finite volume WENO schemes to obtain high-order approximations of non-local scalar conservation laws in one space dimension, where the velocity function depends on a weighted mean of the conserved quantity. The discontinuos Galerkin schemes give the best results, but under a very restrictive Courant-Friedrichs-Lewy (CFL) condition. On the contrary, finite volume WENO schemes can be implemented on larger time steps. In [START_REF] Chiarello | High-order Finite Volume WENO schemes for non-local multi-class traffic flow models[END_REF], the authors extend to 1D systems the finite volume WENO schemes proposed in [START_REF] Chalons | High-order numerical schemes for one-dimensional nonlocal conservation laws[END_REF]. In [START_REF] Friedrich | Maximum principle satisfying CWENO schemes for nonlocal conservation laws[END_REF], central WENO schemes are proposed, which, in contrast with the other high-order schemes for non-local conservation laws, neither require a restrictive CFL condition nor an additional reconstruction step. The paper is organized as follows: Section 2 is devoted to scalar non-local traffic flow models and the corresponding analytical results. A non-local multi-class traffic flow model is recalled in Section 3 and Section 4 presents an upwind numerical scheme, used in Section 5 for numerical tests of the multi-class non-local model with three classes of vehicles.

A class of scalar non-local traffic flow models

In this section, we recall the class of scalar conservation laws with non-local flux arising in traffic modeling studied in [10]

∂ t ρ + ∂ x f (ρ)v(ω η * ρ) = 0, x ∈ R, t > 0, (1)
where

ω η * ρ(t, x) := ∫ x+η x ω η (y -x)ρ(t, y)dy, η > 0. (2)
In ( 1), ( 2), we assume the following hypotheses:

(H) f ∈ C 1 (I; R + ), I = [a, b] ⊆ R + , v ∈ C 2 (I; R + ) s.t. v ≤ 0, ω η ∈ C 1 ([0, η]; R + ) s.t. ω η ≤ 0, ∫ η 0 ω η (x)dx := J 0 , ∀η > 0, lim η→∞ ω η (0) = 0.
In this model, the function ρ(t, x) represents the density of vehicles on the road and η > 0 is proportional to the look-ahead distance of drivers. Indeed, the aim is to describe the reaction of drivers that adapt their speed with respect to a weighted mean of the downstream traffic density. In particular this class of equations includes the traffic flow models presented in [START_REF] Blandin | Well-posedness of a conservation law with non-local flux arising in traffic flow modeling[END_REF][START_REF] Goatin | Well-posedness and finite volume approximations of the LWR traffic flow model with non-local velocity[END_REF][START_REF] Li | Shock formation in a traffic flow model with Arrhenius look-ahead dynamics[END_REF][START_REF] Sopasakis | Stochastic modeling and simulation of traffic flow: asymmetric single exclusion process with Arrhenius look-ahead dynamics[END_REF], which are detailed below. In [START_REF] Sopasakis | Stochastic modeling and simulation of traffic flow: asymmetric single exclusion process with Arrhenius look-ahead dynamics[END_REF], the authors introduce a traffic flow model based on Arrhenius stochastic microscopic dynamics. Using scaling and limit arguments they obtain a macroscopic description of the microscopic dynamics leading to higher-order dispersive partial differential equations. The dynamics includes interactions with other vehicles ahead. The non-local traffic flow model that derives from this stochastic process is

∂ t ρ + ∂ x (ρ(1 -ρ) exp(-ω η * ρ)) = 0, ρ(0, x) = ρ 0 (x), (3)
where the kernel acts only on the space variable x in this way

ω η * ρ(t, x) = ∫ x+η x J 0 η ρ(t, y)dy. (4)
The kernel ω η is an anisotropic short range inter-vehicle interaction potential, η is proportional to the look-ahead distance and J 0 is the interaction strength. Numerical solutions in [START_REF] Kurganov | Non-oscillatory central schemes for a traffic flow model with arrehenius look-ahead dynamics[END_REF] show that solutions of ( 3)-( 4) allow for shock formations in finite time. In [START_REF] Li | Shock formation in a traffic flow model with Arrhenius look-ahead dynamics[END_REF], the authors show that the finite time blowup of the solution of ( 3)-( 4) must occur at the level of the first order derivative and the L p , 1 ≤ p ≤ ∞, norms of the solution are finite close to the blowup time. In this way, the results confirm the formation of shock waves for this non-local model. While pedestrians are likely to react to the presence of people all around them, drivers mostly adapt their velocity to the downstream traffic, assigning a greater importance to closer vehicles. In [START_REF] Blandin | Well-posedness of a conservation law with non-local flux arising in traffic flow modeling[END_REF][START_REF] Goatin | Well-posedness and finite volume approximations of the LWR traffic flow model with non-local velocity[END_REF], the authors consider the following mass conservation equation for traffic flow with non-local velocity depending on a mean downstream density:

∂ t ρ(t, x) + ∂ x ρ(t, x)v ∫ x+η x ρ(t, y)ω η (y -x)dy = 0, for t ∈ R + and x ∈ R, η > 0. The kernel function ω η ∈ C 1 ([0, η]; R +
) is nonincreasing and the support η is proportional to the look-ahead visibility. The mean speed function v is continuous and decreasing. In [START_REF] Berthelin | Regularity results for the solutions of a non-local model of traffic flow[END_REF], the authors also study the regularity results for the solutions of this non-local model proving Sobolev estimates and the convergence of approximate solutions solving a viscous non-local equation.

The Blandin-Goatin model belongs to the class of equations considered in [START_REF] Keimer | Existence, uniqueness and regularity results on nonlocal balance laws[END_REF] and in [START_REF] Chiarello | Global entropy weak solutions for general non-local traffic flow models with anisotropic kernel[END_REF]. In [START_REF] Friedrich | A Godunov type scheme for a class of LWR traffic flow models with non-local flux[END_REF], the authors propose a general non-local vehicular traffic flow model based on a mean downstream traffic velocity. This is the main difference with the models in [START_REF] Blandin | Well-posedness of a conservation law with non-local flux arising in traffic flow modeling[END_REF][START_REF] Chiarello | Global entropy weak solutions for general non-local traffic flow models with anisotropic kernel[END_REF], where drivers adapt their velocity with respect to a mean downstream traffic density. In particular, the model in [START_REF] Friedrich | A Godunov type scheme for a class of LWR traffic flow models with non-local flux[END_REF] allows to capture changes in the velocity function and for this reason, it could be extended to junctions, as in [START_REF] Chiarello | A non-local traffic flow model for 1-to-1 junctions[END_REF].

Analytical results for a class of scalar non-local traffic flow models

Adding an initial condition

ρ(0, x) = ρ 0 (x), x ∈ R, (5)
with ρ 0 ∈ BV(R; I), entropy weak solutions of the Cauchy problem (1), ( 5), are intended in the following sense [START_REF] Amorim | On the numerical integration of scalar nonlocal conservation laws[END_REF][START_REF] Betancourt | On nonlocal conservation laws modelling sedimentation[END_REF][START_REF] Kružkov | First order quasilinear equations with several independent variables[END_REF].

Definition 1 A function ρ ∈ (L 1 ∩ L ∞ ∩ BV)(R + × R; I
) is an entropy weak solution of (1), ( 5), if

∫ +∞ 0 ∫ R | ρ -κ|ϕ t + sgn(ρ -κ)( f (ρ) -f (κ))v(ω η * ρ)ϕ x -sgn(ρ -κ) f (κ)v (ω η * ρ)∂ x (ω η * ρ)ϕ dxdt + ∫ R | ρ 0 (x) -κ|ϕ(0, x)dx ≥ 0 (6)
for all ϕ ∈ C 1 c (R 2 ; R + ) and κ ∈ R. The following results are obtained in [START_REF] Chiarello | Global entropy weak solutions for general non-local traffic flow models with anisotropic kernel[END_REF], approximating the problem (1)-( 5) through a Lax-Friedrichs type numerical scheme and recovering L ∞ and BV estimates for the sequence of approximate solutions. Stability with respect to the initial data is obtained from the entropy condition using the doubling of variable technique.

Theorem 1 (Well-posedness)

Let hypotheses (H) hold and ρ 0 ∈ BV(R; I). Then the Cauchy problem (1), ( 5), admits a unique weak entropy solution ρ η in the sense of Definition 1, such that

min R {ρ 0 } ≤ ρ η (t, x) ≤ max R {ρ 0 }, for a.e. x ∈ R, t > 0. ( 7
)
Moreover, for any T > 0 and τ > 0, the following estimates hold:

TV(ρ η (T, •)) ≤ e C(ω η )T TV(ρ 0 ), (8a)
ρ η (T, •) -ρ η (T -τ, •) L 1 ≤ τe C(ω η )T ( f v + J 0 f v ) TV(ρ 0 ), (8b)
with

C(ω η ) := ω η (0) v ( f ρ 0 + 2 f ) + 7 2 J 0 f v .
Under more regular assumptions, an estimate of the dependence of the solution with respect to the kernel function, the speed and the initial datum ca be proven, see [START_REF] Chiarello | Stability estimates for non-local scalar conservation laws[END_REF]. Above, and in the sequel, we use the compact notation • for • L ∞ .

Corollary 1 (Limit model as η → +∞)

Let hypotheses (H) hold and ρ 0 ∈ BV(R; I). As η → ∞, the solution ρ η of (1), ( 5) converges in the L 1 loc -norm to the unique entropy weak solution of the classical Cauchy problem

∂ t ρ + ∂ x ( f (ρ)v(0)) = 0, x ∈ R, t > 0 ρ(0, x) = ρ 0 (x), x ∈ R. ( 9
)
In particular, we observe that C(ω η ) → 0 in (8a) and (8b), allowing to recover the classical estimates.

We observe that when the look-ahead distance η → ∞, the non-local problem (1)-( 5) becomes a classical transport equation [START_REF] Chiarello | A non-local traffic flow model for 1-to-1 junctions[END_REF]. Besides the mathematical implications, Corollary 1 may give information on connected autonomous vehicle flow characteristics. Indeed, large kernel supports could account for the information range between connected autonomous vehicles. On the other hand, when the visibility η → 0, the non-local problem (1)-( 5) is again reduced formally to a local problem, but the singular local limit is more challenging from the mathematical point of view. The derivation of this limit was initially conjectured in [START_REF] Amorim | On the numerical integration of scalar nonlocal conservation laws[END_REF] for the one-dimensional scalar case motivated by numerical evidence, later corroborated in [START_REF] Blandin | Well-posedness of a conservation law with non-local flux arising in traffic flow modeling[END_REF]. See [START_REF] Colombo | Recent results on the singular local limit for nonlocal conservation laws[END_REF][START_REF] Colombo | On the singular local limit for conservation laws with nonlocal fluxes[END_REF][START_REF] Keimer | On approximation of local conservation laws by nonlocal conservation laws[END_REF][START_REF] Zumbrun | On a nonlocal dispersive equation modeling particle suspensions[END_REF] for more details regarding the analytical derivation of the singular local limit for non-local conservation laws.

Non-local multi-class traffic flow models

In this Section, we consider the following class of non-local systems of M conservation laws in one space dimension, introduced in [START_REF] Chiarello | Non-local multi-class traffic flow models[END_REF]:

∂ t ρ i (t, x) + ∂ x (ρ i (t, x)v i ((r * ω i )(t, x))) = 0, i = 1, ..., M, (10)
where

r(t, x) := M i=1 ρ i (t, x), v i (ξ) := v max i ψ(ξ), ( 11
) (r * ω i )(t, x) := ∫ x+η i x r(t, y)ω i (y -x)dy, (12)
and we assume:

(H)

ω i ∈ C 1 ([0, η i ]; R + ), η i > 0, s.t. ω i ≤ 0, ∫ η i 0 ω i (y)dy = J i . W 0 := max i=1,...,M ω i (0). v max i are the maximal velocities 0 < v max 1 ≤ v max 2 ≤ . . . ≤ v max M . ψ : R + → R + is a smooth non-increasing function s.t. ψ(0) = 1, ψ(r) = 0 for r ≥ 1.
For simplicity, we can consider the function ψ(r) = max{1r, 0}. We couple [START_REF] Chiarello | Global entropy weak solutions for general non-local traffic flow models with anisotropic kernel[END_REF] with an initial datum

ρ i (0, x) = ρ 0 i (x), i = 1, . . . , M. (13)
The model takes into account the distribution of heterogeneous drivers and vehicles characterized by their maximal speeds and look-ahead visibility in a traffic stream.

It is a non-local generalization of the n-populations model for traffic flow described in [START_REF] Benzoni-Gavage | An n-populations model for traffic flow[END_REF] and it is a multi-class version of the one dimensional scalar conservation law with non-local flux proposed in [START_REF] Blandin | Well-posedness of a conservation law with non-local flux arising in traffic flow modeling[END_REF], where ρ i is the density of vehicles belonging to the i-th class, η i is proportional to the look-ahead distance and J i is the interaction strength. The authors consider the following definition of weak solution because solutions to nonlinear conservation laws are in general discontinuous, even if the initial datum is smooth.

Definition 2 A function ρ = (ρ 1 , . . . , ρ M ) ∈ (L 1 ∩ L ∞ )([0,T[ ×R; R M ), T > 0,
is a weak solution of ( 10), ( 13) if

∫ T 0 ∫ ∞ -∞ (ρ i ∂ t ϕ + ρ i v i (r * ω i )∂ x ϕ) (t, x)dx dt + ∫ ∞ -∞ ρ 0 i (x)ϕ(0, x) dx = 0, for all ϕ ∈ C 1 c (] -∞,T[ ×R; R), i = 1, .
. . , M. The main result in [START_REF] Chiarello | Non-local multi-class traffic flow models[END_REF] is the proof of existence of weak solutions to [START_REF] Chiarello | Global entropy weak solutions for general non-local traffic flow models with anisotropic kernel[END_REF], [START_REF] Chiarello | High-order Finite Volume WENO schemes for non-local multi-class traffic flow models[END_REF], locally in time.

Theorem 2 Let ρ 0 i (x) ∈ (BV ∩ L ∞ ) (R; R + ), for i = 1, .
. . , M, and assumptions (H1) -(H3) hold. Then the Cauchy problem [START_REF] Chiarello | Global entropy weak solutions for general non-local traffic flow models with anisotropic kernel[END_REF], ( 13) admits a weak solution on [0,T[ ×R, for some T > 0 sufficiently small.

An upwind numerical scheme

For the numerical simulations in Section 5, we consider the following conservative scheme for the model ( 10)-( 13) introduced in [START_REF] Friedrich | A Godunov type scheme for a class of LWR traffic flow models with non-local flux[END_REF] in the scalar case and generalized in [START_REF] Chiarello | Non-local multi-class traffic flow models[END_REF] to the multi-class case. First of all, we extend ω i (x) = 0 for x > η i . For j ∈ Z and n ∈ N, let x j+1/2 = j∆x be the cell interfaces, x j = ( j -1/2)∆x the cells centers and t n = n∆t the time mesh. We approximate the initial datum ρ 0 i for i = 1, ..., M and the kernel as follows.

ρ 0 i, j = 1 ∆x ∫ x j+1/2 x j-1/2 ρ 0 i (x) dx, j ∈ Z. ω k i := 1 ∆x ∫ (k+1)∆x k∆x ω i (x) dx, k ∈ N, so that ∆x +∞ k=0 ω k i = ∫ η i 0 ω i (x) dx = J i .
We underline that the sum is finite since

ω k i = 0 for some k ≥ N i . Moreover, we set r n j+k = M i=1
ρ n i, j+k for k ∈ N and

V n i, j := v max i ψ ∆x +∞ k=0 ω k i r n j+k , i = 1, . . . , M, j ∈ Z. (14)
We obtain the following upwind scheme:

ρ n+1 i, j = ρ n i, j -λ ρ n i, j V n i, j+1 -ρ n i, j-1 V n i, j , (15)
where we have set λ = ∆t ∆x .

Numerical tests

In this section we show some numerical simulations of the model ( 10)-( 13) for M = 3. The space step is ∆x = 0.001. We use the numerical scheme in Section 4 to approximate the solution of ( 10)-( 13).

Test 1: impact of connected autonomous trucks on a circular road

We study the impact of connected autonomous trucks on road traffic performances on a circular road. For this reason, we consider the space interval [-1, 1] with periodic boundary conditions at x = ±1. In this case model ( 10)-( 13) reads

                         ∂ t ρ 1 (t, x) + ∂ x ρ 1 (t, x)v max 1 ψ((r * ω 1 )(t, x)) = 0, ∂ t ρ 2 (t, x) + ∂ x ρ 2 (t, x)v max 2 ψ((r * ω 2 )(t, x)) = 0, ∂ t ρ 3 (t, x) + ∂ x ρ 3 (t, x)v max 3 ψ((r * ω 3 )(t, x)) = 0, ρ 1 (0, x) = α (0.5 + 0.3 sin(5πx)), ρ 2 (0, x) = β (0.5 + 0.3 sin(5πx)), ρ 2 (0, x) = (1 -α -β) (0.5 + 0.3 sin(5πx)), (16)
with

ω 1 (x) = 1 η 1 , η 1 = 0.3, ω 2 (x) = 2 η 2 1 - x η 2 , η 2 = 0.05, ω 3 (x) = 2 η 3 1 - x η 3 , η 3 = 0.1, ψ(ξ) = max {1 -ξ, 0} , ξ ≥ 0, v max 1 = 0.8, v max 2 = 1.3, v max 3 = 0.8.
Above ρ 1 represents the density of autonomous trucks, ρ 2 the density of humandriven cars and ρ 3 the density of human-driven trucks. α ∈ [0, 1] is the penetration rate of autonomous trucks and β ∈ [0, 1] is the penetration rate of human-driven cars. Figure 1 displays the traffic dynamics in the case α = 0.3 and β = 0.5, and Figure 2 displays the traffic dynamics in the case α = 0 and β = 0.5 (no presence of autonomous trucks). We observe that oscillations are reduced if autonomous trucks are not present on the road, see Figure 3.

Test 2: stretch of straight road with autonomous trucks

In this test case, we propose a stretch of road with human-driven cars and trucks, and autonomous trucks. We consider the space domain [-1, 1], imposing absorbing boundary conditions at the boundaries x = ±1. The initial conditions for the equation [START_REF] Chiarello | Global entropy weak solutions for general non-local traffic flow models with anisotropic kernel[END_REF] with M = 3 are the following. Above ρ 1 represents the density of autonomous trucks, ρ 2 the density of humandriven trucks and ρ 3 the density of human-driven cars. α 1 ∈ [0, 1] is the penetration rate of autonomous trucks. Figure 4 displays the traffic dynamics in the case α 1 = 0.1 and α = 0.4. In this case, the presence of autonomous trucks reduces the vehicle densities during the overtaking phase, see Figure 4.

Fig. 1

 1 Fig. 1 Density profiles corresponding to the non-local problem (16) with α = 0.3 and β = 0.5 at different times.

ρ 1 ( 1 , η 1 = 0. 5 , v max 1 = 0. 8 ,Fig. 2

 1115182 Fig. 2 Density profiles corresponding to the non-local problem (16) with α = 0 and β = 0.5 at different times.

Fig. 3

 3 Fig.3 (t, x)-plots of the total density r(t, x) = ρ 1 (t, x) + ρ 2 (t, x) + ρ 3 (t, x) computed with the upwind scheme[START_REF] Colombo | Recent results on the singular local limit for nonlocal conservation laws[END_REF], corresponding to different penetrations rate of autonomous trucks: (a) mixed autonomous/human-driven traffic, (b) fully human-driven traffic.

ρ 1 (

 1 0, x) = (1 -α 1 ) χ [-0.6,-0.1](x) , ω 2 (x) ) = max {1 -ξ, 0} , ξ ≥ 0.

Fig. 4

 4 Fig. 4 (t, x)-plots of the total density r(t, x) = ρ 1 (t, x) + ρ 2 (t, x) + ρ 3 (t, x) computed with the upwind scheme (15), corresponding to different penetrations rate of autonomous trucks: (a) α 1 = 0.1, (b) α 1 = 0.4.

Acknowledgements

The author wishes to thank Paola Goatin for her suggestions and valuable discussions.