Study on the free corrosion potential at an interface between an Al electrode and an acidic aqueous NaCl solution through density functional theory combined with the reference interaction site model

We investigated the free corrosion potential at an interface between an Al electrode and an aqueous NaCl solution (NaCl(aq)) under acidic conditions via density functional theory combined with the effective screening medium and reference interaction site model (ESM-RISM). Firstly, the electrode potentials for the anodic and cathodic corrosion reactions were obtained from the grand potential profile as a function of the electron chemical potential at the interface. Thereafter, we determined the free corrosion potential using the Tafel extrapolation method. The results of the free corrosion potential were consistent with previous experimental data. By controlling the pH of the NaCl(aq), we determined the pH dependence of the free corrosion potential, and the results agreed well with the experimental results. Our results indicated that the ESM-RISM method duly described the environmental effect of an acidic solution and precisely determined the free corrosion potential. Therefore, the application of this method presents an efficient approach toward calculating the free corrosion potential for various reactions.


Introduction
Protecting materials from corrosion, which entails charge-transfer (CT) reactions at the material/electrolyte interface, plays an essential role in controlling their lifespan [1,2,3,4]. During corrosion reactions, materials react with their environment, and thereafter dissolve into an electrolyte solution. Therefore, owing to the corrosion reaction, the aging degradation of the material gradually progresses. For corrosion protection, the free corrosion potential, which is an equilibrium potential of the corrosion reaction, is a significant physical quantity because it is an indicator of the corrosion reaction from an electrochemical perspective. Consequently, the free corrosion potential of various materials has been extensively measured [5,6,7,8,9,10,11].
Theoretically, the corrosion reaction has been investigated through first-principles density functional theory (DFT) studies [12,13,14,15,16]. However, there are merely a few studies on the free corrosion potential from DFT calculations [15], as it strongly depends on the environmental effects at the material/electrolyte interface, such as temperature, pH , and halogen species of the electrolyte. Therefore, the development of a simulation technique that accurately describes the electrochemical environment at the material/electrolyte interface is essential.
First-principles electronic structure calculations combined with first-principles molecular dynamics (FPMD) is a powerful tool for studying the physical phenomena at the electrochemical interface [17,18,19,20,21]. However, the FPMD approach is expensive as it due to high computational costs because numerous sampling trajectories are needed to describe the average thermal properties of the solution. In the free corrosion potential, the local cell generates a potential difference between the anodic and cathodic reactions. The difference introduces excess electrons on the surface, and the ions in the electrolyte solution screen the electrostatic field generated by the excess electrons, resulting in the formation of an electric double layer (EDL). Therefore, the whole system is neutrally charged, but the origin of the compensating charge on the electrode and in the electrolyte differs. The former is the electron in an electrode, while the latter represents the ions in an electrolyte solution. To introduce these additional charges into the system, we need to use a simulation technique that enables the grand canonical treatment of both electrons and ions. In FPMD simulation, such a particle treatment approach is challenging. Thus, a more flexible simulation technique to describe electrochemical interfaces is required.
To overcome such FPMD limitations, a hybrid solvation model, which is the DFT method combined with the classical liquid or classical DFT, has been developed [22,23,24,25]. In the hybrid solvation model, we explicitly treat an electrode and a few molecules in the electrolyte. Conversely, the remaining electrolyte is replaced by the classical continuum model. The hybrid solvation model has been applied to the electrode/electrolyte interface and satisfactorily describes differential capacitances and electrode potentials [22,24,25].
In this study, the DFT method combined with the effective screening medium (ESM) technique [26,27] and the reference interaction site model (RISM) [28,29,30], termed as ESM-RISM [31,32], is applied to the Al/aqueous NaCl (NaCl(aq)) solution interface. As Al is widely used in our daily lives because of its high anticorrosion properties, it is one of the essential materials applied in the field of corrosion. Here, we focus on the corrosion reaction under acidic conditions. The objective of this study is to calculate the free corrosion potentials at the interface between the Al electrode and NaCl(aq) under acidic conditions using the ESM-RISM approach.
This paper is organized as follows: Firstly, the ESM-RISM method is briefly elucidated in Section 2, and we describe a series of corrosion reactions used in this study. Next, we explain the method for calculating the equilibrium potential from electrochemical and thermodynamic perspectives. Thereafter, the Tafel equation, using ESM-RISM to calculate the free corrosion potential, is discussed. Section 3 presents the results of the calculations. Here, we first present a series of results of the free corrosion potentials at the Al/NaCl(aq) interfaces electrochemically calculated using the ESM-RISM method. Next, the differences in the results between the electrochemical and thermodynamic methods are examined. We thereafter compare the pH dependence of the free corrosion potentials between the present theory and the previous experiment. Finally, the results of the free corrosion potentials depending on the surface orientation are revealed. In Section 4, we discuss the results of the present calculations from the perspectives of methodology, pH, and surface dependence of the corrosion potentials. The last section provides a summary of the present study.

