**2. Computational Details**

In this study, three kinds of Fe models—single crystal, sigma3 or sigma5 twins, and polycrystalline—were used to study the influence of grain boundary on the corrosion process of pure Fe-H2O surface systems. The sigma3 twin boundary is composed of the (111) plane, and the sigma5 twins plane is composed of the (02-1) plane. The structure at the twin grain boundary is stable and the interfacial energy is low. However, the polycrystalline grain boundary structure is unstable, the lattice distortion is large, and there are many defects. The stress in this paper is the shear stress of the Fe substrate on the Y-axis. It is obtained by counting the stress sum of Fe atoms at the grain boundary in the Y direction and dividing it by the corresponding total atomic volume. The volume of polycrystalline grains is 35 nm3. The schematic diagram of the models is shown in Figure 1. The Fe substrate comprises twelve atom layers in the polycrystalline simulation system, including 9215 Fe atoms and existing many defects such as grain boundaries and stacking fault. The water-side part includes the water layers of 40 Å thickness with the density of 0.99 g/cm3 (about 4669 water molecules at 300 K and 1 atm).

**Figure 1.** Schematic diagram of the simulation model for the Fe-water surface system.

All the simulations were based on the open-source Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) code [23] and the atomic interactions were described using the ReaxFF for pure Fe-H2O systems. The original ReaxFF was developed by van Duin et al. [24]. The ReaxFF-MD can simulate the chemical reaction process that includes chemical bonds breaking and forming between atoms; therefore, it is suitable for studying interfacial corrosion. Recently, our group developed a new Fe-H2O ReaxFF potential function [25], which can accurately describe the nature of the defects in Fe substrate and the Fe-water interactions. The interaction between vacancy and hydrogen/oxygen is emphasized to describe the precipitation of hydrogen/oxygen near the vacancy. The energy contributions to the Fe-H2O ReaxFF potential are summarized by the following [25]:

$$E\_{system} = E\_{bond} + E\_{over} + E\_{under} + E\_{val} + E\_{lp} + E\_{Coul} + E\_{H-bond} + E\_{vdW} \tag{1}$$

which includes terms related to bond, angle, over coordination, undercoordination, lone pair, van der Waals, Coulomb, and H-bond energies. The details of the ReaxFF potential have been shown in our recent paper. The essential part of ReaxFF is that the charge equilibration (QEq) method can be used to obtain the charge distribution [26,27]. The charge values were determined at each simulation time step and depended on the system's geometry. This feature made it possible to describe charge transfer in chemical reactions using ReaxFF. Subsequently, this study did not consider an absolute stress value analysis, instead focusing on relations between reactivity and stress states.

The periodic boundary conditions were applied in the *X*- and *Y*-directions in all simulated boxes, and the *Z*-direction uses aperiodic boundary conditions. A reflective wall is added at the upper Z boundary to prevent atoms from passing through the boundary. The Fe substrate and water molecules are relaxed with a Nose/Hoover isothermal–isobaric (NPT) [28] to reach equilibrium. Nose–Hoover thermostat [29,30] is used to maintain the prescribed system temperature during corrosion for the canonical (NVT) ensemble. The MD simulation time step is 0.2 fs. The corrosion process of the pure Fe-H2O system was observed by external electric field [31] acceleration at 300 K and NVT ensemble, which is located at 0~3 Å on the surface of the Fe substrate. We chose the most suitable external electric field to be 325 MV/cm by observing the corrosion phenomenon. All subsequent simulations are performed under this electric field strength. The relationship between

individual atomic charges and the corresponding electrostatic energy *E*(**q**) with time is shown in Equation (2):

$$E(\mathbf{q}) = \sum\_{i=1}^{N} \left[ \chi\_i q\_i + \eta\_i q\_i^2 + \operatorname{Tap} \left( r\_{ij} \right) k\_\varepsilon \frac{q\_i q\_j}{\left( r\_{ij}^3 + \gamma\_{ij}^{-3} \right)} \right] \tag{2}$$

where **q** represents a vector of length *N* containing the charges, *qi* is the charge of ion *i*, *N* is the total number of ions, *kc* is the dielectric constant, *χ<sup>i</sup>* and *η<sup>i</sup>* are represent the electronegativity and the hardness of ion *i*, *Tap*(*r*) is a seventh-order taper function, and *γij* refers to the shielding parameter between the two atoms *I* and *j*.

We reproduced their works to verify the correctness of the ReaxFF potential according to the article of DorMohammadi et al. [19]. The simulation models mainly consist of body-centered cubic (bcc) Fe and water molecules. The corrosion processes of low-index surfaces (100), (110), and (111) were accelerated by applying an external electrical field. The corrosion phenomenon is consistent compared with the works of DorMohammadi et al. The difference is that the diffusion of H atoms into the substrate can be clearly observed during the corrosion process. Therefore, we believe that the ReaxFF potential parameters can simulate Fe-water interfacial corrosion accelerated by an external electric field at room temperature. The relevant verification process is in the supporting materials.
