Computation of Electrical Conductivities of Aqueous Electrolyte Solutions: Two Surfaces, One Property

In this work, we computed electrical conductivities under ambient conditions of aqueous NaCl and KCl solutions by using the Einstein–Helfand equation. Common force fields (charge q = ±1 e) do not reproduce the experimental values of electrical conductivities, viscosities, and diffusion coefficients. Recently, we proposed the idea of using different charges to describe the potential energy surface (PES) and the dipole moment surface (DMS). In this work, we implement this concept. The equilibrium trajectories required to evaluate electrical conductivities (within linear response theory) were obtained by using scaled charges (with the value q = ±0.75 e) to describe the PES. The potential parameters were those of the Madrid-Transport force field, which accurately describe viscosities and diffusion coefficients of these ionic solutions. However, integer charges were used to compute the conductivities (thus describing the DMS). The basic idea is that although the scaled charge describes the ion–water interaction better, the integer charge reflects the value of the charge that is transported due to the electric field. The agreement obtained with experiments is excellent, as for the first time electrical conductivities (and the other transport properties) of NaCl and KCl electrolyte solutions are described with high accuracy for the whole concentration range up to their solubility limit. Finally, we propose an easy way to obtain a rough estimate of the actual electrical conductivity of the potential model under consideration using the approximate Nernst–Einstein equation, which neglects correlations between different ions.