Methods and models
In this section, we describe the theoretical and computational details used in this study. Firstly, we briefly describe the formulation of the ESM-RISM method and grand potential. Secondly, a set of reactions that is assumed to occur under corrosive conditions is described. Next, we explain the calculation methods for the electrode potential, current, and free corrosion potential for the anodic and cathodic reactions. Finally, the computational details are provided.

ESM-RISM method
Firstly, we briefly describe the ESM-RISM method [31,32]. The ESM-RISM is a hybrid solvation model, which is a self-consistent electronic structure calculation that considers the interactions between explicit particles and implicit solvent molecules (and ions) in a solution. The ESM-RISM can treat not only electrons in the explicit system [33] but also RISM components in the implicit solution under grand canonical conditions [31]. Consequently, the well-defined inner potential Φ S , which is the bulk potential of the solution, is used as the reference electrode potential [32].
In the ESM-RISM approach, the total energy of the whole system is described by the Helmholtz free energy A as where, E DFT and ∆µ solv are the total energy of the DFT and the excess chemical potential of the solution determined by the RISM, respectively. ∆µ solv corresponds to the excess free energy by considering the distribution of the implicit solvent at a given temperature, which is determined by solving RISM equations using a closure function [28,29,30]. We define the grand potential Ω using A by: where N e and N 0 e , respectively, denote the total number of electrons for the system with and without an excess charge. In the formulation of Ω, we consider only µ e and ∆N e because ∆µ solv includes a change in the chemical potential and the number of RISM particles by solving the RISM equations. We expect that the Ω profiles are the approximate inverse of the parabolic shape centered at the potential of zero charge (PZC), which is consistent with the constant capacitance model, as shown in previous studies [32,34].
According to a previous study [34], a polynomial model is suitable for describing the profile of Ω. We expand the Ω around the potential at Φ 0 to the third order, as follows: Here, the notation of the first term denotes the Ω(Φ 0 ), while the second, third, and fourth coefficients correspond to the contributions to the Ω from the excess surface charge, capacitance, and change in the capacitance by the change in the potential, respectively. The electrode potential Φ is obtained by using the electron chemical potential µ e as Φ = −µ e /e. We will use this polynomial-fitted Ω in the analysis of the potentials. Figure 1 shows the schematic of the corrosion reaction in acidic aqueous solutions. The corrosion reaction is considered as two CT reactions, i.e., anodic and cathodic reactions. For Al, the anodic reaction for the acidic aqueous solution is written as follows:

Corrosion reaction
Here, we consider that the anodic reaction is Al in the solid state (Al(s)) dissolving into the electrolyte as Al 3+ ion (Al 3+ (sol)) remaining electrons on the Al surface (e − (M)). When the electrolyte solution is acidic NaCl(aq), the Cl − ions participate in the anodic reaction as follows [4]: Al(s) + 4Cl − (sol) where we consider the complex formation resulting from the Al 3+ combined with the Cl − ion of NaCl(aq). The cathodic reaction entails the protons in the aqueous solution receiving electrons at the electrode surface, and then, the H 2 molecule in the gas phase (H 2 (g)) is produced, as shown in figure 1. Thus, the cathode reaction is defined as: According to a previous study [15], the equilibrium potential of the corrosive reactions is thermodynamically calculated by adding the surface, vacancy formation, and adsorption energies of Halogens to the free energy difference for the dissolution of the metal ions from the bulk. Thus, we should consider the corrosion reaction on the Al surface to precisely determine its equilibrium potential. Consequently, as the practical forms of a set of corrosion reactions (4)-(8), we consider that the dissolved Al 3+ ion (or complex) originates from the electrode surface, and the vacancy-type defect (V Al ) is formed on the Al surface by the anodic reaction. In addition, we normalize the number of transferred electrons to one for the corrosion reactions. Therefore, we rewrite the anodic and cathodic reactions as follows: 1/3Al m (s) → 1/3Al + (sol) 1/3Al m (s) + 4/3Cl − (sol) → 1/6Al2Cl 6 (sol) Here, m denotes the number of surface Al atoms. To determine the equilibrium and free corrosion potentials, we use equations (9)- (12) and (13) for the anodic and cathodic reactions, respectively.

Equilibrium potential
To evaluate the free corrosion potential, we determined the equilibrium potentials for the anodic and cathodic reactions using the ESM-RISM method. There are two pathways for calculating the equilibrium potential (E eq ): (i) E eq is determined by µ e at the crossing point of Ω for the left-and right-hand components of the reactions, which gives the equilibrium chemical potential of the electron (µ eq ); and (ii) E eq is calculated by the formulation of electromotive force (EMF) with the difference in Gibbs free energies between before and after the reactions (∆G). Methods (i) and (ii) are, respectively, written as follows: Here, superscripts "a" and "c" indicate physical quantities for anodic and cathodic reactions, respectively, throughout this paper. The subscript "L(R)" represents the left-(right-)hand component of the reaction. Furthermore, n, and F indicate the number of reacted electrons and the Faraday constant, respectively. As using Eq. (14) corresponds to the E eq determined by controlling the bias potential, we refer to method-(i) as the electrochemical method. In contrast to the electrochemical method, we denote method-(ii), using Eqs. (15) and (16), as the thermodynamic method. In the following two subsections, we describe details of the electrochemical and thermodynamic methods used in the present study.

Electrochemical method
Firstly, we describe the practical form of the electrochemical method using reactions (9) and (13) as those for the anode and cathode, respectively. In this method, we calculate the µ eq from µ e at the crossing point of the Ω profiles for the left-and right-hand components of the reactions. The grand potentials for the anodic and cathodic reactions (Ω a and Ω c ) are defined as follows: Subscripts "L" and "R" denote the left-and right-hand components of the anodic and cathodic reactions, respectively. E ZP is the zero-point vibration energy of isolated molecules or ions. Moreover, A[H 2 (g)] represents the sum of the DFT energies for the H 2 molecule in an isolated system and the standard molar for the H 2 gas. As the cathodic reaction indicates the CT reaction via the surface, we include the surface term for the formulation of Ω. The electrochemical method correctly describes the change in the EDL at the interface caused by the electrode potential. Therefore, we naturally describe the change in the interfacial capacitance caused by defect formation at the electrode surface. Nevertheless, because the interaction between the electrode surface (M) and solvated ions (or molecules) (I) is weak, Ω is separated into the M and I parts as follows: Notably, these separation forms are not applicable to desolvated ion problems, such as the charge/discharge process of the lithium-ion battery [36]. We rewrite the formulations for Ω a and Ω c as follows: Here, we use the definition of Ω described by Eq. (2) for deriving the separated formula. Practically, we use the final formula to determine the equilibrium potentials for the anodic and cathodic reactions.

