Computational Study of the Adsorption of Phosphates as Wastewater Pollutant Molecules on Faujasites

: The adsorption of sodium dihydrogen phosphate (NaH 2 PO 4 ) onto X- and Y-type faujasite zeolites was computationally studied using the Density Functional Theory (DFT) method. The structures were capturing pollutant molecules in wastewater. The results are in agreement with the experimental information concerning the inﬂuence of pH on the adsorption activity of phosphate adsorption on zeolites.


Introduction
Phosphorus is a mineral nutrient that is essential for all living species; however, its excessive concentration when dissolved in natural water sources, which is caused by untreated waste from agricultural, domestic, and industrial sources, causes eutrophication problems [1][2][3][4][5]. Common phosphate species with natural pH conditions are H 2 PO 4 − and HPO 4 2− [6]. Recently, the need to remove phosphate from water has become a hot research topic, and a variety of methods have been developed, mainly involving biological, physical, and chemical treatments [7,8]. The disadvantages of traditional treatments are their low selectivity towards phosphate over competing anions and other dissolved organic compounds, their low adsorption capacity in the neutral pH range, the deactivation of their adsorption capacity, and their ineffective regeneration [9][10][11].

Materials and Methods
Density Functional Theory (DFT) calculations were performed using the Vienna Ab initio Simulation Package (VASP) code [13,27]. The Density Functional Theory is based on the Hohenberg-Kohn theorems, which proved that the electron density of the ground state contains the information of an electronic system. DFT is widely used to describe the structural and electronic properties of many-body systems, solving Schrödinger's equation with the Kohn-Sham scheme. In particular, VASP code uses Bloch's theorem to calculate systems with periodic conditions. In this study, the Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional was used for generalized gradient approximation (GGA) [28]. The plane-augmented wave (PAW) method was used to describe the electronion interactions [29]. The plane wave cutoff energy was set at 415 eV. The Kohn-Sham equations were solved by a self-consistent field method until the energy difference was lower than 10 −4 eV. The Brillouin zone integration was sampled using only the Γ-point [30]. The atomic positions were fully relaxed until the forces in each atom were lower than 0.4 eV/Å. Gaussian smearing, with σ = 0.10 eV, was used for the band occupancies and to improve total energy convergence.
A (3 × 3 × 3) cell faujasite supercage was employed as a model; this framework cell was constructed from a cubic cell available in the Materials Studio software, and was later fully optimized. Initially, the cell shape and cell volume of both zeolites models were optimized; then, the atomic positions of the sodium dihydrogen phosphate, of the isolated X-and Y-type faujasite zeolites models, and of the supersystem (adsorbate-zeolite) were fully optimized. The optimized lattice parameters and the volume of the cubic unit cell (96 Si or Al atoms and 384 O atoms) were: a = b = c = 25.028 Å, α = β = γ = 90 • , and V = 15.678 Å 3 . To obtain the FAU-X zeolite model (Si/Al = 1.2), 12 Si atoms were substituted for 12 Al atoms. The negative charge of each of the Al atoms was compensated with a proton (H + ), meaning that 84 H + were added to the structure. The resulting molecular formula of the FAU-X zeolite was (H + ) H 84 A l84 Si 108 O 384 . Similarly, to obtain the FAU-Y zeolite (Si/Al = 2.5), 42 Si atoms were substituted for 42 Al atoms. The negative charge generated by the substitution was compensated by one proton (H + ) so that there were 54 protons in the structure. The resulting molecular formula of the FAU-Y zeolite was H 54 Al 54 Si 138 O 384 . The Al atoms were randomly distributed according to the Lowenstein rule [31]. Figure 1 shows the cubic unit cell of the faujasite zeolite obtained from the Materials Studio software, and Figures 2 and 3 describe the primitive triclinic cell of the FAU-X and FAU-Y zeolites, respectively. Processes 2021, 9, x FOR PEER REVIEW 3 of 12 state contains the information of an electronic system. DFT is widely used to describe the structural and electronic properties of many-body systems, solving Schrödinger's equation with the Kohn-Sham scheme. In particular, VASP code uses Bloch's theorem to calculate systems with periodic conditions. In this study, the Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional was used for generalized gradient approximation (GGA) [28]. The plane-augmented wave (PAW) method was used to describe the electron-ion interactions [29]. The plane wave cutoff energy was set at 415 eV. The Kohn-Sham equations were solved by a self-consistent field method until the energy difference was lower than 10 −4 eV. The Brillouin zone integration was sampled using only the -point [30]. The atomic positions were fully relaxed until the forces in each atom were lower than 0.4 eV/Å . Gaussian smearing, with σ = 0.10 eV, was used for the band occupancies and to improve total energy convergence. A (3 × 3 × 3) cell faujasite supercage was employed as a model; this framework cell was constructed from a cubic cell available in the Materials Studio software, and was later fully optimized. Initially, the cell shape and cell volume of both zeolites models were optimized; then, the atomic positions of the sodium dihydrogen phosphate, of the isolated X-and Y-type faujasite zeolites models, and of the supersystem (adsorbate-zeolite) were fully optimized. The optimized lattice parameters and the volume of the cubic unit cell (96 Si or Al atoms and 384 O atoms) were: a = b = c = 25.028 Å, α = β = γ = 90°, and V = 15.678 Å 3 . To obtain the FAU-X zeolite model (Si/Al = 1.2), 12 Si atoms were substituted for 12 Al atoms. The negative charge of each of the Al atoms was compensated with a proton (H + ), meaning that 84 H + were added to the structure. The resulting molecular formula of the FAU-X zeolite was (H + ) H84Al84Si108O384. Similarly, to obtain the FAU-Y zeolite (Si/Al = 2.5), 42 Si atoms were substituted for 42 Al atoms. The negative charge generated by the substitution was compensated by one proton (H + ) so that there were 54 protons in the structure. The resulting molecular formula of the FAU-Y zeolite was H54Al54Si138O384. The Al atoms were randomly distributed according to the Lowenstein rule [31]. Figure 1 shows the cubic unit cell of the faujasite zeolite obtained from the Materials Studio software, and Figures 2 and 3 describe the primitive triclinic cell of the FAU-X and FAU-Y zeolites, respectively.   Crystallographic data show that, of the possible acid sites, site II is the most occupied [19,32] and is located in the six rings of the sodalite cage from the supercage side [33]. This site is the most accessible to several adsorbates [34,35], therefore, in this study, site II was considered as an adsorption active site in the SDHP molecule.