I. INTRODUCTION
Electrolyte solutions are ubiquitous in nature, where ions can play a key role in many living organisms 1 .Electrolytes are also important in many other fields such as battery technology [2][3][4] and desalination processes 5 .Electrolyte solutions have always been the subject of scientific interest 1,[6][7][8][9] and computer simulations can be a valuable tool for the study of complex phenomena related to these electrolyte solutions in combination with experimental studies.In the 1970s, Heizinger, Vogel, Singer and Sangster [10][11][12][13][14] published the first simulation studies of ionic systems.However, the computational cost of the simulations and the lack of suitable force fields for water and electrolytes did not allow major advances until the recent years.
Computational simulations are a useful tool for studying phenomena elusive to experiments and for predicting properties of interest.Nevertheless, a suitable force field for the studied system is needed.In the case of aqueous electrolyte solutions, force fields for water and ions are necessary 15 .In the case of water, the first force field was proposed by Bernal and Fowler in 1933 15 .Fifty years later, Jorgensen and coworkers started to develop new potential models, such as the TIP3P 16 , the TIP4P 16 and the TIP5P 17 force fields.At the same time, the popular SPC/E force field developed by Berendsen and coworkers 18 .Later, in the 2000s, the knowledge gained from the aforementioned models allowed the development of two of the water models that best reproduce a wide range of properties, which are TIP4P-Ew 19 and TIP4P/2005 20 .In fact, the TIP4P/2005 potential is able to reproduce a variety of properties such as densities, viscosities, and the temperature of the maximum in density (TMD) [21][22][23] .There is no classical model able to reproduce all properties of pure water 21,24 .One option to improve results of the previous mentioned force fields is to use polarizable models, such as the HBP 25 , the MB-Pol 26 or the BK3 27 .However, these force fields are between three and ten times more computationally expensive than the non-polarizable force fields and are also not able to reproduce certain properties simultaneously, such as TMD and melting temperature 24 .Besides water models, ion-ion and ion-water interaction have to be described to study electrolyte solutions.It is in the recent years that a large variety of force fields for salts have been proposed  . Amongthe most popular force fields for ions, we can find the one proposed by Joung and Cheatham 42 (JC), which includes all the alkali halides.These can be used in combination with three different water models-namely; TIP3P, TIP4P-Ew and SPC/E)-the JC-SPC/E being the one that provides overall better results.The other popular force field that is widely used in the literature is the one developed by Smith and Dang 34 (SD) in combination with SPC/E water 18 .Although these force fields have been quite successful in describing many properties of electrolytes (e.g.densities, structure), they fail in describing properties such as solubilities, viscosities and activity coefficients 62 .
In an attempt to overcome the limitations of current force fields for electrolytes, the idea of using scaled charges for the ions was suggested.The concept of scaled charges (i.e., assign a charge smaller than one for monovalent ions) arises from the work of Leontyev and Stuchebrukhov [63][64][65][66][67][68] who proposed a charge of ±0.75 (in electron units) for ions in solution (this was also denoted as Electronic Continuum Correction ECC).The use of scaled charges has also been proposed by Kann and Skinner 69 .However, in this case the value of the scaled charge is selected in such a way that the potential of mean force between ions at infinite dilution and large distances should be the same in experiments and in the force field (so that the model recovers the experimental Debye-Huckel law at infinite dilution).As the potential of mean force depends on the dielectric constant of the water model, so does the value of the scaled charge leading to a value of ±0.85 in the particular case of water when described by the TIP4P/2005 model.The use of scaled charges for ions in solution has undergone a significant expansion in the last years.Different groups have proposed new force fields with scaled charges.including those of Jungwirth and coworkers [70][71][72][73][74] , Barbosa and coworkers 75,76 , Li and Wang 77 and Bruce and van der Vegt 78 .Other authors have proposed the use of different scaled charges for the cation and anion (charging the surrounding water molecules to maintain the system electroneutral) [79][80][81][82][83] .Breton and Joly also studied the effect of introducing scaled charges for studying interfacial properties 84 .In this context, we developed a model for NaCl based on scaled charges 48 .Later, we considered a larger number of salts and we proposed the Madrid-2019 force field 58,59 which includes all the possible alkali halides, some divalent salts (Mg 2+ and Ca 2+ ) and sulphates (SO 2− 4 ).We have shown in previous works that this force field is able to reproduce different properties of interest, such as the salting out effect of methane 85 , the TMD of different salt solutions 86 , the freezing depression of ice in presence of different electrolytes 87 or different properties of seawater 88 .Although scaled charges improve the results in the majority of properties with respect to unit charge models, there is no unique value of the scaled charge that describes all properties correctly.As we have recently shown, the scaled charge can be taken as a fitting parameter depending on the property that one wants to reproduce 61 .Transport properties are among the most interesting properties that can be studied by simulation and that traditional force fields of electrolytes have never been able to reproduce correctly. 89.In our recent work, we proposed the Madrid-Transport force field 60,61 , which uses an scaled charge of q=±0.75 and that is able to reproduce transport properties, such as the viscosities and diffusion coefficients of water and ions in the whole concentration range.This force field has also been able to reproduce transport properties in the presence of hydrogen 90 .It should be mentioned that the introduction of scaled charges, improves a number of properties of electrolytes but deteriorates the value of the free energy of solvation (although it can be corrected via theoretical corrections 85 ).
In this work, we want to analyze in detail the quality of the predictions for electrical conductivities of force fields using either integer or scaled charges.To the best of our knowledge, such a detailed comparison has never been presented before.As we will show in the next section, electrical conductivities can be calculated with the Einstein-Helfand (EH) equation [91][92][93][94][95][96][97] .In a preliminary but pioneering work, Lyubartsev and Laaksonen 98 evaluated the electrical conductivities of NaCl aqueous solutions at different concentrations with the flexible SPC model for water and ions described as charged LJ particles by using the Green-Kubo (GK) equation 99 (which is strictly equivalent to the Einstein-Helfand relation).Due to the computational cost of evaluating the conductivities in this way 100 , many authors tend to calculate the electrical conductivities by using the Nernst-Einstein equation (i.e., neglecting the ion-ion correlations) [101][102][103][104][105][106] .Although this is cheaper from a computational point of view, the results are not exact (as this is an approximation) and overestimate the real conductivities of the model.Electrical conductivities have been accurate and extensively calculated (through GK or EH) by different authors for ionic liquids 91,92,[107][108][109][110] .In the case of aqueous electrolyte solutions, there are a few studies in which the conductivities were properly evaluated 93,111,112 .Martí, Guardia and coworkers 111 calculated electrical conductivities of NaCl solutions at different concentrations using the Smith and Dang 34 ion force field with SPC/E water 18 (SD-SPC/E).Shao et al. 112 also calculated the conductivities of NaCl solutions but using in this case the Joung and Cheatham 42 force field in combination with SPC/E water (JC-SPC/E).In fact, these authors showed interesting results about the existence of finite size effects when calculating conductivities with the Nernst-Einstein relation.There are also semi-empirical fitted models which try to reproduce the electrical conductivities of NaCl solutions in solvent mixtures, such as water-propylene carbonate 113 or water-monoethylene glycol 114 .Other authors have also rigorously computed electrical conductivities for molten salts 115,116 .Nevertheless, there is no comprehensive Molecular Dynamics (MD) study of the performance of different force fields of ions and water for reproducing the electrical conductivities of NaCl solutions in water.
The main purpose of this work is to provide a benchmark to calculate electrical conductivities of aqueous electrolyte solutions, to show that there is a force field able to reproduce the experimental conductivities of NaCl and KCl solutions and finally, from a deeper perspective, we want to demonstrate that to reproduce conductivities the two surfaces present in water have to be simultaneously described.We will properly evaluate the electrical conductivities of different well-known ion force fields in combination with different water force fields.Besides, in this work, we shall introduce a new "conceptual" strategy to determine electrical conductivities.As we have mentioned in the past, water (and aqueous solutions) has two surfaces: the potential energy surface (PES) that describes the energy of each configuration of the system and the dipole moment surface (DMS) that describe the dipole moment of each configuration [117][118][119] .The key idea is that these two surfaces can be described by different fitting parameters (i.e., with different charges in our case).In this work, we use a force field with scaled charges to perform the Molecular Dynamics simulations of the aqueous electrolyte solutions.Only the PES is needed to perform these simulations as we are using the linear response theory, that evaluates transport properties (electrical conductivities in this work) by analyzing the fluctuations of the system when at equilibrium (in absence of an electric field in the case of electrical conductivities).We will use the Madrid-Transport force field (q=±0.75)for the PES as it is able to describe other properties of electrolyte solutions accurately (i.e., viscosities and diffusion coefficients).The obtained trajectories are analyzed employing unit charges to calculate the electrical conductivities (as integer charges instead of partial charges provide better representation of the DMS).In this way, we will show how to evaluate conductivities of aqueous electrolyte solutions with a new methodology which yields to reproduce the experimental conductivities of NaCl and KCl solutions by using a scaled charge force field.

