Effect of Iron Isomorphic Substitution in Mg:Al and Zn:Al-Layered Double-Hydroxide Structures by Means of First Principle Calculations

Layered double hydroxides (LDHs) are important components in terrestrial and extra-terrestrial environments. The presence of iron in these minerals provides them a wide potential application in environmental and materials sciences. In this work, the role of Fe in the crystallographic properties of LDHs M2+:M3+ 2:1 with Mg:(Fe,Al), Mg:Fe, Zn:(Fe,Al), and Zn:Fe is investigated by means of quantum mechanical calculations based on the density functional theory (DFT). Several relative proportions of Fe are studied. The cation ordering of these LDHs has been explored, finding useful insights for experimental synthetic paths of these minerals. The a and b cell parameters increase with the iron concentration. Some diffraction lines at high angle decrease in angle and increase in intensity with the increasing iron concentration. All of them agree with the experimental results. The iron substitutions tend to aggregate.


INTRODUCTION
Layered double hydroxides (LDHs) are metallic hydroxides with the crystal structure of brucite (Mg hydroxide). These hydroxide layers can have isomorphous cation substitutions in the octahedral sheet that in some cases can provoke charge defects, which are compensated with anions intercalated in the interlayer space. This fact gives the main properties of these minerals, which gives their high absorption capacity, variable water content, and anion exchanges. These capabilities and their small particle size make these minerals inorganic membranes with catalytic properties. Moreover, LDHs are present in seabed hydrothermal vents. In this scenario, LDHs provide confined spaces where the adsorption and their catalytic properties can become a possible scenario for the first prebiotic chemical reactions for the origin of life 1−3 in early Earth. 4 Furthermore, Fe hydroxides can also form LDHs with Mg, producing minerals such as pyroaurite or iowaite. 5 These iron minerals can be present in the serpentinization processes of early Earth. The Mg−Fe layers of LDHs are more reactive than no-Fe-bearing LDHs due to the increasing basicity of the hydroxyl groups. In addition, the presence of iron in LDHs also provides redox reactions due to the Fe 3+ −Fe 2+ binary system.
These Fe-bearing LDHs have also great medical and scientific interest. The iron content in LDHs increases the biocompatibility of these solids, helping in tissue regeneration by modulation of collagen production. 6 On the other hand, some Fe-bearing LDH minerals have the Fe 2+ −Fe 3+ pair in the structure, which are known as green rust, being the main mineral the so-called fougerite. 3,5,7 This mineral has been particularly related to the origin of life; therefore, many efforts have been reported for understanding its structure and physicochemical properties. 5,7−9 In addition, the cation ordering in LDHs is becoming an interesting issue, although different researchers have obtained contradictory results. Thus, some studies show that the cations in LDHs are distributed orderly, 10,11 whereas others show that these cations are randomly distributed. 12−14 In all cases, it was found that M 3+ cations avoid close contact. 13,15 Rozov et al. 16 found crystallographic changes with the presence and content of Fe in MgFe-LDH without considering the Fe distribution along the crystal lattice, which has been subsequently demonstrated by Figueiredo et al. 17 The Fe 3+ cation distribution in the octahedral sheet of minerals can affect some properties, such as the nuclear magnetic resonance (NMR) and infrared spectroscopy properties in clay minerals. 18−20 Taking into account the importance of Fe-bearing LDHs and cation ordering in LDHs, one of the aims of this work is to understand, at an atomic scale, the Fe cation distribution in LDHs and the variation of crystallographic properties of Mg:AlFe and Zn:AlFe LDHs with different relative proportions of Fe by means of first principle calculations.