Thermodynamic method
Next, we briefly describe the thermodynamic method. We also use reactions (9) and (13) at the anode and cathode, respectively. Using the ESM-RISM, we calculate the electrode potential E a(c) measured from Φ S with the EMF formulation. The ordinal formulation of EMF and eq. (15) differ from the electrochemical reaction between the neutral and charged states. Thus, the former corresponds to the full-cell reaction, while the latter indicates a half-cell reaction. Previous DFT studies [36,37,38,39] show that the EMF derived through thermodynamics is in good agreement with that of the experimental results. The equilibrium potential obtained by the electrochemical method agrees with that derived by the thermodynamic method, which has been shown by Haruyama et. al [32]. The definitions of the Gibbs free energy for the left-and right-hand components of the anodic and cathodic reactions are: We use the Helmholtz free energy (A) as the approximate Gibbs free energy (G) because A is the same as G by neglecting the volume×pressure term. When we can neglect the change in the details of the EDL, the thermodynamic method gives the same result of E eq as that of the electrochemical. In such a case, the thermodynamic method is flexible compared to the electrochemical approach because the grand potential calculations require many computational demands.

Free corrosion potential
Firstly, we discuss the calculation method for free corrosion potential, which is defined as the equilibrium potential between anodic and cathodic reactions. At this potential, the reaction rate for the anodic reaction is equal to that for the cathodic reaction. Phenomenologically, the Butler-Volmer equation describes the relationship between such currents and electrode potentials in the CT rate-determining step. However, the Butler-Volmer equation requires the simultaneous consideration of the anodic and cathodic reactions. To determine the free corrosion potential, the Tafel equation is useful, which is an approximated formula of the Butler-Volmer equation within the high voltage regime. The Tafel equation for the anodic and cathodic reactions with a single transferred electron is written as: where i a(c) is the current for the anodic (cathodic) reaction, and i 0 is referred to as the exchange current. α a(c) , R, and η are the transfer coefficient, gas constant, and over-potential (η = Φ−E eq ), respectively. Using a logarithmic plot of |i a(c) |, the free corrosion potential is determined by the potential at the crossing point of the anodic and cathodic currents. The log |i (c) |-voltage plot is called the Tafel plot. This approach to determine the free corrosion potential is referred to as the Tafel extrapolation method. Here, we assume that α a = α c = 0.5, as used in the experiments discussed in Ref. [15], and the current for the anodic (cathodic) reaction i a(c) linearly depends on the reaction rates. Based on these assumptions, the dependence of the Tafel equation is written as the difference in grand potentials between the left-and right-hand components of the reaction (∆Ω a(c) ) as follows: where subscripts L and R denote the left-and right-hand components of the reactions, respectively. β indicates the temperature factor. From the logarithmic plot of |i a(c) | as a function of µ e , we obtain the Tafel plot without the exchange current for the anodic and cathodic reactions. As log |i a(c) | directly depends on ∆Ω a(c) , we obtain the free corrosion potential from µ e at the crossing point of ∆Ω a and ∆Ω c . Thus, the definition of the free corrosion potential assumed in this study is the equilibrium potential at the same number of transferred electrons between the anodic and cathodic reactions. Next, we derive a more flexible formula for ∆Ω a(c) , which directly corresponds to the logarithmic plot of the anodic (cathodic) current, using the formulations of the thermodynamic method. For simplicity, we use the anodic reaction (9) as a representative of the reformulation. Substituting Eqs. (22) and (23) to Eq. (32), we obtain: Here, we assume that the curvature of the ground potential (e.g., the interfacial capacitance and its potential change) and the PZC, which is the apex of the ground potential, remain unchanged before and after the reaction. Furthermore, ∆Ω a shows the difference in the Helmholtz free energy before and after the reaction. Therefore, Eq. (33) becomes: ∆Ω a = 1/3A(Al m−1 (s)) Here, we substitute Eqs. (26) and (27) in the second equality, and the EMF formulation [Eqs. (15) and (16)] is used in the last equality. From the same procedure, the ∆Ω c is given as: Using Eqs. (34) and (35), we analytically obtain the chemical potential at the crossing point of ∆Ω a and ∆Ω c . Thus, the Tafel extrapolation is analytically conducted as follows: Here, µ corr is the electron chemical potential for the free corrosion potential. Through this relation, the free corrosion potential E corr is obtained as −µ corr /e. Eq. (36) is also flexible, similar to the thermodynamic method, because we obtain the free corrosion potential without grand potential calculations.
For electronic structure calculations, the cut-off energies for wave functions and augmented charges were 40Ry and 320Ry, respectively. The spin-unpolarized version of the generalized gradient approximation parameterized by Perdew, Burke, and Ernzerhof [52] was used for the electron-electron exchange-correlation function. Furthermore, k-point samplings of 16×16×16, 4×4×1, and only the Γ-point was used to calculate the Al bulk, Al slabs, and ions (or molecules), respectively. The occupation number of electrons was determined by the Gaussian smearing method with a smearing width of 0.01Ry. Four-layered slab models were used to represent Al(100) and Al(111) surfaces. We carried out cell optimization for the Al bulk and obtained an optimized cell parameter of 4.038Å, which was in good agreement with the experimental lattice constant of 4.05Å [53] and used to construct the surface slabs. We will further discuss the details of slab geometries and calculation cells later. For ions and molecules, we carried out structural optimization until all forces acting on the atoms became less than 1.0×10 −3 Ry/Bohr. To calculate the potentials, the E ZP for bonded molecules and ions are determined using the Gaussian16 program package [54] with BLYP/6-311+G(3df, 2pd); the results are listed in Table. 2.
To evaluate the solute, the NaCl(aq) solution at a temperature of 300 K was represented by the RISM [28,29,30]. In the present study, 0.5 M or 1 M NaCl(aq) solutions under the conditions of pH = 1, 2, and 7 were adopted. As pH = 7 represents the neutral condition of the solution, it is out of the range of the present study. However, as a reference, we carried out calculations for neutral solutions. The concentrations of H 2 O solvent and NaCl salt were 1.0 cm 3 /g (55.6 mol/L) and 1mol/L, respectively. For calculations at pH = 1 (or 2), 0.1 (or 0.01)mol/L HCl was added to the NaCl(aq). The calculation of the 1 mol/L NaCl(aq) solution under the condition of pH = 7 [NaCl(aq, pH=7)] was carried out without HCl. To solve the RISM equation, we used the Kovalenko and Hirata type of closure function [30], and the convergence criterion of residuals of correlation functions in the RISM equation was 1×10 −6 . The cut-off energy for the reciprocal representation of the Laue-represented RISM equations was set to 160Ry. For the calculation of ∆µ solv , the Gaussian fluctuation method was used [55]. In the present study, we used the LJ-type classical force field to describe site-site interactions. Table 1 shows the LJ parameters and charges used for the RISM calculations. The Lorentz-Berthelot combination rule [56] was used for the LJ parameters for different types of atoms.
For the ESM-RISM calculations, we used two boundary conditions: vacuum/slab/solvent (VSS) and solvent/ion/solvent (SIS) . Figures 2(a) and (b) show the schematics of unit cells for the ESM-RISM calculation with each boundary condition. For surface calculations, Al(100) and Al(111) slabs had surface areas of 12.114×12.114Å 2 and 8.565×9.891Å 2 , respectively, with a length of 40Å in the z-direction for the DFT Table 1: 12-6 type of Lennard-Jones (LJ) parameters (ε and σ) and charges (Q) for each site of RISM used in the present study are listed. Here, the solvent and solute denote components of the RISM and explicit particles, respectively. An implicit H 2 O molecule is represented by the TIP5P model [40] for precisely determining the hydration free energy for a H 3 O + ion [32]. The LJ parameter of Al 3+ was used to calculate the Al 3+ (sol), AlCl 3 (sol), AlCl − 4 (sol), and Al 2 Cl 6 (sol).  [46] 0.505 4.01 - cell. We carried out calculations of ions and molecules using a cubic cell with dimensions of 20×20×20Å 3 . In the RISM calculations, to obtain the converged distribution function of the solution, an expanded cell with a length of 31.75Å in the z-direction was attached to the right-hand-and both sides of the calculation cells for the VSS and SIS, respectively. For VSS, we take a 10Å of the vacuum region from the left-hand boundary of the calculation cell, whereby the ESM boundary condition of the vacuum is attached to the left-hand side of the unit cell.