II. METHODOLOGY
Conductivities can be calculated using the Einstein-Helfand relation [91][92][93][94][95] : where V is the system volume, k B the Boltzmann constant, T the temperature, ⃗ r i (t) and ⃗ r i (0) the position of the i th particle at time t, and the ⟨[⃗ r i (t) − ⃗ r i (0)] 2 ⟩ term is the mean square displacement (MSD).Taking into account that the dipole moment of the system can be defined as: we can obtain the following equation for the conductivity in function of the mean square dipole displacement of the system: This way, the conductivity can be easily obtained from the slope of the mean square dipole displacement versus time (in the Supporting Information we show this plot for different initial seeds of the MD simulations).
Eq. 1 can be rewritten as: Please note that q i is equivalent to e•z i We define particles between 1 and N 2 as cations, and particles between N 2 + 1 and N as anions.Thus, we can define for a 1:1 electrolyte the Onsager coefficients 1 (Λ ij ) for the different interactions between cations (+) and and anions (-): where N is the total number of ions.Combining Eq. 4 with Eqs.5-8, the electrical conductivity can be calculated from: Eq. 9 is strictly equivalent to Eqs. 1 and 3. Note that in all these equations, we only include the charge and positions of the ions.We do not consider water (solvent) molecules, because it is a neutral molecule which does not contribute to the electrical conductivity.
Since the evaluation of electrical conductivities by the previous methodology can be computationally demanding, some authors use the approximate Nernst-Einstein (NE) equation 120 , which relates the electrical conductivity to the self-diffusion coefficients of the ions.
where q is the charge of the ions, ρ is the number density of the salt and D + and D − are the self-diffusion coefficients of the cation, and anion respectively.The NE approximation assumes that the only terms which contribute to the conductivities are those of the particle with itself.Note that whereas there are N 2 terms in the rigorous expression of the conductivity (Eq.4), there are only 2N terms when using the NE relation (Eq.4).The NE relation, thus, can be obtained from Eq. 1 by assuming that the displacement of ions is independent (i.e., ⟨r i • r j ⟩=0 for i̸ =j).Making this assumption, electrical conductivities can be easily obtained for equimolar salts with Eq. 10.For the evaluation of conductivities with the NE equation, we have calculated the self-diffusion coefficients of both Na + and Cl − by using the Einstein relation (Eq.11): where ⃗ r i (t) and ⃗ r i (0) are the position of the i th particle at time t and 0, and the ⟨[⃗ r i (t) − ⃗ r i (0)] 2 ⟩ term is the MSD.All diffusivities in this work are corrected using the hydrodynamic corrections of Yeh and Hummer 121,122 which are described as where D i is the diffusion coefficient with the applied corrections of Yeh and Hummer, D MD i is the diffusion coefficient initially obtained by simulations, ξ is a dimensionless constant equal to 2.837, η is the computed viscosity at the studied concentration (which is shown to exhibit no finite size effects 121,123,124 ), and L is the length of the simulation box.

III. SIMULATION DETAILS
Electrical conductivity is not a property for which there are many previous results and it is interesting to see whether two groups using different software packages (GROMACS vs LAMMPS), slightly different methodologies (the individual calculation of the Onsager coefficients versus the global expression), and different system sizes yield the same result for a specific model.Since we also aim to use the results of this work as a benchmark for people computing electrical conductivities in the future, and we want to ensure that the good agreement with experiment that is observed here is true regardless of the MD software or post-processing program used for the computation.As one group is located in Madrid (UCM) and the other in Delft (TU Delft), they will be denoted as Madrid and Delft groups.Concentrations of the salts will be given in molality units, so that a solution with a concentration 1 m corresponds to 1 mol of salt per kilogram of water.All the results of this work correspond to room temperature and atmospheric pressure (i.e., 298.15K and 1 bar).Note also that error bars of all the results are calculated by both groups as the standard deviation the property obtained in each run using different initial seeds divided by the square root of the number of runs.