MODELS
The LDH models used in this work were taken from our previous optimized theoretical models. 15 In that paper, the structure of hydrotalcite, 21 Mg 6 Al 2 (CO 3 )(OH) 16 ·4H 2 O, was modified to design an LDH structure with a Mg 2+ :Al 3+ ratio of 2:1, instead of 3:1, chlorine anions in the interlayer space instead of CO 3 2− and removing the water molecules. Moreover, the Zn-LDH structure was also built by replacing Mg 2+ cations with Zn 2+ cations.
From the energetically most favorable optimized structure, a LDH 4 × 3 × 1 supercell was constructed for both Mg:Al-LDH and Zn:Al-LDH. 15 The effect of Fe 3+ cations in the LDH structure was studied by substituting different amounts of aluminum cations in the supercell Mg 24 (Al,Fe) 12 . These substitutions were always done in pairs up to a maximum of 12Fe 3+ cations per supercell, which indicates the maximum amount of substitutions of trivalent cations in the supercell with a general formula Mg 24 Al 12−x Fe x LDH, named in short as Mg:Al 12−x Fe x . The same was done for Zn:Al-LDH with Zn 24 Al 12−x Fe x LDH models, named Zn:Al 12−x Fe x . All our models are dry, without any water molecules in the structure, to optimize the computational facilities and consider that this aspect does not alter the relative energy calculations between structures.

METHODS
The crystal structure and energies of the different LDH models were calculated using density functional theory (DFT). All these electronic calculations were carried out with the CASTEP code, implemented in the Materials Studio package, 22 with the generalized gradient approximation (GGA) functional with the exchange and correlation functional of Perdew−Burke−Ernzerhof (PBE). 23 All optimization calculations were performed using On-The-Fly-Generated (OTFG) ultrasoft pseudopotentials including the Koelling−Harmon relativistic treatment, 24 Grimme dispersion correction, 25 and an energy cutoff of 571.4 eV. The spin polarization was applied in two different ways, optimizing the spin state, or maintaining the minimal spin state. The Hubbard U correction 26,27 was also taken into account with U = 2.5 eV in the d orbitals of Fe of the Mg 24 Al 12−x Fe x LDH and Zn 24 Al 12−x Fe x LDH crystals (x ≤ 12). The results, with and without these corrections, were compared in fully optimized crystal structures (i.e., atomic positions and cell parameters). The convergence tolerance parameters were set to 5 × 10 −6 eV/atom for energy, 0.0001 Å for maximum displacement, 0.02 GPa for maximum stress, and 0.01 eV/Å for the maximum force. The convergence gradient for the self-consistent field (SCF) was 1 × 10 −7 eV/atom in the density matrix. All of the crystal models shown in this work have been plotted using the Vesta software. 28 Powder X-ray diffraction (XRD) patterns of the calculated structures were simulated with the REFLEX code included in Materials Studio. 22 The simulation was carried out considering a copper X-ray source. Moreover, the samples were considered to have a crystallite size of 35 nm, which produces diffractograms more similar to those obtained experimentally in previous works.