Results
In this section, we present the results obtained through the ESM-RISM method. Firstly, we describe the Ω profiles at the Al(100) and Al(111)/NaCl(aq, pH=1) interfaces with and without V Al . Secondly, we discuss the electrode potential results for the anodic and cathodic reactions, and then present the results of the Tafel plots for the anodic and cathodic reactions. Thirdly, we compare results of free corrosion potential between the electrochemical and thermodynamic methods. Thereafter, we appraise the pH-dependent behavior of the free corrosion potential relative to the previous experiment. Finally, the surface dependence of the free corrosion potentials and defect formation energies are presented. The results of the grand potential profile evaluation as a function of µ e [vs. Φ S ] for the Al(100) and Al(111)/1 M NaCl(aq, pH = 1) interfaces with and without V Al are shown in Fig. 3. The Ω profiles are obtained by introducing excess charge (∆N e ) to the electrode surface, where ∆N e changes from -1.0e/cell to +1.0e/cell at 0.2e/cell intervals. Here, we neglect the specific adsorption of halogen ions (Cl − ) on the surfaces. To compare Ω profiles between the surfaces, we select the maximum values of Ω (Ω PZC ) as the zero energy for each Ω. As expected, all the Ω profile results show the approximate inverse of the parabolashape centered at the point ∆N e = 0, which is consistent with the constant capacitance model. The µ e [vs.Φ S ] at ∆N e = 0 for clean Al(100) and Al(111)/1M NaCl(aq, pH=1) interfaces are −4.44eV and −4.10eV, respectively. Thus, the PZC for electrodes depends on the surface orientation, which is consistent with the findings of a previous experiment [57]. For defected surfaces, the curvatures and PZC of Ω are similar to those for clean surfaces. This result indicates that the change in the details of the EDL by introducing V Al is relatively small. Thus, we expect that the difference in the results of equilibrium potentials between Eqs. (14) and (15) is rather small.