A. Madrid Group
All MD simulations were performed with GROMACS 125,126 (version 4.6.7).The leap-frog integrator algorithm 127 with a time step of 2 fs was used.We have also employed periodic boundary conditions in all directions.Temperature and pressure were kept constant using the Nosé-Hoover thermostat 128,129 and Parrinello-Rahman barostat 130 , both with a coupling constant of 2 ps.For electrostatics and van der Waals interactions, the cut-off radii were fixed at 1.0 nm and long-range corrections in the energy and pressure were applied to the Lennard-Jones part of the potential.The smooth PME method 131 was used to account for long-range electrostatic forces.Water geometry was maintained using the LINCS algorithm 132,133 .To compute conductivities, we have simulated systems of 4440 water molecules and the corresponding number of ions for the desired concentration (e.g., 80 NaCl molecules for a concentration of 1 m as shown in Tables 1 and 2).We have performed an initial N pT simulation of 20 ns to accurately calculate the volume of the system.After that, using the average volume obtained in the N pT simulation we have carried out five independent runs in N V T ensemble.Runs of 200 ns were performed for the lowest concentration (i.e., 1 m) and of 120 ns for the higher concentrations (2, 4, and 6 m).Thus, typically around 600 ns (5x120) or 1 µs (5x200) are needed to compute electrical conductivities of each model and thermodynamic state.
The electrical conductivities were obtained from fitting the mean square dipole displacement (Eq.2) versus time between 50 and 1000 ps as shown in Eq. 3 (in the SI we also provide the results from fitting the data between 50 and 2000 ps).

B. Delft Group
MD simulations are carried out on the Delft-Blue super-computer at TU Delft 134 with the Largescale Atomic/Molecular Massively Parallel Simulator (LAMMPS: version August 2018) 135 .Periodic boundary conditions are imposed in all directions, and the velocity-Verlet algorithm is used with a time step of 2 fs.The Nosé-Hoover thermostat and barostat [128][129][130] are set with coupling constants of 0.1 and 1 ps, respectively.The SHAKE algorithm is used to fix the bond lengths and angles of water 135,136 .A cut-off of 1.0 nm is used for both the Lennard-Jones and electrostatic potentials.Long-range electrostatic interactions are modelled using the particle-particle particle-mesh (PPPM) method 137,138 with a relative error 139 of 10 −5 .Analytic tail corrections 137 are applied to the Lennard-Jones interactions for both energies and pressures.Initial configurations are created using PACKMOL (v20.3.1) 140.To compute electrical conductivities, self-diffusivities, and shear viscosities, the OCTP plugin 141 is used.In this plugin, the Einstein relations are used in combination with the order-n algorithm 137,142 to compute transport properties.All details on the OCTP plugin can be found in Ref. 141 The approach to evaluate conductivities is based on the computation of the Onsager coefficients (Λ ij ) for the cation-cation, anion-anion and cation-anion interactions independently as shown in Eqs.4-9.These equations are used to compute the exact electrical conductivities accounting for ion-ion correlations.In the SI we have collected the results for the individual contributions to the electrical conductivities (i.e., σ ++ ,σ +− , σ −+ , and σ −− ) computed from each individual Onsager coefficient.To evaluate the electrical conductivities, the system sizes were of 555 or 1000 water molecules (the corresponding number of ion molecules is dictated from the molality of each system) as we have listed in Tables 1 and 2. To accurately compute the average volume of the simulation box, simulations of 20 ns in the N pT ensemble were initially carried out (10 ns equilibration runs followed by 10 ns production runs).The self-diffusivities, viscosities, and Onsager coefficients are calculated from production runs of 200 ns in the N V T ensemble.Three different simulations were carried out with different initial velocities for all molalities to obtain statistics.Thus, a total simulation time of ca.600 ns is required to compute electrical conductivities of each model and thermodynamic state.

A. Electrical conductivities of popular force fields
In previous studies, we have shown that popular force fields (that use integer charges for the ions) are not able to reproduce transport properties, such as viscosities, diffusion coefficients of water in salt solutions, and selfdiffusion coefficients of ions 58,59,61 .Here we investigate if these models also fail in accurately predicting another important transport property, namely the electrical conductivity.In Fig. 1, results of the electrical conductivity of NaCl obtained in this work for two popular force fields (i.e., Joung and Cheatham 42 and Smith and Dang 34 ) that use integer charges combined with the SPC/E water model are compared to experiments.Note that contrary to previous studies that only focused on concentrations up to 4 m, in this work we have evaluated the conductivities in the whole concentration range, i.e., up to the experimental solubility limit of NaCl (6.1 m).Our results at 4 m for the JC-SPC/E are in excellent agreement with those obtained by Shao et al. 112 using the Green Kubo formalism.Results presented in Fig. 1 were obtained by the Madrid group.It is clear that electrical conductivities are underestimated with respect to experiments by the models using integer values for the charge.The SD-SPC/E is slightly more accurate compared to the JC-SPC/E, but also underestimates the electrical conductivities.Thus, it is evident that these two popular force fields for NaCl are not able to reproduce electrical conductivities.These models overestimate the experimental viscosities of NaCl solutions 61 .It is not surprising that electrical conductivities are underestimated, as intuitively one would expect that an overestimate of the viscosity would lead to an underestimate of the diffusion coefficient of the ions and, therefore, to an underestimate of the the electrical conductivity.