Adsorption Energy
The adsorption energy of SDHP on FAU-X and FAU-Y zeolites was calculated using the following Equation [36][37][38]: where: EFAU-SHDP: The energy of the zeolite molecular supersystem.  Crystallographic data show that, of the possible acid sites, site II is the most occupied [19,32] and is located in the six rings of the sodalite cage from the supercage side [33]. This site is the most accessible to several adsorbates [34,35], therefore, in this study, site II was considered as an adsorption active site in the SDHP molecule.

Adsorption Energy
The adsorption energy of SDHP on FAU-X and FAU-Y zeolites was calculated using the following Equation [36][37][38]: where: EFAU-SHDP: The energy of the zeolite molecular supersystem. Crystallographic data show that, of the possible acid sites, site II is the most occupied [19,32] and is located in the six rings of the sodalite cage from the supercage side [33]. This site is the most accessible to several adsorbates [34,35], therefore, in this study, site II was considered as an adsorption active site in the SDHP molecule.

Adsorption Energy
The adsorption energy of SDHP on FAU-X and FAU-Y zeolites was calculated using the following Equation [36][37][38]: where: E FAU-SHDP : The energy of the zeolite molecular supersystem.

Molecular Optimization
The P-O bond distances were 1.5 Å and 1.6 Å in the optimized molecular structure of sodium dihydrogen phosphate (NaH 2 PO 4 ), the P-OH bond distance was 1.7 Å, and the O-H bond distance was 1.0 Å. The distance between the Na atom and the P atom was 2.7 Å. The angle formed between the oxygen atoms (121 • ) was greater than the angles formed by the oxygen atoms that make up the hydroxyl groups (109 • and 107 • ). Figure 4 shows the optimized structure of SDHP.