Tafel plot from the electrochemical method
Next, we fit the results of Ω by a third-order polynomial function, which is shown as solid lines in Fig. 3. The fitted curves duly described the original data obtained by the ESM-RISM calculations. We calculated Ω for all combinations of interfaces between Al surfaces and solutions and performed polynomial fitting for the results of Ω. To simplify further analysis, hereafter, we will use the Ω profiles obtained by the polynomial fitting. Figure 4 shows the results of the Ω profiles at the Al(100)/1 M NaCl(aq, pH = 1) interface for leftand right-hand components of anodic and cathodic reactions. Here, we only show the results of Ω using reaction (9) as a representative of the anodic reactions. The equilibrium potentials for the anodic (cathodic) reaction were obtained by the chemical potentials at the crossing points of Ω a(c) L and Ω a(c) R (µ a(c) eq ). From the Ω profiles, the anodic reaction is reduced to a lower µ e than µ a eq because in such a reaction, Ω R becomes lower than Ω L . Conversely, the cathodic reaction occurs at a higher µ e than µ a(c) eq . Thus, the corrosion reaction occurs in the range of µ e from µ c eq to µ a eq . The result of µ c eq corresponds to the µ e of the reversible hydrogen electrode (RHE) potential at pH = 1, when it is measured from the inner potential. The electrode potential measured from the inner potential (Φ [vs. Φ S ]) was converted to that measured from RHE (E [vs. RHE]) by the RHE potential relative to the inner potential (E RHE [vs. Φ S ]) as follows: Here, µ RHE [vs.Φ S ] denotes µ eq [vs. Φ S ] corresponding to the RHE. Notably, µ RHE [vs. Φ S ] is equal to the equilibrium potential of the cathodic reaction with respect to the inner potential of various solutions. In addition, the RHE potential (E RHE ) measured from the inner potential at temperature T can be represented using the standard hydrogen electrode (SHE) potential E SHE as follows: Here, k B is the Boltzmann constant, and the coefficient of the second term in Eq. (38) is 0.059 V at a temperature of 298 K. In the ESM-RISM calculation, we can determine only E RHE [vs. Φ S ] as a relative potential because we usually choose the inner potential as the reference potential level. Consequently, we cannot directly compare the electrode potentials using different solutions without the contact (or Volta) potential difference [32]. However, we can convert the RHE potential to the SHE potential using eq. (38). This conversion allows us to compare electrode potentials with different solutions obtained by the ESM-RISM, using the unique reference potential. Figure 5 shows the results of the logarithmic plot of i a(c) at the Al(100)/1 M NaCl(aq, pH = 1) interface as a function of voltage vs. RHE obtained by the ESM-RISM. Here, we neglect the constant values i.e., β and α a(c) , in Eq. (1) for simplicity. The voltages at the crossing points of ∆Ω for the cathodic and anodic reactions denote the free corrosion potentials for each anodic reaction, which are measured from the RHE potential. In addition, the voltages at ∆Ω = 0 correspond to the equilibrium potentials measured from the RHE potential for each reaction. Thus, the anodic and cathodic reactions do not occur for voltages in the ∆Ω > 0 range. From this plot for various pH conditions of NaCl(aq), we can obtain the pH dependence of the free corrosion potentials measured from the RHE potential, and thereafter calibrate the reference electrode potential from the RHE to the SHE using Eq. (38).  Here, we compare the results of the free corrosion potentials by the electrochemical and thermodynamic methods, as shown in Table 3 with the products of the anodic reactions. When the change in the interfacial capacitance is negligibly small, the free corrosion potentials obtained by both methods have the same value. We found that the results of free corrosion potentials by the electrochemical method (E corr (Ω)) were in reasonable agreement with those obtained using the thermodynamic method (E corr (G)). The differences between E corr (Ω) and E corr (G) originate from the deviation of ∆Ω from the linear function. This difference in ∆Ω is regarded as the change in the interfacial capacitance due to the formation of V Al . Thus, the change in the details of EDL slightly affects the free corrosion potential (the absolute value of the difference in the free corrosion potential between the methods is an order of error within 0.01 V). Figure 6 presents the results of free corrosion potential vs. SHE for various conditions of the NaCl(aq) solution, as well as the experimental data for the 1 M, 0.5 M, and 0.1 M NaCl(aq) solutions, which are, respectively, obtained from Refs. 6, 5 and 8. The free corrosion potentials for Al 3+ , AlCl 3 , and Al 2 Cl 6 are in reasonable agreement with the experimental results. Conversely, the results of E corr for AlCl − 4 are lower than those for other products. This result suggests that the main corrosion products of the anodic reaction due to H 2 creation are Al 3+ , AlCl 3 , and Al 2 Cl 6 . The calculated E corr decreases with increasing pH. This tendency is also in reasonable agreement with the experimental data under acidic conditions. Therefore, the present model of the corrosive reaction duly reproduces the experimental data at the Al electrode in the acidic NaCl(aq) solution. Although the results of E corr for the Al/NaCl(aq, pH = 7) interface are reference data, we compared their results to those of the experiment. Compared to the acidic region, the calculated E corr for the neutral condition of NaCl(aq) deviates from the experimental data. Thus, the present model is limited in the acidic condition of the NaCl(aq) solution.  Here, we show the results of E corr depending on the surface orientation, which are presented in Fig. 7. We found that the results of E corr for Al(100) are approximately the same as those for Al(111). The main