B. Electrical Conductivities of Scaled Charge models
Unit charge force fields cannot reproduce electrical conductivities of aqueous solutions, but in previous works we have demonstrated that scaled charge models significantly improve the description of transport properties as viscosities and self-diffusion coefficients 60,61 .Thus, we have now employed the Madrid-2019 force field to observe whether we can improve the results of unit charge force fields.This force field uses scaled charges for the ions (in particular 0.85).In Fig. 1, we show the results for the conductivities of the Madrid-2019 model (red empty squares).Surprisingly we do not obtain better results than those of unit charge models.This was also observed by Gullbrekken et al. 144 in their recent work in which they studied the electrical conductivities by using an integer charge force field but then when scaling the charges  the do not observe differences.Regarding Madrid-2019 results, one may wonder how is it possible than a force field which yields better results for viscosities and diffusion coefficients of aqueous solutions does not perform equally well for electrical conductivity.We shall provide now a possible explanation for this puzzling behavior.
As discussed in detail previously 117 , water has two different surfaces, the Potential Energy Surface (PES) and the Dipole Moment Surface (DMS).In the absence of an external macroscopic electric field, all properties of a system can be determined in computer simulations from the PES.The DMS surface is not needed to determine any property of a system when not applying an external electric field.However, certain properties describe the response of a system to an external electric field.
In particular, the dielectric constant and the electrical conductivity are response functions of this type.Obviously, these properties are relevant only when an external electric field is applied to the system.Therefore, to determine these response functions in computer simulations, it is required to describe both the PES and the DMS.A clarification is in order now.The PES simply gives the energy of a system provided the positions of all the nuclei of the system.The PES only depends on the position of the atoms, and does not depend on any macroscopic property as for instance the viscosity, or the dielectric constant.Sometimes it is stated that the dielectric constant enters in the description of the PES in the case of electrolytes.This is not correct.Even for electrolytes one simply needs to know the position of the atoms to determine the energy of the system, and the value of the dielectric constant is not needed.The origin of this confusion arises from the fact that at infinite dilution, and when r tend to ∞, the potential of mean force w(r) between two ions obtained defined from the radial distribution function g(r) using the expression: can be obtained from the knowledge of the dielectric constant ϵ r as follows: But the potential of mean force is not the PES , and besides this expression is only valid in the Debye-Huckel limit (i.e., at infinite dilution of electrolyte and infinitely large distances).The summary is that the dielectric constant does not enter in the description of the PES, and in the absence of an external electric field all properties of an aqueous electrolyte solution can be obtained from the PES and the knowledge of the DMS is not needed.
Imagine that a model is able to describe correctly the viscosities and the diffusion coefficients.This is an indicator that the PES is described correctly.Imagine now that the electrical conductivity is not well described.How to solve this paradox?The answer is rather simple: the PES is well described but the DMS is not well described.Often partial charges are used to describe the PES and these charges are also used to describe the DMS.We have suggested sometime ago that the charges that are good to describe the PES may not be suitable to describe the DMS 117 .We have provided some indications by analyzing the behavior of the dielectric constant of water.Jorge and coworkers 119,145 have shown that the same is true for the dielectric constant of alcohols, and Bowman and coworkers followed up on this idea 146,147 .
Here, we shall use different charges to describe the PES and the DMS of electrolytes in water.In particular, we shall use scaled charges for the ions in water for the PES while we shall use non-scaled charges (i.e., integer charges) to describe the DMS.The way to implement this idea is rather simple.Since we are using linear response theory (so that the electrical conductivity is computed by simulating the system in the absence of the electric field), we shall perform MD simulations to obtain the trajectories using the scaled charge of the ions (i.e., ± 0.85 for the particular case of the Madrid-2019 force field).In sharp contrast, when the trajectory is analyzed using Eqs.1-9, integer charges are used for the ions.The results from this approach are presented for the Madrid-2019 force field in Fig. 1 (orange empty squares).As clearly shown, the computed conductivitiesare closer to the experimental data.These results showcase that a better description is obtained when simultaneously describing both surfaces the PES and DMS.However, we do not obtain a perfect agreement with the experimental electrical conductivities of the aqueous solutions.This is not entirely surprising, as the Madrid-2019 force field (that uses a scaled charge of ±0.85) improves the description of transport proper-ties of electrolytes in water, but is not able to yield quantitative agreement with experiments.To this end, we have recently proposed a force field for NaCl and KCl (denoted as Madrid-Transport) that is able to predict transport properties of these electrolyte solutions with excellent agreement to experiment 60,61 .Therefore, we shall compute the electrical conductivities of NaCl and KCl solutions using the Madrid-Transport force field (that uses an scaled charge of ±0.75).We shall implement the main idea of this work, namely, to use scaled charges to obtain the trajectories and integer charges to describe the DMS (i.e., using integer charges in Eqs.1-9).In Fig. 2, we present both properties (viscosities and self-diffusion coefficients of water) in the whole concentration range of each salt up to the experimental solubility limit.Results are independently calculated by two different research groups: Madrid (blue) and Delft (red) and are consistent within the error bars.The results are in good agreement with experiments showing that these two transport properties that are obtained from the PES are described satisfactorily by the Madrid-Transport force field.
In Fig. 3 the results for the electrical conductivities of NaCl and KCl solutions using the Madrid-Transport force field are shown.Note that these results are computed with our novel approach (i.e., using scaled charges (q=±0.75) to describe the trajectories of the system and unit charges to compute the electrical conductivities).We show the results obtained from two different research groups: Madrid (blue) and Delft (red).Both groups have adopted the same approach to calculate conductivities (i.e., employing the EH equations) with some minor differences.Madrid has evaluated the conductivity from the mean square dipole displacement, taking into account all the interactions between ions (cation-cation, anion-anion and cation-anion).Delft has evaluated the Onsager coefficients independently for the cation-cation, anion-anion, and cation-anion interactions, and then summed up all the contributions.Both approaches are equivalent and, thus, should yield identical results.It is important to note that both groups have used different software, constraint algorithms, system sizes, simulation times, and fitting methods (see Simulation Details).Even so, the results for the conductivities obtained by both groups are equal within the error bars.Note also that despite the system sizes used by Delft and Madrid groups are different (between 555 and 4440 water molecules), the electrical conductivities are in agreement, showing that no finite size effects in the computation of electrical conductivities are observed.Nevertheless, since both self and collective (Maxwell Stefan and Fick) diffusivities exhibit significant finite size effects 121,122,150,151 , a thorough investigation for electric conductivities should be also performed.This is particularly important for small concentrations (i.e., below 1 m).This was also previously studied by Shao et al. 112 for the JC-SPC/E model concluding that there were not finite size effects when using the Green-Kubo equation.In Fig. 3, we can observe that the conductivities obtained by both groups for the Madrid-Transport force field reproduce the experimental conductivities of both NaCl and KCl aqueous solutions for the whole concentration range.Note that the electrical conductivity of KCl is significantly larger than that of NaCl at the same concentration.This difference is correctly described by the Madrid-Transport force field.Thus, the Madrid-Transport force field and the use of scaled charges for the trajectories and integer charges for computing conductivities (i.e., describing the PES with scaled charges and the DMS with integer charges) allows for the first time to correctly reproduce experimental conductivities (and other transport properties as viscosities and water diffusion coefficients).In Tables 1 and 2, we have collected the computed conductivities for each system (by both groups) along with the conductivities obtained by using the NE equation (Eq.10) with and without applying the finite-size corrections to the self-diffusion coefficients.

