Theoretical Study on ORR/OER Bifunctional Catalytic Activity of Axial Functionalized Iron Polyphthalocyanine

Designing efficient ORR/OER bifunctional electrocatalysts is very significant for reducing energy consumption and environmental protection. Hence, we studied the ORR/OER bifunctional catalytic activity of iron polyphthalocyanine (FePPc) coordinated by a series of axial ligands which has different electronegative coordination atom (FePPc-L) (L = -CN, -SH, -SCH3, -SC2H5, -I, -Br, -NH2, -Cl, -OCH3, -OH, and -F) in alkaline medium by DFT calculations. Among all FePPc-L, FePPc-CN, FePPc-SH, FePPc-SCH3, and FePPc-SC2H5 exhibit excellent ORR/OER bifunctional catalytic activities. Their ORR/OER overpotential is 0.256 V/0.234 V, 0.278 V/0.256 V, 0.280 V/0.329 V, and 0.290 V/0.316 V, respectively, which are much lower than that of the FePPc (0.483 V/0.834 V). The analysis of the electronic structure of the above catalysts shows that the electronegativity of the coordination atoms in the axial ligand is small, resulting in less distribution of dz2, dyz, and dxz orbitals near Ef, weak orbital polarization, small charge and magnetic moment of the central Fe atom, and weak adsorption strength for *OH. All these prove that the introduction of axial ligands with appropriate electronegativity coordinating atoms can adjust the adsorption of catalyst to intermediates and modify the ORR/OER bifunctional catalytic activities. This is an effective strategy for designing efficient ORR/OER bifunctional electrocatalysts.


Introduction
Due to the increase in fossil fuel emissions, the energy crisis and global environmental pollution problems have become increasingly prominent.At present, fuel cells are considered an effective substitute for fossil fuels, because they can protect natural resources and the environment [1][2][3].Oxygen reduction reaction (ORR) [4] and oxygen evolution reaction (OER) [5] are two important reactions in fuel cells.Up to now, it has been proven that Pt-based precious metals are the most efficient electrocatalysts for ORR [6,7].IrO 2 and RuO 2 are classic OER electrocatalysts [8,9].However, their high cost and low reserves hinder wide application.Thus, it is necessary to explore other non-precious catalysts to replace precious metal catalysts.People have been studying all kinds of non-precious catalysts of ORR and OER [10][11][12][13][14]. Transition metal phthalocyanine (TMPc)-based electrocatalysts have been of broad concern due to their low preparation cost, high conjugated structure, high thermal stability, chemical stability, and excellent catalytic performance [15].Since 1964, cobalt phthalocyanine molecule was used as an ORR catalyst [16], TMPcs have been extensively studied as superior ORR [17][18][19], OER [20][21][22], and ORR/OER [23][24][25][26] difunctional catalysts.Transition metal polyphthalocynine complexes (TMPPc) formed by multiple TMPcs show better conductivity, stability, and electrocatalytic activity than that of TMPc monomers.Previous research has proved that TMPPc can be used as promising ORR [27][28][29][30], OER [31,32], and ORR/OER difunctional electrocatalysts [33].For in-stance, a conjugated aromatic network two-dimensional material composed of TMPPc (TM = Fe, Co, Fe/Co) with a large-conjugated plane, highly exposed active sites, and strong conductivity exhibited outstanding ORR performance [27].It was reported that a thinner edge-anhydride-functionalized CoPPc with more exposed active sites and wider interlayer spacing showed better ORR catalytic performance [28].Meanwhile, polymeric azo linkage cobalt phthalocyanine was synthesized by the diazotization method.It was mixed with the benchmark catalyst IrO 2 and showed a good OER catalytic performance [31].In addition, it was reported a battery of polymeric cobalt phthalocyanines containing S or SO 2 linkers coated on various substrates showed good OER activity [32].Especially, FePPc has been recognized as a promising ORR electrocatalyst both experimentally and theoretically [27,29,34].However, FePPc is not a good OER electrocatalyst on account of the high OER overpotential [33].It is always the goal of people to explore or design a lowcost and efficient ORR/OER bifunctional electrocatalyst.Therefore, it is very meaningful to find an efficient strategy to make FePPc become an excellent ORR/OER bifunctional electrocatalyst.
Many studies have reported that introducing axial ligands to the active site is a very useful way to improve the electrocatalytic activity of MN 4 complexes.The axial ligands can change the electron density of the central metal to influence the binding energy between the central metal atom and the reaction intermediate [35,36].For example, a biomimetic iron phthalocyanine with pyridine functionalized carbon nanotubes as an axial ligand on the Fe atom was synthesized and exhibited better ORR catalytic activity compared with Pt/C catalyst.The higher ORR activity was mainly due to the adjustment of the electronic structure of the Fe atom by axial ligands [37].Three different penta-coordinated FePc electrocatalysts were successfully prepared and showed better ORR catalytic performance than the FePc [38].The influence of axial ligands on the ORR catalytic performance of FePPc was studied by introducing a series of axial ligands (L = -F, -Cl, -Br, -I, -OH et al.) to Fe atom of FePPc to form FePPc-L in the experiment.It was found that FePPc-I showed good stability and ORR catalytic performance [35].It was also proved that appropriate axial ligands can promote the ORR catalytic activity of FeN 4 [36].
Thus, it can be seen that the introduction of axial ligands with different electronegative coordination atoms on Fe atoms in FePPc can affect the electronic structure of the FeN 4 active center of FePPc in different degrees, thereby affecting the catalytic activity of FePPc.Therefore, in this study, based on DFT calculations, we introduced a series of axial ligands (L = -CN, -SH, -SCH 3 , -SC 2 H 5 , -I, -Br, -NH 2 , -Cl, -OCH 3 , -OH, and -F) with different electronegative coordination atoms on the Fe atom of FePPc to form FePPc-L and studied their ORR/OER catalytic activity.It is expected to understand how the electronegativity of coordination atoms in axial ligands affects the catalytic activity of the FePPc and design effective ORR/OER bifunctional catalysts.We look forward to the research that will be instructive to explore the efficient ORR/OER bifunctional electrocatalysts.