RESULTS AND DISCUSSION
As stated above in the model section, the Mg:Al and Zn:Al-LDH crystal structures used in this work were previously fully optimized by Pimentel et al. 15 using the same theoretical procedure as described above. However, other parameters have to be taken into account in systems with Fe 3+ , i.e., spin polarization and Hubbard correction. Therefore, the best method of calculation had to be chosen prior to the in-depth study of the LDH structures with Fe 3+ . Then, we used the Mg-LDH with 12Fe 3+ cations per supercell (Mg 24 Fe 12 ) as a model structure (Mg:Fe 12 ) ( Table 1). Other approaches not shown in the table were used to optimize the LDH structure, such as no spin polarization and using free spin polarization without the Hubbard correction. However, after several trials, no convergence was obtained. Therefore, the best approach for calculating the LDH with Fe 3+ is to use Hubbard correction with the optimized spin state (SP+U) and with the minimal net spin state equal to zero (SP0+U), whose optimized cell parameters are in good agreement with the previously reported experimental results 17 (Table 1). The SP+U method yields the closest cell parameters to the experimental values. In addition, the crystal structure optimized with SP+U has 0.078 eV/uc (i.e., eV per unit cell) lower energy than SP0+U, indicating that the high spin state is slightly more stable than the paramagnetic low spin state. Although both methods have been used in this work, the main method used is SP+U, with which a larger amount of possible LDH structures with Fe 3+ cations have been explored (including all of the Mg 24 Al 12−x Fe x LDH and Zn 24 Al 12−x Fe x LDH structures, Table  1). These results contrast with those previously obtained for Mg 24 Al 12−x Fe x LDH bearing carbonate anions, in which the use of Hubbard correction was found to yield a slight overestimation of the lattice parameters. 29 The largest difference between the theoretical and experimental results is found in the c parameter, which could be due to the different amounts of water in the structure (our theoretical structure being anhydrous). This difference has also been observed in previous works. 30 −32 Once the best calculation method has been determined, the Mg:Al 10 Fe 2 LDH structures with different arrangements of two Fe 3+ were studied ( Figure 1): placed in different layers (Mg:Al 10 Fe 2 DL), placed separated in the same layer (Mg:Al 10 Fe 2 SL), or placed together (Mg:Al 10 Fe 2 SLTa). After a full optimization of the geometry (i.e., lattice parameters and atomic positions), the geometry and the energy between the different structures were compared ( Figure 2 and Table 2). No significant energy difference was observed between these structures. Considering the low energy differences per unit cell between these three models, it can be assumed that all of them could be equally found in experimental LDHs, although the distribution of iron cations in the LDH unit cell has not been studied experimentally so far.
The next step to study the Mg-LDH structures with two Fe 3+ cations was to determine the effect of the orientation of the Fe 3+ cation pairs in the Mg:Al 10 Table 2). The difference in energies between the three arrangements is higher than those between DL and SL. Although the energy differences between these arrangements are not significant enough, the lowest-energy structure (Mg:Al 10 Fe 2 SLTa) is chosen as the default cation arrangement to build and calculate the structures with a higher amount of Fe 3+ cations.
In the Mg:Al 10 Fe 2 SLTa model, two Al 3+ cations were replaced with additional two Fe 3+ cations, being together in the next layer with respect to the previous Fe 3+ cation pair (Mg:Al 8 Fe 4 DL) or being all Fe 3+ cations together in the same layer (Mg:Al 8 Fe 4 SL). In the same way, two additional Fe 3+ cations were substituted, placing them together but in another layer (Mg:Al 6 Fe 6 ), where two Fe 3+ cations formed pairs in each layer in the 2−2−2 sequence along the c-axis. From this structure, two additional Al 3+ cations of the central layer were substituted by two Fe 3+ cations (Mg:Al 4 Fe 8 ), with the 2−4−2 sequence of the Fe 3+ cations in the layers of the supercell along the c-axis. Finally, all Al 3+ cations were substituted by Fe 3+ cations in the structure described above (Mg:Fe 12 ). All these LDH crystal structures were fully optimized, i.e., atomic positions and lattice parameters (Figure 3 and Table 3). The model with 4Fe 3+ cations together in the same layer (Mg:Al 8 Fe 4 SL) is 0.061 eV/uc less stable than that with the Fe 3+ cation pairs in different layers (Mg:Al 8 Fe 4 DL).
Hypothetical isodesmic reactions can be designed to explore some significant substitutions. The isodesmic reactions have the advantage of having the same structure, same bonds, and same number and type of atoms either in reactants or products. All errors due to the electron basis set, correlation exchange functionals, pseudopotentials, and different approximations for