C. Electrical Conductivities by using the Nernst-Einstein equation
Although the correct way to calculate conductivities is to use the Green-Kubo or the Einstein-Helfand equations.However, the NE equation (Eq.10) is widely used due to its simplicity.In Tables 1 and 2, the electrical conductivities predicted by the approximate NE formula are also presented.As can be seen, the NE relation overestimates the true conductivity of the model as obtained from the EH equations.Not surprisingly, the NE equation does not provide the exact value of the electrical conductivity of the model (see the Supporting Information).The reason for that is that it neglects correlations between different ions.The deviations increase with the concentration but they are clearly visible even at a concentration of 1 m.Our results suggest that correlations between different ions tend to decrease the value of the electrical conductivity.The conclusion from the results of Tables 1 and 2 is that NE should not be used to estimate electrical conductivities as it provides incorrect results.If one wants to compute the true conductivity of a force field, either the GK or EH expressions should be used.
Despite being inaccurate, the NE expression is often used in many papers to obtain electrical conductivities.There are two main reasons for this.The first one is that the calculation of self-diffusion coefficients is relatively easy and computationally cheap.Thus, if the NE formalism is used, the electrical conductivities are obtained with no additional computational cost.Diffusion coefficients have good statistics as one can accumulate the statistics of each individual ion, while the EH is expensive as it is a  for different molalities (m in units of mol salt kgwater -1 ) and system sizes.All simulations were performed at 1 bar and 298.15 K, using the Madrid-Transport model.The number of water molecules (nW) and NaCl molecules (ns), the corresponding densities (ρ in units of kg m −3 ) and viscosities (η in units of mPa•s) are shown for all molalities.Additional electrical conductivities computed using the Nernst-Einstein with (σNE+YH in units of S•m −1 ) and without (σNE in units of S•m −1 ) Yeh-Hummer finite-size corrections 121,122,124  global property and each configuration contributes a single value of the correlation function.In short, diffusion coefficients of ions can be obtained with good accuracy in runs of 20 ns, whereas ca.600 ns are needed to have reasonable statistics of the electrical conductivity.Additionally, the computation of the diffusion coefficient is implemented in many MD packages, whereas this is not the case of the electrical conductivity.The second reason why NE is so popular is because most of the force fields tend to underestimate significantly the electrical conductivity compared to the experiments when computed rigorously from the EH formalism.However, since NE electrical conductivities are much larger they tend to be in better agreement with the experimental results.This creates the paradox that despite NE tends to overestimate the electric conductivities, in many cases is closer to experimental data, and thus, there is a resistance to abandon its use.Obviously, this apparent agreement arises TABLE 2. Computed electrical conductivities (σ in units of S•m −1 ) of aqueous KCl solutions from the EH relations (Eq. 3) for different molalities (m in units of mol salt kgwater -1 ) and system sizes.All simulations were performed at 1 bar and 298.15 K, using the Madrid-Transport model.The number of water molecules (nW) and KCl molecules (ns), the corresponding densities (ρ in units of kg m −3 ) and viscosities (η in units of mPa•s) are shown for all molalities.Additional electrical conductivities computed using the Nernst-Einstein with (σNE+YH in units of S•m −1 ) and without (σNE in units of S•m −1 ) Yeh-Hummer finite-size corrections 121,122,124  FIG. 3. Electrical conductivities as a function of NaCl concentration obtained with the Madrid-Transport force field in independent research groups: Madrid (blue) using Eq. 3 and Delft (red) using Eq. 9, for NaCl (circles) and KCl (squares) aqueous solutions at temperature T = 298.15K, pressure p = 1 bar.Experimental results have been taken from Ref. 143 from a cancellation of errors (i.e., a poor force field along with a poor way of computing the actual conductivity of the force field can provide good agreement with experiments).This is illustrated in Fig. 4, where the electrical conductivities of the JC-SPC/E and SD-SPC/E obtained from NE and from EH are compared to experiment.As it can be seen, the agreement is better with NE.As we have shown NE does not describe correctly the electrical conductivity of the force field so that the apparent improvement is obtained from fortuitous cancellation of two errors (the force field and the way to compute the electrical conductivity).Another "apparent" advantage of the NE formalism is that since diffusion coefficients are quite sensitive to the size of the system, one can often find a system size for which the agreement with experiment is excellent.This adds another degree of freedom for "finetuning" the final value of electric conductivity.However, this should not be accepted when striving for accurate computation of properties.This becomes even more pronounced by the fact that the system size dependency of the electrical conductivity when computed from the EH formalism is quite small as it has been shown in this work (i.e., compare the results from Madrid and from Delft) and also shown in other works 112 .
After this work, we strongly advise against using NE when the one can actually compute the conductivity using rigorous ways.NE can be used when one has no access to such a computation or when a quick estimation of the order of magnitude of the conductivity is needed.However, even in such cases, the researcher should keep in mind that the NE value is overestimated respect to the real conductivity.However, we will suggest an approximate (and not rigorous) way of at least correcting NE results.The idea is simple.We have analyzed the typical ratio between the electrical conductivities obtained rigorously (i.e., EH relation) and those obtained from the approximate NE formalism, and have analyzed whether it is approximately constant for different force fields and concentrations.When computing the electrical conductivity from NE, we have used the finite size corrected diffusion coefficients.In this way, the NE conductivities computed here are determined using one of the most correct approaches to estimate diffusion coefficients in the thermodynamic limit.In Fig. 5, we show the ratio between the conductivities obtained by EH divided by the conductivities calculated with NE as a function of the concentration.Notice that for models with different charges (JC-SPC/E, Madrid-2019 and Madrid-Transport), for all concentrations and even for both studied salts (NaCl and KCl), the conductivities using the EH equation are about 30% lower than when using the NE equation.Thus, a rough approximation of the EH conductivities can be obtained simply by employing: where σ EH is the conductivity rigorously computed by using the EH relation, and σ N E+Y H is the approximate conductivity calculated with the NE equation (after including Yeh-Hummer corrections to the values of the diffusion coefficients).If one wants to roughly estimate the correct conductivities of a model (i.e., evaluated by the EH relation), one can only employ the NE equation and then apply our rule described in Eq. 15.In this way, one can obtain an approximate (with a typical error of about 5-8%) but still reasonable estimate of the true conductivity of the force field under consideration from the initial guess provided by the NE relation.In fact, regarding the recent work from Gullbrekken et al. 144 , this rule also works properly for their results.Note also that this scaling factor in Eq. 15 (0.7) is empirical and is not related to the charge of ions in the force field.In any case, we recommend to evaluate properly the conductivities by using the GK or EH relations without any approximation to obtain rigorously the correct conductivity of the force field.