Molecular Optimization
The P-O bond distances were 1.5 Å and 1.6 Å in the optimized molecular structure of sodium dihydrogen phosphate (NaH2PO4), the P-OH bond distance was 1.7 Å , and the O-H bond distance was 1.0 Å . The distance between the Na atom and the P atom was 2.7 Å . The angle formed between the oxygen atoms (121°) was greater than the angles formed by the oxygen atoms that make up the hydroxyl groups (109° and 107°). Figure 4 shows the optimized structure of SDHP.

Optimized X-and Y-Type Faujasite Zeolites
The framework model was a unit cell of approximately 25 Å in each side, and was formed by a 192 XO4 (X = Si or Al) tetrahedral linked by oxygen atoms, in agreement with bibliographic information [39,40]. The cell parameter and volume before and after a full optimization of the faujasite structure are shown in Table 1. The unit cell of the Y-type faujasite zeolite included 54 Al atoms and 138 Si atoms, and the X-type faujasite zeolite included 84 Al atoms and 108 Si atoms; the optimized cell volume of FAU-X was 15.286 Å 3 , and that of FAU-Y was 14.975 Å 3 .
Therefore, as the Si atoms were substituted for Al atoms, the cell's structure and dimensions were modified. The higher the Si/Al ratio, the more the unit cell's dimensions were increased and the more uniform the structure. Both effects were observed in experimental studies [41,42]. The diameter of the central pore was considered as the length between two oxygen atoms located at opposite extremes of the pore of the optimized framework cell, as shown in Figure 5. The diameter for FAU-X was 12.41 Å , and for FAU-Y it was 12.47 Å ; these values were reasonably comparable with experimental data, which indicated that the pore diameters vary between 6 Å and 12 Å [43,44]. An O-H bond distance of 0.97 Å was

Optimized X-and Y-Type Faujasite Zeolites
The framework model was a unit cell of approximately 25 Å in each side, and was formed by a 192 XO 4 (X = Si or Al) tetrahedral linked by oxygen atoms, in agreement with bibliographic information [39,40]. The cell parameter and volume before and after a full optimization of the faujasite structure are shown in Table 1. The unit cell of the Y-type faujasite zeolite included 54 Al atoms and 138 Si atoms, and the X-type faujasite zeolite included 84 Al atoms and 108 Si atoms; the optimized cell volume of FAU-X was 15.286 Å 3 , and that of FAU-Y was 14.975 Å 3 . Therefore, as the Si atoms were substituted for Al atoms, the cell's structure and dimensions were modified. The higher the Si/Al ratio, the more the unit cell's dimensions were increased and the more uniform the structure. Both effects were observed in experimental studies [41,42].
The diameter of the central pore was considered as the length between two oxygen atoms located at opposite extremes of the pore of the optimized framework cell, as shown in Figure 5. The diameter for FAU-X was 12.41 Å, and for FAU-Y it was 12.47 Å; these values were reasonably comparable with experimental data, which indicated that the pore diameters vary between 6 Å and 12 Å [43,44]. An O-H bond distance of 0.97 Å was obtained in the Brønsted acid sites of the optimized zeolite models; this value was very close to the experimental and computational data [45,46]. Processes 2021, 9, x FOR PEER REVIEW 6 of obtained in the Brønsted acid sites of the optimized zeolite models; this value was ver close to the experimental and computational data [45,46].
(a) (b) The sodalite cages in the faujasite-type zeolites were connected to the tetrahedr through hexagonal prisms. These sodalite units formed so-called β cavities, which we accessible through cavities with a diameter of 2.0 Å to 2.5 Å , and formed by the oxyge atoms of the unshared hexagonal faces [43,47,48]. In the FAU-X structure, the distan between the oxygen atoms of the hexagonal prism was 3.80 Å ; likewise, the diameter the largest sodalite pore was 4.70 Å . The pore dimensions of the FAU-Y structure in th hexagonal prism and the sodalite were 3.7 Å and 5.2 Å , respectively, as shown in Figu 6.  The sodalite cages in the faujasite-type zeolites were connected to the tetrahedral through hexagonal prisms. These sodalite units formed so-called β cavities, which were accessible through cavities with a diameter of 2.0 Å to 2.5 Å, and formed by the oxygen atoms of the unshared hexagonal faces [43,47,48]. In the FAU-X structure, the distance between the oxygen atoms of the hexagonal prism was 3.80 Å; likewise, the diameter of the largest sodalite pore was 4.70 Å. The pore dimensions of the FAU-Y structure in the hexagonal prism and the sodalite were 3.7 Å and 5.2 Å, respectively, as shown in Figure 6.
Processes 2021, 9, x FOR PEER REVIEW 6 of 12 obtained in the Brønsted acid sites of the optimized zeolite models; this value was very close to the experimental and computational data [45,46]. The sodalite cages in the faujasite-type zeolites were connected to the tetrahedral through hexagonal prisms. These sodalite units formed so-called β cavities, which were accessible through cavities with a diameter of 2.0 Å to 2.5 Å , and formed by the oxygen atoms of the unshared hexagonal faces [43,47,48]. In the FAU-X structure, the distance between the oxygen atoms of the hexagonal prism was 3.80 Å ; likewise, the diameter of the largest sodalite pore was 4.70 Å . The pore dimensions of the FAU-Y structure in the hexagonal prism and the sodalite were 3.7 Å and 5.2 Å , respectively, as shown in Figure  6.