Structure and Stability
The top and side view of the optimized geometrical configuration of the unit cells of FePPc-L (L = -CN, -SH, -SCH 3 , -SC 2 H 5 , -I, -Br, -NH 2 , -Cl, -OCH 3 , -OH, and -F) are shown in Figure 1a,b.In the optimized unit cell of FePPc as shown in Figure S1 in the Supporting Information (SI), the lattice parameters and the average bond length of Fe-N bond are 10.66 Å and 1.93 Å, respectively, which is concordant with prior results [29,34].The lattice parameters of various FePPc-L are not much different from that of FePPc.However, average bond length of Fe-N (d Fe-N ) is slightly longer than that of FePPc as shown in Table S1.This is because the interaction between axial ligands and Fe atom weakens the Fe-N bond.The thermodynamic stability of catalysts usually can be assessed by the binding energy (E b ).
The negative E b means that the catalyst is thermodynamically stable, and it can be easily prepared in experiments.The E b of FePPc-L (L = -CN, -SH, -SCH 3 , -SC 2 H 5 , -I, -Br, -NH 2 , -Cl, -OCH 3 , -OH, and -F) are calculated according to the following formula: where E FePPc-L , E FePPc , and E L are the energy of FePPc-L, FePPc, and isolated axial ligand (L) in the vacuum calculated by DFT calculation.As seen in Figure 1c, the thermodynamic stability of all FePPc-L is good on account of their negative E b .Moreover, ORR and OER are performed in an alkaline medium, and there is competition between OH − and axial ligands in the adsorption with the central Fe atom [39].The E b calculated for FePPc-I is −1.74 eV, which is obviously larger than that of FePPc-OH (−2.50 eV) as seen in Figure 1c.However, FePPc-I was successfully prepared and exhibited higher stability in an alkaline medium and showed higher ORR catalytic activity than that of FePPc-OH [35].The reason may be that it is not easy for OH − to replace the axial ligand -I due to the comprehensive control of steric hindrance, reaction kinetics, thermodynamic equilibrium, etc. in real experiments, and FePPc-I showed good stability and high ORR catalytic activity in alkaline medium.This can be used as an analogy.Among all FePPc-L considered, FePPc-I had the largest E b , other ligands have stronger interaction with Fe atoms than that of FePPc-I.If FePPc-I can stably exist, so may other FePPc-L.
Molecules 2024, 29, x FOR PEER REVIEW 3 of 17 energy (Eb).The negative Eb means that the catalyst is thermodynamically stable, and it can be easily prepared in experiments.The Eb of FePPc-L (L = -CN, -SH, -SCH3, -SC2H5, -I, -Br, -NH2, -Cl, -OCH3, -OH, and -F) are calculated according to the following formula: where EFePPc-L, EFePPc, and EL are the energy of FePPc-L, FePPc, and isolated axial ligand (L) in the vacuum calculated by DFT calculation.As seen in Figure 1c, the thermodynamic stability of all FePPc-L is good on account of their negative Eb.Moreover, ORR and OER are performed in an alkaline medium, and there is competition between OH − and axial ligands in the adsorption with the central Fe atom [39].The Eb calculated for FePPc-I is −1.74 eV, which is obviously larger than that of FePPc-OH (−2.50 eV) as seen in Figure 1c.However, FePPc-I was successfully prepared and exhibited higher stability in an alkaline medium and showed higher ORR catalytic activity than that of FePPc-OH [35].The reason may be that it is not easy for OH − to replace the axial ligand -I due to the comprehensive control of steric hindrance, reaction kinetics, thermodynamic equilibrium, etc. in real experiments, and FePPc-I showed good stability and high ORR catalytic activity in alkaline medium.This can be used as an analogy.Among all FePPc-L considered, FePPc-I had the largest Eb, other ligands have stronger interaction with Fe atoms than that of FePPc-I.If FePPc-I can stably exist, so may other FePPc-L.

