Analysis of Local and Global Aromaticity in Si 3 C 5 and Si 4 C 8 Clusters. Aromatic Species Containing Planar Tetracoordinate Carbon

: The minimum energy structures of the Si 3 C 5 and Si 4 C 8 clusters are planar and contain planar tetracoordinate carbons (ptCs). These species have been classiﬁed, qualitatively, as global ( π ) and local ( σ ) aromatics according to the adaptive natural density partitioning (AdNDP) method, which is an orbital localization method. This work evaluates these species’ aromaticity, focusing on conﬁrming and quantifying their global and local aromatic character. For this purpose, we use an orbital localization method based on the partitioning of the molecular space according to the topology of the electronic localization function (LOC-ELF). In addition, the magnetically induced current density is analyzed. The LOC-ELF-based analysis coincides with the AdNDP study (double aromaticity, global, and local). Moreover, the current density analysis detects global and local ring currents. The strength of the global and local current circuit is signiﬁcant, involving 4n + 2 π - and σ -electrons, respectively. The latter implicates the Si-ptC-Si fragment, which would be related to the 3c-2e σ -bond detected by the orbital localization methods in this fragment.

Given the NICS problems mentioned above, in this work we analyze the local and global ring currents in the Si3C5 and Si4C8 systems. In addition, we analyzed the chemical bonding with the ELF-LOC method, an orbital localization scheme in the domains of the electronic localization function (ELF). The ELF-LOC, like AdNDP, provides information related to electronic delocalization, given its flexibility to identify orbitals distributed in more than two atomic centers. Our results confirm the presence of two delocalization circuits, one global (π) and one local (σ). More importantly, the global and local ring currents are diatropic and significant (according to the ring current strength, RCS, calculations). These results strongly support these systems' double aromatic character (local and global) and their role in their stability.

Computational Methods
Geometry optimizations were performed at B3LYP [65]/6-311+G* [66,67] level. Vibrational frequencies were evaluated at the same level to confirm the optimized structures as true minima on their potential energy surface using the Gaussian16 program [68]. Cartesian coordinates of the optimized structures are shown in Table S1.
A detailed description of the ELF-LOC algorithm can be found elsewhere [69][70][71][72]. In this work, all numerical determinations were performed at the B3LYP level within density functional theory (DFT), using the atomic STO-3G basis sets. The overlap integrals over the ELF regions, required to calculate the localized natural orbitals, were obtained from the GAMESS computational package [73] and a modified version of the ToPMoD program [74]. The orbital localization was performed using our codes [69,70].
Current densities were computed with the GIMIC program [75,76] using the gaugeincluding atomic orbital (GIAO) [77] method. In the calculations, the magnetic field is directed along the z-axis, i.e., perpendicular to the molecular plane. The unit for current susceptibility is nA T −1 and the results are, therefore, independent of the magnitude of the magnetic field. For a qualitative analysis, vector plots of the current density in a plane placed 0.0 and 0.5 Å above the molecular plane were generated. Diatropic (aromatic) and paratropic (antiaromatic) currents are assumed to circle clockwise and counterclockwise, respectively. Current pathways are visualized using Paraview [78,79]. The ring current strengths (RCS), a measure of the net current intensity around the molecular ring, were obtained after considering different integration planes (see Figures S1-S3). The integration planes correspond to cut-off planes perpendicular to the chosen bonds of the molecule and extend horizontally along the ring's plane in 3.6 Å, with 2.6 Å above and below the Scheme 1. Design strategy of aromatic ptC systems employed in reference [62].