V. CONCLUSIONS
In this work, we have evaluated the electrical conductivities of NaCl and KCl aqueous solutions up to their solubility limit by using the Einstein-Helfand equation.We have computed the electrical conductivities of four different NaCl force fields (i.e., JC-SPC/E, SD-SPC/E, Madrid-2019 and Madrid-Transport).We have shown that force fields with integer charges are not able to reproduce electrical conductivities, as they underestimate them considerably.Note that for these models integer charges are used to describe both the PES and the DMS.We tested a new approach that is based on the more general idea that different charges should be used to reproduce correctly the PES and the DMS.To implement this idea, we used the Madrid-Transport force field for NaCl and KCl (with scaled charges of ±0.75 for the ions) to obtain the trajectories of the system in the absence of an external electric field.The EH equation is used with integer charges for the ions.The basic idea is that scaled charges describe the PES better, and integer charges describe the DMS more accurately, as was suggested some time ago. 117,119By using this approach, we have shown that the Madrid-Transport force field (a model which provides excellent results for transport properties such as viscosities and diffusion coefficients) reproduces the electrical conductivities of NaCl and KCl for the whole concentration range.Certainly, one can argue that it is possible to use a model with q=±1 to describe the PES and with q higher than ±1 for describing the DMS.By using this idea (developed in this work) the agreement with experiments for electrical conductivities will also improve for models that use integer charges for the PES.However, the cost of that is that the transported charge will not be 1 e , a result that has no physical meaning.In addition, the force field with that charge would not reproduce the viscosities and self-diffusion coefficients of the water as in the rest of unit charge models.It seems therefore that to reproduce the DMS one should use the integer charge that actually is transported (1 e) but in contrast, for describing the PES, one can use scaled charges to better describe the relative weight of some configurations with respect to others.What quantum chemistry says about the value of the charge that is indeed transported?Is it correct to assume formal integer charges when computing electrical conductivities?This issue has been discussed in two important papers 115,152 .For simplicity we shall discuss this issue using the Green-Kubo formalism 153,154 but one could also use the totally equivalent Helfand-Einstein relation.The electrictical conductivity can be rigorously computed from ab-initio calculations using the expression: where ⃗ j i is the current density vector which is defined as: Apparently everything seems normal, but now comes the first surprise.The charge Q * i (t) is not a scalar but a time dependent tensor with components given by: where α, β = x, y, z and r i,β refers to the component β of the vector defining the position of ion i and M α is the α component of the dipole moment of the system.It should not come as a surprise that the value of the charge of the ion for the calculation of the conductivity comes from derivatives of the dipole moment of the system which is well defined and not from arbitrary schemes partitioning the electronic density among the different ions of the system (as for instance the Bader method 155,156 ).However, it has been shown 115,152 that one can obtain the rigorous value of the conductivities using a non-time dependent scalar for the charge of ion i (usually denoted as the topological charge q i,top ).This is summarized in the following equations: The conclusion is that although the charge of an ion is a time dependent tensor, there is a non-time dependent scalar value of the charge (denoted as the topological charge) that leads to the correct value of the electrical conductivity.Finally, Grasselli and Baroni 115 , and French et al. 152 have shown that the value of the topological charge is just the value of the formal charge.Quantum chemistry suports the use of the formal integer charge (i.e., 1 e for monovalent ions) in the calculation of the electrical conductivity.One should not use different values of the formal charge when computing electrical conductivities.That is what we have done in this work.However, quantum chemistry does not say anything about the charge that better fits the PES as it should be regarded as a fitting parameter of the force field, and we have shown that for aqueous electrolyte solutions the scaled charge with value 0.75 e provides an excellent description of transport properties.
Finally, we have shown that when using the Nernst-Einstein equation the electrical conductivities are incorrectly calculated due to the neglect of the correlation between different ions, showing discrepancies even at low concentrations.We propose a rule of thumb to obtain a rough estimate of the EH conductivities.The recipe is simple.It consists in multiplying by 0.7 the conductivity obtained from the NE equation (after applying YH corrections to the values of the diffusion coefficients of the individual ions).Nevertheless, this is an approximate correction, and our advice is to use the correct expression (i.e., EH or GK) to obtain rigorously the electrical conductivity of a certain force field.
The results presented in this work independently obtained by two research groups could be useful in the future as benchmark results to be reproduced by groups interested in computing electrical conductivities of electrolyte solutions.The success of the Madrid-Transport force field in reproducing the transport properties of NaCl and KCl solutions can be regarded as workcase example showing how fruitful the idea of using different charges to describe the PES and DMS can be in the future.The community performing classical simulations should benefit from this "mental" flexibility.In fact, the community performing ab-initio calculations is already using in a effective way as they are using different fitting parameters when developing neural networks for the PES and the DMS 157,158 .Why should we not realize that we can do the same with our force fields?