Adsorption of Phosphate on X-and Y-Type Faujasites
Two models were used to assess the effect of the substitution of Si atoms for Al atoms on phosphate adsorption activity: a model with a Si/Al ratio of 1.2 for FAU-X, and a model with a Si/Al ratio of 2.5 for FAU-Y. The calculated adsorption energy values for phosphate on both zeolites are shown in Table 3.   Table 2 shows the calculated structural parameters of three computational method Experimental X-ray crystallographic data are included [45].

Adsorption of Phosphate on X-and Y-Type Faujasites
Two models were used to assess the effect of the substitution of Si atoms for Al atom on phosphate adsorption activity: a model with a Si/Al ratio of 1.2 for FAU-X, and a mod with a Si/Al ratio of 2.5 for FAU-Y. The calculated adsorption energy values for phospha on both zeolites are shown in Table 3. Table 3. Adsorption energy of NaH2PO4 on FAU-X and FAU-Y.

Structure
Adsorption Energy   Table 2 shows the calculated structural parameters of three computational methods. Experimental X-ray crystallographic data are included [45].

Adsorption on the FAU-X Model
In the most favored position of the phosphate on FAU-X, with the lowest adsorption energy, the oxygen atoms of the molecule were close to the hydrogen atoms of the pore, forming a bidentate complex, as shown in Figure 8. The bond distance was 1.5 Å. This indicated a possible electrostatic attraction between FAU-X and NaH 2 PO 4 . This interaction plays an important role in the catalytic activity [49,50]. In the most favored position of the phosphate on FAU-X, with the lowest adsorption energy, the oxygen atoms of the molecule were close to the hydrogen atoms of the pore, forming a bidentate complex, as shown in Figure 8. The bond distance was 1.5 Ǻ. This indicated a possible electrostatic attraction between FAU-X and NaH2PO4. This interaction plays an important role in the catalytic activity [49,50]. The adsorption energy can be used as a criterion to describe the interaction mechanism between two systems [4]. In this case, an adsorption energy value of 31.64 kJ/mol was obtained, suggesting an exothermic reaction. Similar DFT studies of H2PO4 -adsorption as bidentate complexes report an energy of 94.4 kJ/mol at low pH conditions on hydrated Fe oxide [51], and of 22.5 kJ/mol in modified polyether sulfone [1].
Experimentally, it is known that positively charged adsorbents can interact effectively with negatively charged H2PO4 -through electrostatic attractions, concluding that, in acidic conditions, phosphate adsorption activity is favored. Adsorption tests in the pH range from 2 to 9 indicate that adsorption activity is greater in the pH range between 4 and 6. In the strongly acidic pH range, phosphate adsorption activity decreases; this could be attributed to the formation of a weak hydrophosphate salt [1,[51][52][53][54]. In fact, the FAU-X structure is more acidic than FAU-Y, as the latter contains more protons that compensate for the negative charge resulting from the substitution of silicon for aluminum. Therefore, the adsorption energy described for the FAU-X system is different than that for FAU-Y.
It is also known that, at pH values of above 7, there is a reduction in the adsorption of H2PO4 -, which is attributed to the fact that the surface charge of the adsorbent becomes more negative at higher pH values, causing greater electrostatic repulsion towards the negatively charged phosphate anions [55][56][57][58].

Adsorption on FAU-Y Model
The more Si atoms are substituted by Al atoms, the more the protons compensate for the negative charge, therefore, the FAU-Y model has fewer H + that compensate for the negative charge generated by the Al atoms in the system. In this case, an oxygen atom of the molecule was placed on a hydrogen atom at site II of the zeolite; this conformation allowed for the placing of two of the hydrogen atoms of the phosphate onto the oxygen atoms of the zeolite model, as shown in Figure 9. The adsorption energy can be used as a criterion to describe the interaction mechanism between two systems [4]. In this case, an adsorption energy value of 31.64 kJ/mol was obtained, suggesting an exothermic reaction. Similar DFT studies of H 2 PO 4 − adsorption as bidentate complexes report an energy of 94.4 kJ/mol at low pH conditions on hydrated Fe oxide [51], and of 22.5 kJ/mol in modified polyether sulfone [1].
Experimentally, it is known that positively charged adsorbents can interact effectively with negatively charged H 2 PO 4 − through electrostatic attractions, concluding that, in acidic conditions, phosphate adsorption activity is favored. Adsorption tests in the pH range from 2 to 9 indicate that adsorption activity is greater in the pH range between 4 and 6. In the strongly acidic pH range, phosphate adsorption activity decreases; this could be attributed to the formation of a weak hydrophosphate salt [1,[51][52][53][54]. In fact, the FAU-X structure is more acidic than FAU-Y, as the latter contains more protons that compensate for the negative charge resulting from the substitution of silicon for aluminum. Therefore, the adsorption energy described for the FAU-X system is different than that for FAU-Y.
It is also known that, at pH values of above 7, there is a reduction in the adsorption of H 2 PO 4 − , which is attributed to the fact that the surface charge of the adsorbent becomes more negative at higher pH values, causing greater electrostatic repulsion towards the negatively charged phosphate anions [55][56][57][58].

Adsorption on FAU-Y Model
The more Si atoms are substituted by Al atoms, the more the protons compensate for the negative charge, therefore, the FAU-Y model has fewer H + that compensate for the negative charge generated by the Al atoms in the system. In this case, an oxygen atom of the molecule was placed on a hydrogen atom at site II of the zeolite; this conformation allowed for the placing of two of the hydrogen atoms of the phosphate onto the oxygen atoms of the zeolite model, as shown in Figure 9.
After full optimization, the O(SHDP)-H(FAU-Y) bond distance was 1.4 Å, and the hydrogen atoms of the molecule were located at a distance of 2.3 Å and 2.1 Å, respectively. In addition, the protons of the hydroxyl groups of the molecule interacted with the nearest oxygen atoms of FAU-Y, and one of the oxygen atoms was vertically located on an H atom of this zeolite. Likewise, the O-H bond distance of FAU-Y, closer to the NaH 2 PO 4 molecule, was 1.1 Å, which suggested an interaction between an O atom of the NaH 2 PO 4 molecule and an H atom of FAU-Y. A similar computational study of H 3 PO 4 adsorption on ZSM-5 showed analogous behavior [3,59].
Processes 2021, 9, x FOR PEER REVIEW 9 of 12 Figure 9. Optimized structure of the shaped system for the adsorbent (FAU-Y) and the adsorbate (NaH2PO4).
After full optimization, the O(SHDP)-H(FAU-Y) bond distance was 1.4 Ǻ, and the hydrogen atoms of the molecule were located at a distance of 2.3 Ǻ and 2.1 Ǻ, respectively. In addition, the protons of the hydroxyl groups of the molecule interacted with the nearest oxygen atoms of FAU-Y, and one of the oxygen atoms was vertically located on an H atom of this zeolite. Likewise, the O-H bond distance of FAU-Y, closer to the NaH2PO4 molecule, was 1.1 Ǻ, which suggested an interaction between an O atom of the NaH2PO4 molecule and an H atom of FAU-Y. A similar computational study of H3PO4 adsorption on ZSM-5 showed analogous behavior [3,59].
The hydroxyl groups of sodium dihydrogen phosphate form bonds between the hydrogen atoms and the oxygen atoms of the structure, which act as alkaline sites, and these bonds have an important impact on the stability of the adsorption complexes. In another study, an adsorption energy of 161 kJ/mol was obtained from Mg/Ca-modified biochar structures [4,60]. Furthermore, the adsorption energies of H2PO4 -vary in the range of 265.35-317.66 kJ/mol in experiments with hydroxylated alfa-SiO2 [61].

Conclusions
The Y-type faujasite zeolite shows greater activity for the adsorption of NaH2PO4 molecules, with an adsorption energy of 161 kJ/mol, compared to the X-type faujasite zeolite system, with an adsorption energy of 31.64 kJ/mol. Both Brønsted acidic and Lewis alkaline sites favor the adsorption activity of zeolites. The optimized phosphate molecular structure and unit cell geometrical parameters of the zeolites are in agreement with experimental data.
The terminal hydrogen atoms of zeolites are found at a distance of 0.97 Å from the oxygen atoms. This value increases to 1.1 Å for hydroxides interacting with phosphate. The electrostatic interactions of the phosphate and the oxygen atoms of the zeolite provokes changes in the arrangement of the hydrogen atoms of the molecule.
The supercage or the central pore of zeolites are accessible sites for the interaction of the phosphate molecule.
The corroboration of the obtained computational results with experimental data will be very interesting, as will the complete analysis of all the variables involved in the complex process of decreasing the pollutants of wastewater at the theoretical and experimental level. This paper only reports the results obtained from dihydrogen phosphate; however, the study of the adsorption of other phosphate species on these models, and also on other adsorbents, is very interesting and promising. The hydroxyl groups of sodium dihydrogen phosphate form bonds between the hydrogen atoms and the oxygen atoms of the structure, which act as alkaline sites, and these bonds have an important impact on the stability of the adsorption complexes. In another study, an adsorption energy of 161 kJ/mol was obtained from Mg/Ca-modified biochar structures [4,60]. Furthermore, the adsorption energies of H 2 PO 4 − vary in the range of 265.35-317.66 kJ/mol in experiments with hydroxylated alfa-SiO 2 [61].

Conclusions
The Y-type faujasite zeolite shows greater activity for the adsorption of NaH 2 PO 4 molecules, with an adsorption energy of 161 kJ/mol, compared to the X-type faujasite zeolite system, with an adsorption energy of 31.64 kJ/mol. Both Brønsted acidic and Lewis alkaline sites favor the adsorption activity of zeolites. The optimized phosphate molecular structure and unit cell geometrical parameters of the zeolites are in agreement with experimental data.
The terminal hydrogen atoms of zeolites are found at a distance of 0.97 Å from the oxygen atoms. This value increases to 1.1 Å for hydroxides interacting with phosphate. The electrostatic interactions of the phosphate and the oxygen atoms of the zeolite provokes changes in the arrangement of the hydrogen atoms of the molecule.
The supercage or the central pore of zeolites are accessible sites for the interaction of the phosphate molecule.
The corroboration of the obtained computational results with experimental data will be very interesting, as will the complete analysis of all the variables involved in the complex process of decreasing the pollutants of wastewater at the theoretical and experimental level. This paper only reports the results obtained from dihydrogen phosphate; however, the study of the adsorption of other phosphate species on these models, and also on other adsorbents, is very interesting and promising.
The results of this study contribute to the objective of previous experimental and computational studies that have provided interesting information to improve the understanding of the trapping of phosphates on zeolites as adsorbents.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.