Computational Methods
Geometry optimizations were performed at B3LYP [65]/6-311+G* [66,67] level. Vibrational frequencies were evaluated at the same level to confirm the optimized structures as true minima on their potential energy surface using the Gaussian16 program [68]. Cartesian coordinates of the optimized structures are shown in Table S1.
A detailed description of the ELF-LOC algorithm can be found elsewhere [69][70][71][72]. In this work, all numerical determinations were performed at the B3LYP level within density functional theory (DFT), using the atomic STO-3G basis sets. The overlap integrals over the ELF regions, required to calculate the localized natural orbitals, were obtained from the GAMESS computational package [73] and a modified version of the ToPMoD program [74]. The orbital localization was performed using our codes [69,70].
Current densities were computed with the GIMIC program [75,76] using the gaugeincluding atomic orbital (GIAO) [77] method. In the calculations, the magnetic field is directed along the z-axis, i.e., perpendicular to the molecular plane. The unit for current susceptibility is nA T −1 and the results are, therefore, independent of the magnitude of the magnetic field. For a qualitative analysis, vector plots of the current density in a plane placed 0.0 and 0.5 Å above the molecular plane were generated. Diatropic (aromatic) and paratropic (antiaromatic) currents are assumed to circle clockwise and counterclockwise, respectively. Current pathways are visualized using Paraview [78,79]. The ring current strengths (RCS), a measure of the net current intensity around the molecular ring, were obtained after considering different integration planes (see Figures S1-S3). The integration planes correspond to cut-off planes perpendicular to the chosen bonds of the molecule and extend horizontally along the ring's plane in 3.6 Å, with 2.6 Å above and below the ring. The two-dimensional Gauss-Lobatto algorithm [76,80] was used to integrate the current passing through an integration plane.
Vector plot visualization of the current density in the plane of the molecule and 0.5Å above are reported in Figures S1-S5. It is essential to mention that the negative (diatropic) and the positive (paratropic) NICS at the center of the molecules are associated with aromaticity and anti-aromaticity. In contrast, for the RCS, a positive (diatropic) and negative (paratropic) sign correspond to aromatic and anti-aromatic molecules, respectively. For both the NICS and the RCS, values close to zero suggest non-aromatic behavior [81]. Both the NICS and the current density analysis were performed at B3LYP/6-311+G* level.
The NICS were computed using the gauge-including atomic orbital (GIAO) [77] method and dissected into their core, σ, and π contributions, using the natural chemical shielding (NCS) [82] analysis as implemented in the NBO 6.0 program [83]. To evaluate the NICS, we used NICSall, a simple program developed in our group, which is interfaced with the Gaussian program. NICSall helps to prepare the inputs and submits the calculations to generate the data according to the user's requirements. In this work, we computed the NICS in 2D. To do this, we first estimated the box size, which in this case was defined by the sides equal to 1.5 times the length and width of the molecule (centered and placed in the XY plane) by the height (z-axis), taking their lowest value. NICSall prepares the inputs to fill the grid with a step size of 0.2 Å (this is a default value; it could be modified by the user). Finally, NICSall delivers the outputs: text files with the properties (scans, FiPC-NICS) or cube files to plot maps and isosurfaces. The NICS plots were performed with the VisIt 3.0.2 program [84].
It is essential to mention that both the ring currents analysis [76] and the ELF analysis [85,86] are not very dependent on the method or the basis set used in their calculation.

Results
For the sake of clarity, our analysis will be divided into two parts: the chemical bonding analysis according to the ELF-LOC method and the global and local aromaticity analysis according to magnetic criteria.

Chemical Bonding Analysis according to the ELF-LOC Method
The orbital localization provided by the ELF-LOC method reveals a chemical bonding pattern similar to that described by AdNDP, as shown in Figures 1 and 2. In Si 3 C 5 and Si 4 C 8 , the C 5 and C 8 rings are connected by 2-center 2-electron C-C σ-bonds (2c-2e). For Si 4 C 8 , one C-C σ-bond (2c-2e) that splits the C 8 ring into two pentagons is also detected. Additionally, delocalized bonds are detected, three π-bonds in the case of Si 3 C 5 and five in Si 4 C 8 . Note that in this set, the π-orbitals are distributed over the entire molecular structure. Moreover, delocalized σ-bonds (3c-2e) are also detected in each Si-ptC-Si triangular fragment. These results support global πand local σ-delocalization, suggesting possible global and local aromaticity according to Hückel's 4n + 2 rule [87][88][89].