FIG. 2 .
FIG. 2. Transport properties of NaCl and KCl aqueous solutions at different concentrations obtained with the Madrid-Transport force field by Madrid (blue) and Delft (red) groups at at temperature T = 298.15K and pressure p = 1 bar.(a) Viscosities and(b) Self-diffusion coefficients of water (corrected for system size effects).Experimental results have been taken from Refs.148,149

FIG. 4 .
FIG. 4. Electrical conductivities as a function of NaCl concentration obtained with different models studied in this work by two different methodologies (i.e., by using the EH relation and the Nernst-Einstein equation) at T = 298.15K, p = 1 bar.(a) JC-SPC/E (b) SD-SPC/E.Experimental results have been taken from Ref. 143
in an ab initio machine-learning model of water.Proceedings of the National Academy of Sciences 2022, 119, e2207294119.
143.1.Electrical conductivities (computed by the Madrid group) as a function of NaCl concentration obtained with the different models studied in this work using Eq. 3 at temperature T = 298.15K,and pressure p = 1 bar.Note that Madrid-2019 PES (red empty squares) uses scaled charges in both MD simulations and for computing the conductivities, but Madrid-2019 PES+DMS uses scaled charges for the dynamics and integer charges for computing the conductivities.Experimental results have been taken from Ref.143

TABLE 1 .
Computed electrical conductivities (σ in units of S•m −1 ) of aqueous NaCl solutions from the EH relations (Eq. 3) are reported as well.Numbers in parentheses are the uncertainty in the last digit of the results.
are reported as well.Numbers in parentheses are the uncertainty in the last digit of the results.