ACS Earth and Space Chemistry
http://pubs.acs.org/journal/aesccq Article calculating the electronic structure and optimal atomic distances and the cell parameters would be canceled. This approach has been applied previously in octahedral cation substitutions in phyllosilicates 33 and organic reactions. 34 Three reactions have been designed by simplifying the names of components (Fe 12    layer, leaving two layers up and down without any substitution in our model. As a consequence, certain distortion is created between the layers, introducing an additional destabilization in the structure. A similar tendency is found in the Zn-LDH structures bearing 2 (Zn:Al 10 Fe 2 ), 6 (Zn:Al 6 Fe 6 ), and 12 (Zn:Fe 12 ) Fe 3+ cations per supercell (Figure 4). In such structures, the energy difference between Zn:Fe 12 and Zn:Al 6 Fe 6 was only 0.0629 eV/uc lower than that between Mg:Fe 12 and Mg:Al 6 Fe 6 . In the same way, the energy difference between Zn:Fe 12 and Zn:Al 10 Fe 2 is slightly higher, 0.0014 eV/uc, than that between Mg:Fe 12 and Mg:Al 10 Fe 2 .
The size of the lattice parameters increases with the increasing amount of Fe 3+ cations in the Zn-LDH structures (Tables 3 and 4), following linear relationships ( Figure 5), similar to Mg-LDH. This behavior of the LDH structures has also been observed in Mg-and Zn-LDH crystals synthesized experimentally in the laboratory. 16,17 By comparing both theoretical and experimental results ( Figure 5), we find that the computational method used here not only provides crystal structures with cell parameters in good agreement with those previously reported but also reproduces a slight change in the a/b cell axes measured in synthesis experiments with different iron contents. This behavior was also observed experimentally   Figure 6). Small differences can be observed due to the nature of the cation joined to them. In general, the ZnZnZn OH groups are longer than the MgMgMg ones. In the X:Fe 12 (X = Mg or Zn) models, the d(O−H) bond length follows the sequence XXX < XFeX < FeXFe, although the differences are smaller in the models with Zn than those with Mg. However, in the X:Al 10 Fe 2 SLTa models, the differences are diffused and not so sharp, indicating that other factors are involved, such as nonbonding interactions with the Cl − anions, and orientation of the OH bonds in the interlayer space. A similar phenomenon was observed previously in the octahedral sheet of phyllosilicates, which produces a shift to lower infrared frequencies in the vibration modes of these OH groups when the Fe 3+ are joined to these OH groups. 35,36 For each calculated structure, the diffractograms were also simulated ( Figure 7). The diffractograms are similar to those   experimentally obtained in the previous works on the LDH structures. 16,17,32 Taking into account that no symmetry was imposed in these calculations, some small peaks appear in the calculated diffractograms, in addition to the main LDH peaks that cannot be assigned to any LDH plane, which are spurious planes produced for considering a small number of unit cells. In our theoretical diffractograms (Figure 7), the positions of the (003) and (006) peaks are quite close for the different amounts of Fe cations considered; for example, in the case of the (003) peak, the positions are around 11.9°for both MgFe-LDH and ZnFeLDH structures. The main differences in the diffractograms are those related to the (110) and (113) peaks. In these cases, structures with higher iron content have lower 2θ angles (Figure 7). In MgFe-LDH structures, the (110) peak is placed at 60.5°for Mg-LDH (without any iron), and its position moves to lower spaces up to 59.4°for Mg:Fe 12 , while the position of the (113) peak changes from 62 to 60.8°. In the ZnFeLDH structures, the same differences have been observed, changing the positions of the peaks from 59.2 to 58.3°for (110) and from 61.3 to 59.7°for (113) by increasing the Fe content. Figueiredo et al. 17 described the same behavior of the (110) and (113) peaks with the increasing amount of iron in the structure, i.e., their position moved to lower 2θ angles. This result also confirms the goodness of our calculations, since in addition to reproducing the structure, the same behavior has been observed with the increase of the iron content in the LDH structure, as described in previous works, i.e., the increasing the parameter a and b axes and the displacement toward lower angles of the peaks (110) and (113) in the diffractograms.

CONCLUSIONS
Our DFT calculations reproduce the crystal structures of Mg:Al 12−x Fe x -LDH and Zn:Al 12−x Fe x -LDH with different amounts of Fe 3+ cation. The high spin states are the most probable in these systems. A certain preference for clustering of The first column indicates the nature of the three cations to which O is bound. X indicates Mg or Zn cations, depending on the structure considered. Distances are given in Å. Fe 3+ cations can be observed in these systems; however, this tendency is not clear probably due to the high content of Mg. Similar behavior was observed in phyllosilicates with a high proportion of Mg. 20 These DFT calculations have reproduced the linear relationships of the cell parameters increasing with the relative proportion of Fe 3+ cations observed experimentally in both systems, Mg-LDH and Zn-LDH. The in-depth study of the structural parameters of the calculated LDHs has allowed us to determine that the presence of iron increases the length of the bonds to which it is linked, affecting both M−O and O−H bonds.
Moreover, the displacement of the (110) and (113) peaks to lower 2θ angles with the increasing amount of iron in the LDH structure have been also observed, reproducing the experimental results previously reported.