ORR/OER Mechanism of FePPc and FePPc-L
Previous studies [29,35] have proven that the ORR mechanism on FePPc and FePPc-L is the four-electron in an alkaline medium.The OER is the reversible process of the ORR, so we studied the ORR/OER bifunctional catalytic activity of FePPc and FePPc-L (L = -CN, -SH, -SCH3, -SC2H5, -I, -Br, -NH2, -Cl, -OCH3, -OH, and -F) based on the four-electron mechanism in the alkaline medium in this work.
The ORR in alkaline medium can be described as: The OER is the reversible process of the ORR, and the reaction process is described as:

ORR/OER Mechanism of FePPc and FePPc-L
Previous studies [29,35] have proven that the ORR mechanism on FePPc and FePPc-L is the four-electron in an alkaline medium.The OER is the reversible process of the ORR, so we studied the ORR/OER bifunctional catalytic activity of FePPc and FePPc-L (L = -CN, -SH, -SCH 3 , -SC 2 H 5 , -I, -Br, -NH 2 , -Cl, -OCH 3 , -OH, and -F) based on the four-electron mechanism in the alkaline medium in this work.
The ORR in alkaline medium can be described as: The OER is the reversible process of the ORR, and the reaction process is described as: R2': R4': where * denotes the adsorption site of the catalyst, and *OOH, *O, and *OH are the adsorbed intermediates.The free energy change of the ORR (∆Gx:x = 1-4) and of OER (∆Gx':x = 1-4) can be described by Equations ( 1)-( 8) and calculated based on computational hydrogen electrode (CHE) [40] according to the following equation: where ∆E, ∆ZPE, and ∆S are the changes in total energy, the zero point energy, and the entropy calculated by DFT.∆G U is the free energy change due to the applied electrode potential U, it was calculated according to the formula ∆G U = −neU, where e, n, and U are the transferred electron, the number of e, and the applied potential, respectively.∆G pH is the free energy change due to a change of H + concentration, it was obtained according to the formula ∆G pH = 2.303 K B T × pH (pH = 14 and T = 298.15K).
The potential-determination step (PDS) in the ORR and OER is defined as the reaction step with the maximum the free energy change [41].The ORR overpotential (η ORR ) and OER overpotential (η OER ) are obtained by Equations ( 10) and ( 11): The OER theoretical overpotential (η OER ) is defined as: where 0.401 V is the theoretical equilibrium potential when pH = 14 calculated based on the computational hydrogen electrode (CHE) model [40]. .The free energy of H 2 O(l) is derived as G H2O (l) = G H2O(g) + RT × ln(p/p 0 ) since only G H2O (g) can be directly obtained by DFT calculations, where R is the ideal gas constant, T = 298.15K, p = 0.035 bar, and p 0 = 1 bar.The free energy of O 2 (g) has been estimated as G O2 = 2G H2O(l) -2G H2 + 4.92 eV.The free energy of OH − was derived as G OH− = G H2O(l) -G H+ , G OH − = G H2O(l) -½G H2 -2.303K B T × pH.When pH = 14, ∆G = −1.605eV, ∆G' = 1.605 eV, the theoretical equilibrium potential is 0.401 V.
The Implicit solvation calculations were performed adopting a continuum solvation model of vaspsol [42] to consider the OER/ORR catalytic activity and electronic structure analysis of FePPc and FePPc-L.The free energy diagrams of ORR and OER of FePPc and FePPc-L are drawn as shown in Figures 2 and 3, respectively.

Adsorption Properties and OER/ORR Activity Descriptor
Adsorption-free energy can reflect the adsorption strength of different catalysts to the same substance.For good electrocatalysts, the adsorption strength of the reaction intermediates is moderate, which is conducive to the activation and desorption of the reaction intermediates [41].In ORR and OER, the adsorption-free energies of *OOH, *O, or *OH (∆G *OOH , ∆G *O , or ∆G *OH ) usually are used as the catalytic activity descriptor of the catalyst [46,47].∆G *OOH , ∆G *O , and ∆G *OH of the FePPc and FePPc-L were calculated by the following equations: In order to find the OER/ORR bifunctional catalytic activity descriptor of FePPc-L, the linear relationships between ∆G *OOH and ∆G *OH , ∆G *O , and ∆G *OH were calculated, and they were described as: ∆G *OOH = 0.88∆G *OH + 1.32 (R 2 = 0.98), and ∆G *O = 1.12∆G *OH + 0.52 (R 2 = 0.86), as shown in Figure 4a 4c and 4d, respectively.From Figure 4c, it can be found that the data does not follow a typical volcano curve, but a plateau.On the left branch of the curve, η ORR decreases with the decrease of adsorption strength of FePPc-L to intermediate *OH, the PDS is R4 (*OH + e − → * + OH − ).On the right branch of the curve, η ORR has no clear change trend with the increase of adsorption strength of FePPc-L to intermediate *OH, which is limited by the number of catalysts.The PDS is R1 (* + O 2 + H 2 O + e − → *OOH + OH − ).When ∆G *OH = 0.15 eV, this corresponds to the optimal ORR overpotentials of 0.25 V. FePPc-CN (η ORR = 0.256 V), FePPc-SH (η ORR = 0.278 V), FePPc-SCH 3 (η ORR = 0.280 V), and FePPc-SC 2 H 5 (η ORR = 0.278 V) are near the top of the volcanic curve and show better ORR catalytic performance.From Figure 4d, it can be seen that there is a similar volcano curve relationship between η

Electronic Properties Analysis
To deeply explore the origin of the OER/ORR bifunctional catalytic activity, we analyzed the electronic structure of FePPc-L, especially for FePPc-CN, FePPc-SH, FePPc-SCH3, and FePPc-SC2H5.
Firstly, we analyzed the Bader charge of the Fe atoms (QFe), the total Bader charge of four N atoms connected to Fe atom in FePPc and FePPc-L, and the total Bader charge of axial ligands (QL) in FePPc-L as shown in Table 2.It is found that the Fe atoms in FePPc-L showed smaller positive Bader charge in most of FePPc-L than that in FePPc (1.345|e|), especially for FePPc-CN, FePPc-SH, FePPc-SCH3, and FePPc-SC2H5.All QL of FePPc-L are negative, which indicates that the axial ligands exhibit an electron-withdrawing inductive effect.It is following the fact that the electronegativity of Fe (1.83) is smaller than that of coordinating atoms in axial ligands of FePPc-L: C (2.550) < S (2.580) < I (2.660) < Br (2.960) < N (3.040) < Cl (3.160) < O (3.440) < F (3.980) [49].Furthermore, we analyzed the spin magnetic moment of the Fe atom (MFe) in FePPc and FePPc-L.It is found that the MFe changes greatly in FePPc-L compared with that in FePPc (1.866 μB).For FePPc-CN, FePPc-SH, FePPc-SCH3, and FePPc-SC2H5, their smaller electronegativity of coordinating atoms in axial ligands corresponds to smaller MFe.Therefore, the electronegativity of coordinating atoms in axial ligands has a great influ-

Electronic Properties Analysis
To deeply explore the origin of the OER/ORR bifunctional catalytic activity, we analyzed the electronic structure of FePPc-L, especially for FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 .
Firstly, we analyzed the Bader charge of the Fe atoms (Q Fe ), the total Bader charge of four N atoms connected to Fe atom in FePPc and FePPc-L, and the total Bader charge of axial ligands (Q L ) in FePPc-L as shown in Table 2.It is found that the Fe atoms in FePPc-L showed smaller positive Bader charge in most of FePPc-L than that in FePPc (1.345|e|), especially for FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 .All Q L of FePPc-L are negative, which indicates that the axial ligands exhibit an electron-withdrawing inductive effect.It is following the fact that the electronegativity of Fe (1.83) is smaller than that of coordinating atoms in axial ligands of FePPc-L: C (2.550) < S (2.580) < I (2.660) < Br (2.960) < N (3.040) < Cl (3.160) < O (3.440) < F (3.980) [49].Furthermore, we analyzed the spin magnetic moment of the Fe atom (M Fe ) in FePPc and FePPc-L.It is found that the M Fe changes greatly in FePPc-L compared with that in FePPc (1.866 µ B ).For FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 , their smaller electronegativity of coordinating atoms in axial ligands corresponds to smaller M Fe .Therefore, the electronegativity of coordinating atoms in axial ligands has a great influence on the charge and magnetic moment of Fe atoms, which may affect the OER/ORR catalytic activity.
Table 2. Q Fe , Q N , and Q L represent the Bader charge of the Fe atom, the total Bader charge of four N atoms connected to Fe atom, and the total Bader charge of axial ligands in FePPc-L.The positive and negative bader charges present the charge accumulation and depletion, respectively.M Fe represents the magnetic moment of Fe atoms in FePPc and FePPc-L.We also analyzed the projected density of states (PDOS) of the Fe 3d orbitals (dx 2 -y 2 , dz 2 , dxz, dyz, and dxy) in FePPc and FePPc-L as shown in Figure 5. Compared with FePPc, for FePPc-L (L = -I, -Br, -OCH 3 , and -F), the polarization of the 3d orbitals is enhanced, which is consistent with bigger electronegativity of coordinating atoms in axial ligands, and the bigger magnetic moment of Fe atom (M Fe ).Polarization of the 3d orbitals is weaker for other FePPc-L (L = -CN, -SH, -SCH 3 , and -SC 2 H 5 ), which is consistent with smaller electronegativity of coordinating atoms in axial ligands and the smaller M Fe as shown in Figure 6a and Table 2.For FePPc-L, the dx 2 -y 2 , dz 2 , dyz, and dxz orbitals all split and became narrower.dz 2 , dyz, and dxz orbitals shift to a higher energy level relative to the fermi level (E f ) in different degrees.The distribution Polarization of dz 2 , dyz, and dxz orbitals near E f also decreases, especially for FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 .The ORR/OER intermediates OOH, O, and OH are adsorbed on the iron atom of FePPc-L from the side opposite to the axial ligand through the interaction between the P orbit of oxygen atom and the dz 2 , dyz and dxz orbits of iron atom.It is well known that the dz 2 , dyz and dxz orbitals of Fe atom extend in the axial direction in FePPc-L.The Axial ligands are bonded to Fe atoms by interaction with dz 2 , dyz, and dxz orbitals of the Fe atom, which affects the distribution of the dz 2 , dyz and dxz orbitals of the Fe atom.Axial ligands with different electronegative coordination atoms have different effects on dz 2 , dyz, and dxz orbitals of the Fe atom.The axial ligands with a smaller electronegativity of the coordination atom (-CN, -SH, -SCH 3 , and -SC 2 H 5 ) reduce the distribution of dz 2 , dyz, and dxz orbitals near E f , thereby weakening the adsorption strength of the catalyst to the intermediate and enhancing the catalytic activity of catalysts.This means that the dz 2 , dyz, and dxz orbitals play a major role in OER/ORR catalytic activity [50].Therefore, dz 2 , dyz, and dxz orbitals can be used as a descriptor of the catalytic activity of FePPc-L.What needs to be noted is that the ΔG*OH of FePPc-Cl and FePPc-OCH3 is smaller than that of other FePPc-L, which means the adsorption strength for *OH is stronger.However, their ICOHP(Fe-O) of *OH are not very big, which is probably because of their stronger solvation effect.As shown in Figure S2a-c, the ΔG*OOH, ΔG*OH, and ΔG*OH of FePPc and FePPc-L in aqueous solution all are smaller than that in vacuum, particularly for FePPc-Cl and FePPc-OCH3, all ΔG*OOH, ΔG*O, and ΔG*OH descend too much.This may be the reason why the η ORR of FePPc-Cl and FePPc-OCH3 are higher in aqueous solution as in Figure S2d.Mayer bond order is a very effective chemical bond analysis method [52,53].The Mayer bond order of the Fe-O bond of the *OH of all FePPc-L were also calculated.The growth trend of the Mayer bond order of the Fe-O bond and the electronegativity of coordinating atoms in axial ligands of FePPc-L is broadly consistent as shown in Figure 6a,d.It means that the small the electronegativity of coordinating atoms in axial ligands What needs to be noted is that the ∆G *OH of FePPc-Cl and FePPc-OCH 3 is smaller than that of other FePPc-L, which means the adsorption strength for *OH is stronger.However, their ICOHP(Fe-O) of *OH are not very big, which is probably because of their stronger solvation effect.As shown in Figure S2a-c, the ∆G *OOH , ∆G *OH , and ∆G *OH of FePPc and FePPc-L in aqueous solution all are smaller than that in vacuum, particularly for FePPc-Cl and FePPc-OCH 3 , all ∆G *OOH , ∆G *O , and ∆G *OH descend too much.This may be the reason why the η ORR of FePPc-Cl and FePPc-OCH 3 are higher in aqueous solution as in Figure S2d.
Mayer bond order is a very effective chemical bond analysis method [52,53].The Mayer bond order of the Fe-O bond of the *OH of all FePPc-L were also calculated.The growth trend of the Mayer bond order of the Fe-O bond and the electronegativity of coordinating atoms in axial ligands of FePPc-L is broadly consistent as shown in Figure 6a,d.It means that the small the electronegativity of coordinating atoms in axial ligands of FePPc-L corresponds the weak Fe-O bond of *OH, and a small adsorption strength for *OH.For FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 , the small Mayer bond order of the Fe-O bond of the *OH corresponds to the small electronegativity of coordinating atoms in axial ligands.

Computational Methods
Spin-polarized density functional theory calculations are performed using the Vienna ab initio simulation package (VASP5.4.4) [54,55] combined with the projector augmented wave (PAW) method [56].The electron-correction interactions were described by the generalized gradient approximation (GGA) of Perdew-Burke-Ernzerhof (GGA-PBE) [57].The van der Waals (vdW) interactions were treated via Grimme's DFT-D3 correction method [58].The plane wave cutoff energy was 500 eV.In the process of geometric optimization, the convergence criterion of force is 0.01 eVÅ −1 , energy is 10 −5 eV, and the vibrational frequency calculations are 10 −7 eV.A vacuum layer is set as 20 Å to prevent the interaction between adjacent layers.The k-point mesh based on the Monkhorst-Pack scheme [59] was set to 5 × 5 × 1 and 10 × 10 × 1 to optimize the structure and calculate the electronic structure, respectively.The implicit solvation calculations were performed adopting a continuum solvation model of vaspsol [42].The interaction between two atoms is analyzed by the Crystal orbital Hamilton population (COHP) [51] and the Mayer bond order [52,53].

Conclusions
In summary, through performing DFT calculations, we studied the ORR/OER bifunctional catalytic activity of iron polyphthalocyanine with a series of axial ligands (L = -CN, -SH, -SCH 3 , -SC 2 H 5 , -I, -Br, -NH 2 , -Cl, -OCH 3 , -OH, and -F) with different electronegative coordination atoms on the Fe atom (FePPc-L) in alkaline media.It is expected to understand how the electronegativity of coordination atoms in axial ligands affacts the catalytic activity of the FePPc-L and design effective ORR/OER bifunctional catalysts.The important results are as follows: Compared with other FePPc-L, FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 exhibit more excellent ORR/OER bifunctional catalytic activities.Their ORR/OER overpotential are 0.256 V/0.234 V, 0.278 V/0.256 V, 0.280 V/0.329 V, and 0.290 V/0.316 V, respectively, which are much lower than that of the FePPc (0.483 V/0.834 V).The analysis of the electronic structure of the above catalysts shows that the electronegativity of the coordination atoms in the axial ligand is small, resulting in less distribution of dz 2 , dyz, and dxz orbitals near E f , weak orbital polarization, small charge and magnetic moment of the central Fe atom, and weak adsorption strength for *OH.All these prove that the axial ligands with appropriate electronegative coordination atoms can affect the electronic structure of FeN 4 in FePPc, adjust the adsorption of catalyst to intermediates, and modify the ORR/OER bifunctional catalytic activities.The conclusions are generally in good agreement with our study goal.So this is an effective strategy of designing ORR/OER bifunctional electrocatalysts with low cost and high efficiency.

Supplementary Materials:
The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/molecules29010210/s1, Figure S1: The top and side view of the unit cell of FePPc.The green, gray, blue, and pink present H, C, N, and Fe atoms, respectively.

Figure 1 .
Figure 1.(a,b) The top and side view of the unit cell of FePPc-L, the green, gray, blue, and pink present H, C, N, and Fe atoms, respectively, and wine red presents the axial ligands.(c) Binding energy (Eb/eV) of FePPc-L.

Figure 1 .
Figure 1.(a,b) The top and side view of the unit cell of FePPc-L, the green, gray, blue, and pink present H, C, N, and Fe atoms, respectively, and wine red presents the axial ligands.(c) Binding energy (E b /eV) of FePPc-L.

Figure 3 .
Figure 3. Free energy diagrams of OER on FePPc and FePPc-L at U = 0 V and U = 0.401 V in alkaline solution.The potential-determination step (PDS) of the OER is marked with an orange upward arrow.The stable geometrical configurations of the intermediates (*OOH, *O, *OH) are also exhibited.The green, gray, blue, red, cyan, yellow, bright green, orange, purple, and pink present H, C, N, O, F, S, Cl, Br, I, and Fe atoms, respectively.(a-k) Free energy diagrams of OER on FePPc-L (L = -CN, -SH, -SCH 3 , -SC 2 H 5 , -I, -Br, -NH 2 , -Cl, -OCH 3 , -OH, and -F).(l) Free energy diagram of OER on FePPc.

17 Figure 4 .
Figure 4. (a) The linear relationship of ΔG*OOH vs. ΔG*OH.(b) The linear relationship of ΔG*O vs. ΔG*OH.(c) The relationship curve of η ORR and ΔG*OH.The blue vertical dotted line represents the most moderate ΔG*OH corresponding to the optimal ORR overpotential.(d) The relationship curve of η OER and ΔG*O-ΔG*OH.The blue vertical dotted line represents the most moderate ΔG*O-ΔG*OH corresponding to the optimal OER overpotential.

Figure 4 .
Figure 4. (a) The linear relationship of ∆G *OOH vs. ∆G *OH .(b) The linear relationship of ∆G *O vs. ∆G *OH .(c) The relationship curve of η ORR and ∆G *OH .The blue vertical dotted line represents the most moderate ∆G *OH corresponding to the optimal ORR overpotential.(d) The relationship curve of η OER and ∆G *O − ∆G *OH .The blue vertical dotted line represents the most moderate ∆G *O − ∆G *OH corresponding to the optimal OER overpotential.

Figure 6 .
Figure 6.(a) Electronegativity of coordinating atoms in axial ligands of various FePPc-L.(b) ΔG*OH of FePPc-L.(c) Integrated crystal orbital Hamilton population (ICOHP) of Fe-O bond in *OH of FePPc-L.(d) Mayer bond order of Fe-O bond in *OH of FePPc-L.In addition, the increasing trend of the electronegativity of coordinating atoms in axial ligands of FePPc-L is roughly consistent with the decreasing trend of the ΔG*OH of FePPc-L as shown in Figure6a,b.For FePPc-CN, FePPc-SH, FePPc-SCH3, and FePPc-SC2H5, their small electronegativity of coordinating atoms in axial ligands corresponds to their small adsorption strength for *OH.We further analyzed the bond strength of Fe-O bond in *OH of FePPc and FePPc-L based on crystal orbital Hamiltonian (COHP)[51] as shown in Figure 7.It is found that the bonding states distribution below Ef of Fe-O, dyz (Fe)-O, dz 2 (Fe)-O, and dxz (Fe)-O in *OH of FePPc-L are less than that in FePPc, the anti-bonding states below Ef of Fe-O bond, dyz (Fe)-O, dz 2 (Fe)-O, and dxz (Fe)-O in *OH are more than that in FePPc, especially for FePPc-CN, FePPc-SH, FePPc-SCH3, and FePPc-SC2H5.It is also illustrated by the integral COHP of the Fe-O [ICOHP (Fe-O)] in *OH of FePPc and FePPc-L as shown in Figure 6c.The increasing trend of the electronegativity of coordinating atoms in axial ligands of FePPc-L is roughly consistent with the increasing trend of the ICOHP (Fe-O) of *OH in FePPc-L as shown in Figure 6a,c.For FePPc-CN, FePPc-SH, FePPc-SCH3, and FePPc-SC2H5, the small electronegativity of coordinating atoms in axial ligands corresponds to the weak adsorption strength of *OH.It once again shows that the dz 2 , dyz, and dxz orbitals play a major role in OER/ORR catalytic activity.The electronegativity of the coordination atom in axial ligands can affect the distribution of the dz 2 , dyz, and dxz orbitals of the Fe atom, thereby weakening the adsorption strength of the catalyst to the intermediate and enhancing OER/ORR catalytic activity of catalysts.What needs to be noted is that the ΔG*OH of FePPc-Cl and FePPc-OCH3 is smaller than that of other FePPc-L, which means the adsorption strength for *OH is stronger.However, their ICOHP(Fe-O) of *OH are not very big, which is probably because of their

Figure 6 .
Figure 6.(a) Electronegativity of coordinating atoms in axial ligands of various FePPc-L.(b) ∆G *OH of FePPc-L.(c) Integrated crystal orbital Hamilton population (ICOHP) of Fe-O bond in *OH of FePPc-L.(d) Mayer bond order of Fe-O bond in *OH of FePPc-L.In addition, the increasing trend of the electronegativity of coordinating atoms in axial ligands of FePPc-L is roughly consistent with the decreasing trend of the ∆G *OH of FePPc-L as shown in Figure6a,b.For FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 , their small electronegativity of coordinating atoms in axial ligands corresponds to their small adsorption strength for *OH.We further analyzed the bond strength of Fe-O bond in *OH of FePPc and FePPc-L based on crystal orbital Hamiltonian (COHP)[51] as shown in Figure 7.It is found that the bonding states distribution below E f of Fe-O, dyz (Fe)-O, dz 2 (Fe)-O, and dxz (Fe)-O in *OH of FePPc-L are less than that in FePPc, the anti-bonding states below E f of Fe-O bond, dyz (Fe)-O, dz 2 (Fe)-O, and dxz (Fe)-O in *OH are more than that in FePPc, especially for FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 .It is also illustrated by the integral COHP of the Fe-O [ICOHP (Fe-O)] in *OH of FePPc and FePPc-L as shown in Figure 6c.The increasing trend of the electronegativity of coordinating atoms in axial ligands of FePPc-L is roughly consistent with the increasing trend of the ICOHP (Fe-O) of *OH in FePPc-L as shown in Figure 6a,c.For FePPc-CN, FePPc-SH, FePPc-SCH 3 , and FePPc-SC 2 H 5 , the small electronegativity of coordinating atoms in axial ligands corresponds to the weak adsorption strength of *OH.It once again shows that the dz 2 , dyz, and dxz orbitals play a major role in OER/ORR catalytic activity.The electronegativity of the coordination atom in axial ligands can affect the distribution of the dz 2 , dyz, and dxz orbitals of the Fe atom, thereby weakening the adsorption strength of the catalyst to the intermediate and enhancing OER/ORR catalytic activity of catalysts.

Figure 7 .
Figure 7. Crystal orbital Hamilton population (COHP) of the Fe-O bond in *OH of FePPc and FePPc-L.The horizontal red dotted line indicates the Ef.Positive contribution to the right represents the bonding states, while negative contribution to the left represents the anti-bonding states.

Figure 7 .
Figure 7. Crystal orbital Hamilton population (COHP) of the Fe-O bond in *OH of FePPc and FePPc-L.The horizontal red dotted line indicates the E f .Positive contribution to the right represents the bonding states, while negative contribution to the left represents the anti-bonding states.
Figure S2: (a-c) ∆G *OOH , ∆G *O , and ∆G *OH of FePPc and FePPc-L in vacuum and in aqueous solution.(d-f) η ORR , η OER , and η ORR + ηOER of FePPc and FePPc-L in vacuum and in aqueous solution.Table S1: The lattice parameters (in Å) of cell of FePPc and FePPc-L.Average bond length of Fe-N bond (d Fe-N , in Å) of FePPc and FePPc-L.

Table 1 .
where G * , G *OOH , G *O , and G *OH are the free energies of catalyst, adsorbed intermediates (*OOH, *O and *OH) in alkaline solution.G H 2 and G H 2 O are the free energies of H 2 and H Adsorption free energy intermediates (∆G *OOH , ∆G *O , and ∆G *OH , eV) in alkaline solution.