Interface
Al(100) Al(111) Vacuum 0.54eV 0.58eV NaCl(aq, pH=1) 0.60eV 0.65eV NaCl(aq, pH=2) 0.61eV 0.65eV NaCl(aq, pH=7) 0.61eV 0.65eV contribution of the surface dependence of E corr is attributed to the stability of the Al atoms on the outermost surfaces. Here, we show the results of the vacancy formation energy for the outermost Al atom (E(V Al )), calculated as: The results of E(V Al ) are listed in table 4. The results of E(V Al ) for the vacuum/Al interfaces agree with those of previous DFT calculations [58,59,60]. For all interfaces, E(V Al ) at Al(100) is slightly smaller than that at Al(111). Thus, the stability of the Al atoms, depending on the surface orientation, is weak. This explains the small difference in the results of E corr between Al(100) and Al(111).

Discussions
In this section, we first discuss the comparison between the electrochemical and thermodynamic methods. Secondly, the pH dependence of the free corrosion potential and corrosive reactions is discussed. Finally, we elucidate the surface-dependent corrosion.
If the changes in capacitance and the PZC potential before and after the reaction are negligible, the thermodynamic method is an effective way to reduce the computational cost. This is because the calculation of the SIS is much cheaper than that of the VSS. Therefore, this method is useful for a systematic study of E corr with varying parameters, such as pH and temperature. However, it is necessary to simulate the change in the EDL, which plays a central role, using electrochemical methods, such as studies on the adsorption of ions (or molecules) [34], and the desorption/insertion processes of ions [36].
Based on the agreement between the present theoretical and experimental data, we emphasize that the ESM-RISM method satisfactorily reproduces the nature of the solid/liquid interface under acidic conditions. In contrast, the difference in the results for the neutral NaCl(aq) solution between the theoretical and experimental data originates from the formation of a passive layer. According to previous studies [4,61], the Al surface is covered with a passive layer of oxide or hydroxide under the conditions of pH>4. The passive layers ensue from the dissolved O 2 participating in the cathode reaction; thus, the passive layer alters the details of the corrosion reactions and the free corrosion potential. Therefore, the model of corrosion used in the present study does not reproduce the free corrosion potential for neutral NaCl(aq) solutions. To determine the whole behavior of the pH dependence of E corr , we should systematically consider other anodic and cathodic reactions. We must investigate such systematic considerations on this corrosion model in future studies.
Finally, we briefly discuss the surface dependence of the corrosion potential. The present results of free corrosion potentials weakly depend on the surface orientations. However, a previous experiment for single-crystal Al shows that the order of the difference in the pitting corrosion potential between the surface orientations is approximately 10 mV [62], which is higher than the present results of E corr . This difference between the present and previous studies is mainly due to the absence of Cl − adsorption in the present study. The previous experiment suggests that pitting corrosion occurs via halogen adsorption onto the Al surface [62]. Thus, to further understand the surface dependence of the corrosion reaction, we need to consider the halogen adsorption on individual Al surfaces.