Current Density Analysis
According to the magnetic criteria, in the presence of an external magnetic field, aromatic (antiaromatic) molecules sustain diatropic (paratropic) currents. In contrast, in nonaromatic molecules, the currents in one direction or the other cancel out, giving a resultant current strength close to zero [53][54][55]76,81,90,91]. In this work, we define "local ring current" as the current circuits distributed in a local molecular ring. In contrast, we describe the "global ring current" circuit as distributed around the whole molecule. This assignment was introduced by Sundholm et al. [61] and recently used by our group to highlight the shortcomings of the NICS in assessing the aromaticity of polycyclic systems [92]. Note that Aihara introduced a similar concept to distinguish between different current density pathways [93].

Current Density Analysis
According to the magnetic criteria, in the presence of an external magnetic field, aromatic (antiaromatic) molecules sustain diatropic (paratropic) currents. In contrast, in nonaromatic molecules, the currents in one direction or the other cancel out, giving a resultant current strength close to zero [53][54][55]76,81,90,91]. In this work, we define "local ring current" as the current circuits distributed in a local molecular ring. In contrast, we describe the "global ring current" circuit as distributed around the whole molecule. This assignment was introduced by Sundholm et al. [61] and recently used by our group to highlight the shortcomings of the NICS in assessing the aromaticity of polycyclic systems [92]. Note that Aihara introduced a similar concept to distinguish between different current density pathways [93].
We can glimpse the patterns of current flows by inspecting current density plots. Hence, the magnetically induced current density of Si3C5 and Si4C8 systems, calculated in a plane located 0.5 Å above the molecular plane, are shown in Figures 3 and 4 (part a). These plots guide our selection of the integration planes. The different contributions (bond, atomic and ring currents) are defined and quantified by analyzing the integration profiles along the integration planes, see Figures S1 and S2. In this way, it is possible to computationally determine the local and global character of the induced currents. The We can glimpse the patterns of current flows by inspecting current density plots. Hence, the magnetically induced current density of Si 3 C 5 and Si 4 C 8 systems, calculated in a plane located 0.5 Å above the molecular plane, are shown in Figures 3 and 4 (part a). These plots guide our selection of the integration planes. The different contributions (bond, atomic and ring currents) are defined and quantified by analyzing the integration profiles along the integration planes, see Figures S1 and S2. In this way, it is possible to computationally determine the local and global character of the induced currents. The strength of the local currents (diatropic) is computed and then subtracted from the strength of the diatropic current connecting the local rings to determine the global RCS [61,92]. Accordingly, it has been possible to identify the different delocalization circuits in the studied systems (part b of Figures 3 and 4). In addition, the intensity of each current is also shown in these figures. For the Si 3 C 5 system, an intense and paratropic current is detected inside the C 5 (RCS = −8.06 nA/T). This ring current is more paratropic than that exhibited by the cyclopentadienyl anion (−3.5 nA T −1 ), the results of which are shown in Figure S3. A global diatropic ring current with a strength of 12.2 nA T −1 is also detected. The latter result is from the delocalization of the six π-electrons and is weaker than that of the cyclopentadienyl anion (14.5 nA T −1 ). The differences may be due to the effect of the polarization of the π-electron cloud toward the silicon atoms. This hypothesis is supported by analyzing the Si 2 C 5 H 2 system, with one less bridged silicon atom exhibiting a paratropic (RCS = −3.6 nA T −1 ) and diatropic (RCS = 13.4 nA T −1 ) current intensity closer to those of the cyclopentadienyl anion (see Figures S3 and S5). Interestingly, a local diatropic ring current involving the Si-ptC-Si fragment is also detected, in complete agreement with that predicted by AdNDP and the ELF-LOC, which indicated the presence of a σ-delocalized 3c-2e bond. Moreover, the RCS value (5.2 nA/T) suggests that this species has a moderate local σ-aromatic character.
For the case of the Si 4 C 8 system (Figure 4), the current density analysis detects two paratropic and local ring currents (within each C 5 ring). Similar to Si 3 C 5 , these currents are more paratropic (RCS = −4.7 nA T −1 ) than those of the pentalene dianion (RCS = −4.4 nA T −1 ), whose current density analysis is shown in Figure S4. We also detect a weak global paratropic ring current (RCS = −1.6 nA T −1 ) distributed around the C 8 fragment in addition to an intense global diatropic current (RCS = 10.4 A T −1 ). The global diatropic current is less intense than those of the pentalene dianion (RCS = 10.4 A T −1 ), presumably because of the polarization of the π-electron cloud toward the silicon. Finally, two local sigma currents are detected, involving the Si-ptC-Si circuits, similar to that exhibited by Si 3 C 5 system. Moreover, these currents are of moderate intensity (RCS = 5.8 nA T −1 ), leaving in evidence the local sigma aromatic character of these species.
those of the cyclopentadienyl anion (see Figures S3 and S5). Interestingly, a local diatropic ring current involving the Si-ptC-Si fragment is also detected, in complete agreement with that predicted by AdNDP and the ELF-LOC, which indicated the presence of a σ-delocalized 3c-2e bond. Moreover, the RCS value (5.2 nA/T) suggests that this species has a moderate local σ-aromatic character.
For the case of the Si4C8 system (Figure 4), the current density analysis detects two paratropic and local ring currents (within each C5 ring). Similar to Si3C5, these currents are more paratropic (RCS = −4.7 nA T −1 ) than those of the pentalene dianion (RCS = −4.4 nA T −1 ), whose current density analysis is shown in Figure S4. We also detect a weak global paratropic ring current (RCS = −1.6 nA T −1 ) distributed around the C8 fragment in addition to an intense global diatropic current (RCS = 10.4 A T −1 ). The global diatropic current is less intense than those of the pentalene dianion (RCS = 10.4 A T −1 ), presumably because of the polarization of the π-electron cloud toward the silicon. Finally, two local sigma currents are detected, involving the Si-ptC-Si circuits, similar to that exhibited by Si3C5 system. Moreover, these currents are of moderate intensity (RCS = 5.8 nA T −1 ), leaving in evidence the local sigma aromatic character of these species.    Figure 5 shows the σ-and π-components of the NICSzz computed in both the molecular plane and a plane perpendicular to the molecular plane for the Si3C5 and Si4C8 systems. These plots are in complete agreement with the ring current analysis. The σ-component clearly depicts paratropic regions at the center of the C5 rings in both Si3C5 and Si4C8, in accord with the presence of paratropic ring currents inside these rings. In addition, an intense diatropic region (long-range) centered on the local Si-ptC-Si ring is observed supporting the presence of a local diatropic current in this region. The π-component exhibits a strong diatropic region around the whole Si3C5 and Si4C8 ring (long-range), in agreement with the presence of a global diatropic ring current. As indicated previously [61,92], the NICS analysis does not allow the identification of all the ring current circuits. However, it is a suitable complement to understand the magnetic behavior of these species and the interpretation of their aromaticity under the magnetic criterion. Figure 6 shows a classical NICSzz analysis, i.e., discrete measurements at the center of the ring and 1.0 A above. The values measured in the molecular plane are pretty large (presumably because of the coupling of bond contributions). In contrast, the values above the plane agree with the current strength values measured in the induced current density analysis. In addition, Figures S6-S8 show an identical NICSzz analysis for Si2C5H2, the cyclopentadienyl anion (C5H5 − ), and the pentalene dianion (C8H6 2− ). As seen with the current density analysis, diatropicity due to π-delocalization decreases for species with ptCs, presumably by the electron cloud polarization towards silicon atoms.  Figure 5 shows the σand π-components of the NICS zz computed in both the molecular plane and a plane perpendicular to the molecular plane for the Si 3 C 5 and Si 4 C 8 systems. These plots are in complete agreement with the ring current analysis. The σ-component clearly depicts paratropic regions at the center of the C 5 rings in both Si 3 C 5 and Si 4 C 8 , in accord with the presence of paratropic ring currents inside these rings. In addition, an intense diatropic region (long-range) centered on the local Si-ptC-Si ring is observed supporting the presence of a local diatropic current in this region. The π-component exhibits a strong diatropic region around the whole Si 3 C 5 and Si 4 C 8 ring (long-range), in agreement with the presence of a global diatropic ring current. As indicated previously [61,92], the NICS analysis does not allow the identification of all the ring current circuits. However, it is a suitable complement to understand the magnetic behavior of these species and the interpretation of their aromaticity under the magnetic criterion. Figure 6 shows a classical NICS zz analysis, i.e., discrete measurements at the center of the ring and 1.0 A above. The values measured in the molecular plane are pretty large (presumably because of the coupling of bond contributions). In contrast, the values above the plane agree with the current strength values measured in the induced current density analysis. In addition, Figures S6-S8 show an identical NICS zz analysis for Si 2 C 5 H 2 , the cyclopentadienyl anion (C 5 H 5 − ), and the pentalene dianion (C 8 H 6 2− ). As seen with the current density analysis, diatropicity due to π-delocalization decreases for species with ptCs, presumably by the electron cloud polarization towards silicon atoms.