Summary
We investigated the free corrosion potential at the Al/NaCl(aq) interface using the ESM-RISM method. Firstly, we considered a set of corrosion reactions and formulated a method to determine the free corrosion potentials using the ESM-RISM. The equilibrium potentials for the anodic and cathodic reactions were determined from the results of the grand potentials. Thereafter, we obtained the Tafel plot by evaluating the difference in the grand potentials between the right-and left-hand components of the reactions. The results of the free corrosion potentials between the electrochemical and thermodynamic methods was compared. The results of the free corrosion potentials between these two methods agreed with each other within the order of 0.01 V. Therefore, we concluded that the change in the EDL before and after the reaction was inconsiderable. In such a case, we submit that the thermodynamic method is more suitable for obtaining the free corrosion potential owing to a lower calculation cost.
The pH dependence of the obtained free corrosion potentials was compared to the experimental results with the reference potential of the SHE. The pH dependence of the free corrosion potentials under acidic conditions is in good agreement with the experimental results. From this result, we emphasize that the ESM-RISM correctly describes the environmental effects of NaCl(aq) solutions. Nonetheless, the results for neutral NaCl(aq) solutions deviated from the obtained experimental data. This phenomenon suggests that we need to consider the formation of the passive layer on the electrode in the case of neutral NaCl(aq) solutions.
Finally, we examined the dependence of the free corrosion potential on the surface orientations. The difference in the free corrosion potentials between Al(100) and Al(111) is small. This result is consistent with the results of the vacancy formation energies. From these results, we conclude that the areal dependence of the free corrosion potential at the Al/NaCl(aq) interface is small. However, by considering the pitting corrosion, the corrosion potential depends on the surface orientations, which may be caused by Cl − adsorption. To better understand the concept of corrosion reactions, we additionally considered halogen adsorption on the surface.