NICS Analysis
Chemistry 2021, 3, FOR PEER REVIEW 8 Figure 5. Isolines of the σ-and π-components of NICSzz for the studied molecules (at the GIAO-B3LYP/6-311+G*//B3LYP/6-311+G* level). The isolines are plotted in both the molecular plane (left) and a plane perpendicular to the molecular plane (right). The color scale at the bottom is in ppm.

Conclusions
The global (π) and local (σ) aromaticity of Si 3 C 5 and Si 4 C 8 clusters have been reviewed, employing an orbital localization method (ELF-LOC) and magnetically induced current density analysis.
The ELF-LOC investigation leads to bond descriptions like those reported previously, based on the adaptive natural partitioning analysis (AdNDP) method. These systems show σand π-delocalization in agreement with Hückel's rule of 4n + 2 electrons.
The induced current density analysis identifies global and local aromatic current circuits. In the case of Si 3 C 5 , a local paratropic current (of moderate-intensity, RCS = −8.1 nA T −1 ) is detected inside the C 5 ring. However, a diatropic current of higher intensity (RCS = 12.16 nA T −1 ) is also present, which provides a global aromatic character to this species. A local σ-diatropic current involving the Si-ptC-Si fragment is also identified, with moderate intensity (RCS = 5.2 nA T −1 ). These results highlight the double aromatic character of this species. In the case of Si 4 C 8 , a weak global paratropic current (RCS = −1.5 nA T −1 ) and an intense global diatropic current (RCS = 10.4 nA T −1 ) are identified. In addition, two local diatropic currents around the Si-ptC-Si rings (RCS = 5.8 nA T −1 ) are also identified. This study highlights the double aromatic character of these clusters and the importance of this delocalization pattern in the stabilization of these species.
The ptC structures of the Si 3 C 5 and Si 4 C 8 clusters correspond to the global minima and, together with their double aromatic character, suggest that gas-phase experiments could detect these species.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/ 10.3390/chemistry3040080/s1, Table S1: Cartesian coordinates of the studied systems, Figure S1:  Figure S6: (a) Computed values (in ppm) of σand π-components (left and right sides) of NICS zz at and above the ring centers (local and global) of Si 2 C 5 H 2 (at the GIAO-B3LYP/6-311+G*//B3LYP/6-311+G* level). The blue/red dots denote diatropic/paratropic character, and the dot size is in line with the NICS zz magnitude. (b) Isolines of the σand π-components of NICSzz for Si 2 C 5 H 2 (at the GIAO-B3LYP/6-311+G*//B3LYP/6-311+G* level). The isolines are plotted in both the molecular plane (left) and a plane perpendicular to the molecular plane (right). The color scale at the bottom is in ppm, Figure S7: Isolines of the σand π-components of NICSzz for both the cyclopentadienyl anion (C 5 H 5 − ) and the pentalene dianion (C 8 H 6 2− ) (at the GIAO-B3LYP/6-311+G*//B3LYP/6-311+G* level). The isolines are plotted in both the molecular plane (left) and a plane perpendicular to the molecular plane (right). The color scale at the bottom is in ppm, Figure S8: Computed values (in ppm) of σand π-components (left and right sides) of NICS zz at and above the ring centers (local and global) of both the cyclopentadienyl anion (C 5 H 5 − ) and the pentalene dianion (C 8 H 6 2− ) (at the GIAO-B3LYP/6-311+G*//B3LYP/6-311+G* level). The blue/red dots denote diatropic/paratropic character, and the dot size is in line with the NICS zz magnitude.