Open Access Article
This Open Access Article is licensed under a
Creative Commons Attribution 3.0 Unported Licence

Orbital contributions to magnetically induced current densities using gauge-including atomic orbitals

Rinat T. Nasibullin *, Maria Dimitrova , Rashid R. Valiev and Dage Sundholm *
Department of Chemistry, Faculty of Science, University of Helsinki, P. O. Box 55 (A. I. Virtanens plats 1), FIN-00014, Finland. E-mail: Dage.Sundholm@helsinki.fi; Rinat.Nasibullin@helsinki.fi

Received 23rd January 2025 , Accepted 31st March 2025

First published on 1st April 2025


Abstract

We have developed a method to calculate orbital contributions to magnetically induced current density (MICD) susceptibilities in molecules using gauge-including atomic orbitals (GIAO). The methods implemented in the GIMIC program have been used for analyzing orbital contributions to magnetically induced ring-current (MIRC) strengths. We have studied five aromatic, one nonaromatic, and four antiaromatic molecules. We show here that the contributions to the MIRC strength of all orbitals belonging to a given irreducible representation of the molecular point group in the presence of an external magnetic field are divergence free, whereas the MICD susceptibility of the individual orbitals are generally not divergence free. The largest contribution to the MIRC strength of antiaromatic molecules originates from the transition between the highest occupied molecular orbital (HOMO) and the lowest unoccupied molecular orbital (LUMO), whereas aromatic molecules have significant contributions involving many occupied orbitals. The MIRC contributions of σ orbitals are significant for planar molecules with a strained molecular structure.


1 Introduction

The aromaticity concept was originally introduced to characterize molecules whose electron-delocalization properties resemble those of benzene.1,2 The electron delocalization of aromatic molecules leads to enhanced thermodynamic stability, reduced bond-length alternation, particular reactivity, and typical magnetic and spectroscopic properties.3 The aromatic nature according to the magnetic criterion can be determined by measuring or calculating proton nuclear magnetic resonance (1H NMR) chemical shifts or by calculating various magnetic response properties such as nucleus-independent chemical shifts (NICS),4 magnetically induced current density (MICD) susceptibilities,5–11 and magnetically induced ring currents (MIRC).11–15

These computational approaches have to some extent been used for estimating contributions to the aromatic nature of individual molecular orbitals.16–23 A general and reliable method for calculating orbital contributions to MICD susceptibilities and molecular properties would be useful, since decomposition into orbital contributions provides insights into the electronic structure and can be used for predicting changes due to substitution, structural variations, or charge redistribution. Common approaches to qualitatively understand the role of the frontier molecular orbitals (MO) in aromatic bonds comprise analyses of orbital energies and spatial distributions of orbitals as well as using symmetry selection rules.24,25

Symmetry-based selection rules of the transitions between occupied and virtual orbitals have been used for assessing which occupied orbitals contribute to the MICD susceptibilities and magnetic properties. When the electronic transition between the highest-occupied molecular orbital (HOMO) and the lowest-unoccupied molecular orbital (LUMO) is electric-dipole (translationally) allowed, ring-shaped molecules may sustain a diatropic MIRC in the presence of an external magnetic field. Aromatic molecules can be identified using group theory. Molecules can be aromatic when the product of the irreducible representations of the HOMO, the LUMO and the translation operator contains the total symmetric irreducible representation. Analogously, molecules can be antiaromatic when the irreducible representation of the rotational operator multiplied by the irreducible representations of the HOMO and LUMO contains the total symmetric irreducible representation. The transition is then magnetically dipole (rotationally) allowed, which eventually leads to a paratropic MIRC.26–28 Nonaromatic molecules are in this context conjugated molecular ring that are neither aromatic nor antiaromatic. Even though a transition is formally translation or rotationally allowed, a molecule may still be nonaromatic due to a large spatial separation of the occupied and virtual orbitals or because the energy separation between them is large.

The symmetry analysis provides a qualitative interpretation of molecular aromatic properties, whereas the MIRC strengths provide quantitative information about the aromatic nature because molecules sustaining a net diatropic MIRC are aromatic, whereas antiaromatic molecules have a net paratropic MIRC.29 The MIRC strength can be used for determining how strongly aromatic or antiaromatic the molecule is.11

Orbital contributions to the MICD susceptibility have previously been calculated using the continuous transformation of the origin of the current density (ipsocentric) approach.7,10,30 In the ipsocentric approach, the origin of the vector potential of the external magnetic field is the position where the MICD susceptibility is calculated. Steiner and Fowler performed the first calculations of orbital contributions to the MICD susceptibility using the few-orbital model.24,31 They studied the MICD susceptibility of benzene with the ipsocentric approach using an early version of the SYSMOIC package.14 The ipsocentric approach employs atomic basis functions that do not depend on the magnetic perturbation, whereas gauge-including atomic orbitals (GIAO) also known as London atomic orbitals incorporate the gauge origin in the basis functions32,33 leading to a fast basis-set convergence of the MICD susceptibility and NMR shielding constants.34,35

The gauge-including magnetically induced currents (GIMIC) method employs GIAOs.11–13,15,28 GIMIC calculations need the one-electron density matrix and the three magnetically perturbed density matrices as input. The density matrices can be calculated at various levels of theory using common quantum chemistry program packages. Cartesian coordinates of the molecular structure and basis-set information are also needed as input data. The GIMIC method has been successfully applied in MICD susceptibility studies of aromatic and antiaromatic molecules including planar molecular rings,29 Möbius-twisted molecular rings,36,37 toroidal molecules,38 organometallic molecules,39–41 and inorganic molecular clusters.42,43

MICD susceptibilities calculated with GIMIC suffer from a tiny charge leakage, i.e., the divergence of the MICD susceptibility does not vanish when using finite basis sets. Divergence-free MIRC strengths can be calculated by performing numerical line integration of the induced magnetic field along the stagnation line of the MIRC vortex using the Ampère-Maxwell law.44 However, it has been shown that the charge-leakage error of GIMIC calculations is negligibly small already when using commonly used moderate-size basis sets and that it practically vanishes when using large basis sets.11 Recently, an expression for correcting the divergence was derived enabling calculations of divergence-free MICD susceptibilities using any basis set.45

NMR shielding tensors have been partitioned into contributions of various orbitals by performing natural chemical shielding (NCS) analyses.46 However, unitary transformations of the orbitals lead to spurious magnetic shielding contributions,25 which can be identified as unphysical core and valence-σ contributions to the MICD susceptibility and to the NMR shielding tensors due to the mixing of occupied orbitals belonging to the same irreducible representation.43,47 Orbital contributions to the MIRC strengths can also be estimated by partitioning the MICD susceptibility.25,48 Orbital contributions and atomic orbital contributions to NMR shielding tensors can then be calculated by numerical integration.49,50

In this work, we report an extension of the GIMIC method that enables calculations of orbital contributions to the MICD susceptibility. This is achieved by transforming the density matrices in the basis-function (atomic) basis to the canonical orbital basis where molecular orbital contributions can be identified. The selected orbital contributions to the density matrices are transformed back to the basis-function basis. The orbital contributions to the MICD susceptibility can be analyzed and visualized in the same way as the total MICD susceptibility.

We have studied eleven planar molecules and one nearly planar monocyclic molecule with different types of aromatic character. Six of them are aromatic (benzene (D6h), borazine (D3h), 1,5-dibora-2,4-diazabenzene (C2v), porphin (D2h), cyclopropane (D3h) and the cyclopropenium cation (D3h)), one is nonaromatic (1,4-cyclohexadiene, (D2h)), and four are antiaromatic (cyclooctatetraene (D4h and D2d), hexadehydro[12]annulene (D3h), tetraoxa-isophlorin (D2h) and cyclobutadiene (D2h)), where the molecular point group is given within parenthesis.

Cyclobutadiene and the cyclopropenium cation were studied because they have bond angles that significantly deviate from the standard bond angles of sp2 hybridized hydrocarbons leading to ring strains. These strained molecules, along with molecules like cyclopropane, cyclobutane, and cubane, sustain an MIRC in the valence σ orbitals leading to what is called σ aromaticity or σ antiaromaticity. The direction (tropicity) of the MIRC is determined by the number of electrons in the σ orbitals of the C–C bonds.29,51–59

The results for some of the molecules are compared to available data calculated using the ipsocentric approach. Orbital contributions to the MICD susceptibility of most of the studied molecules are reported here for the first time to the best of our knowledge.

The method we use for calculating orbital contributions to the MICD susceptibility is described in Section 2. We describe the employed computational levels in the Section 3. The obtained results are presented in Section 4. The main results are summarized in Section 5, which also contains the most important conclusions of our study.

2 Methodology

2.1 Decomposition of density matrices

The one-electron density matrix in the atomic-orbital basis (DAO) consist of contributions of the individual occupied orbitals
 
image file: d5sc00627a-t1.tif(1)
where Nocc is the number of occupied orbitals and ni is the occupation number of MO i, and Ci contains the coefficients of MO i. The orbital contributions to DAO (DAOi) can be identified by transforming it to the canonical-orbital basis DMO. The basis transformation is done by multiplying from right with SC, where S is the overlap matrix and C is matrix containing the MO coefficients, and with its transpose CTS from left
 
DMO = CTSDAOSC.(2)

The first-order magnetically perturbed density matrices in the AO basis for each Cartesian direction of the magnetic field (Bλ, λ ∈ {x, y, z}) can be obtained from the MO coefficients and the first-order change of the MO coefficients (Cλ) due to the external perturbation Bλ

 
image file: d5sc00627a-t2.tif(3)
Cλ can be expressed as a unitary (Uλ) transformation of the unperturbed occupied MO coefficients
 
image file: d5sc00627a-t3.tif(4)
which is the orbital response due to the external magnetic field Bλ in the limit of vanishing perturbation.60–62 The occupied–virtual block of Uλ (Uλ,iv) is obtained by solving the orbital response equations. The virtual–virtual block of Uλ does not contribute to the perturbed density matrices.63 The occupied–occupied block of Uλ (Uλ,ij) does not vanish when perturbation-dependent GIAOs are used.64 However, it can be determined i.e., by differentiating the orthonormalization condition60,61
 
image file: d5sc00627a-t4.tif(5)

Using the relation in eqn (4), the final expression for the unitary transformation matrix can be written as65

 
image file: d5sc00627a-t5.tif(6)
where SBλ is the first derivative of the overlap matrix with respect to the external magnetic field (Bλ) in the limit of vanishing perturbation. The contributions of orbital i to the perturbed density matrices can also be obtained by performing the same unitary transformation and orbital selection procedure as applied on the one-electron density matrix. The orbital contribution to the perturbed density matrices in the MO basis consists of one row and one column of the matrices. The common elements of two orbitals in the occupied–occupied block are shared equally between the two orbitals.

The orbital contributions to the density matrices are transformed to the AO basis. The dimension of the AO density matrices is Nbasis2, where Nbasis is the number of basis functions. The orbital contributions to the MICD susceptibility can then be calculated with GIMIC in the same way as the total MICD susceptibility by using the corresponding orbital contributions to the density matrices DAOλ,i and DAOi as input. Contributions to the MIRC strength of the ith MO are obtained by modifying the occupation matrix ni and performing the linear transformations and thereby avoiding time-consuming calculations.

The computational costs for calculating orbital contributions to the MICD susceptibility are comparable to those of all-electron GIMIC calculations because the computational time for preparing the density matrices in the AO basis for selected orbitals is much smaller than for the GIMIC calculation, whose computational time is independent of whether it is the total MICD susceptibility or orbital contributions to it.

The perturbed density matrices contain information about the response of all orbitals through the mixing of occupied and virtual orbitals, whereas in calculations of orbital contributions using the few-orbital model, only selected virtual orbitals are considered.24,25,31 The selection might lead to wrong conclusions because many small contributions from pairs of occupied and virtual orbitals may affect the interpretation. We can perform calculations using the few-orbital model by identifying the orbital response involving selected occupied and virtual orbitals, whereas the rest of the response terms are set to zero. We have not implemented the approximate few-orbital model because the computational costs are the same as when considering the response of all orbitals.

2.2 Group theoretical analysis

The external magnetic field reduces the symmetry of the molecular point group when symmetry operations of the molecular point group do not map the magnetic field onto itself.66,67 The resulting molecular point group in the presence of the external magnetic field is a subgroup of the molecular point group in the absence of the field and of C∞h,68 which is the point group of the external magnetic field. The order of the subgroup is usually larger than one, when the magnetic field coincides with a symmetry axis of the molecular point group.66 The point group of the molecular structure in the presence of an external magnetic field can conveniently be determined by using the flow chart reported in ref. 66.

Orbital contributions to the MICD susceptibility are not in general divergence free, whereas the orbital contributions of all orbitals in a given irreducible representation of the combined point group are divergence free. Thus, the sum of the orbital contributions to the MIRC of all orbitals in each irreducible representation is then charge conserving, which leads to an MIRC whose strength is independent of the orientation of the integration plane.

The point groups of the studied molecules in the presence of an external magnetic field along the main symmetry axis, which is perpendicular to the molecular ring, are C6h (benzene), C4h (planar cyclooctatetraene), S4 (bent cyclooctatetraene), C3h (borazine, cyclopropenium cation, and hexadehydro[12]annulene) , C2h (porphin, 1,4-cyclohexadiene, tetraoxa-isophlorin and cyclobutadiene), and Cs (1,5-dibora-2,4-diazabenzene).

2.3 Spatial averaging of the orbital contributions

The MICD (J(r)) is obtained by contracting the MICD susceptibility tensor with the direction vector of the external magnetic field. For a given direction of the external magnetic field, the circulation direction (tropicity) of the MIRC vortex can be in the classical diatropic (aromatic) direction, which is here the positive direction or in the negative paratropic (antiaromatic) direction. The tropicity can be assessed by visualizing the vector field. The profile of the MIRC and its strength (I) can be obtained by numerically integrating J(r) that passes through a surface (S) cutting the molecular ring.
 
image file: d5sc00627a-t6.tif(7)
where the dot represents the scalar product of J(r) with the normal of the infinitesimal cross-section area dS. The two-dimensional integration domain begins in the center of the MIRC vortex and extends to a sufficiently large distance from the molecule in the other spatial directions. Some inaccuracy is introduced when the stagnation line of the MIRC vortex and the edge of the integration plane do not coincide. However, the studied molecules have high symmetry implying that the magnetic field, the stagnation line of the MIRC vortex and the main symmetry axis coincide. The integration plane is usually positioned so that it cuts the ring in the middle of a chemical bond. However, since charge conservation is almost fulfilled when using GIAOs and standard basis-set sizes, I is only weakly dependent on the position of S.69 The tropicity of J(r) can be identified and separated into diatropic and paratropic contributions by following the vector field around the vortices of J(r) using e.g., the Runge–Kutta algorithm.70

Since contributions to the MIRC strength of individual orbitals (Ii) may depend on the position of S, unique orbital contributions to the MIRC strength can be obtained by calculating the average of the MIRC strength for various orientations of S around the ring. This can be done by rotating S around the center of the ring with one of its edges along the MIRC vortex as shown in Fig. 1A. The average MIRC strength 〈Ii〉 is then obtained as

 
image file: d5sc00627a-t7.tif(8)
which corresponds to numerical integration around a circle using the trapezoidal rule. We use an integration step of Δα = 1° to ensure high accuracy.


image file: d5sc00627a-f1.tif
Fig. 1 The orientations of the integration plane for calculating the angular dependence of the MIRC strength of benzene (A). The angular dependence of the MIRC strength of benzene for all orbitals in each irreducible representation of the C6h point group (B). The angular dependence of the MIRC strength of the core, valence σ, and π orbitals of benzene (C).

3 Computational details

The molecular structures were optimized with Turbomole71–73 at the density functional theory (DFT) level using the B3LYP functional,74,75 the empirical D4 dispersion correction76 and the def2-TZVP77 basis sets. The one-electron density matrix and the three perturbed density matrices were obtained by calculating the nuclear magnetic resonance (NMR) shielding constants at the DFT level78,79 using the B3LYP74,75 functional, which is a popular functional that has been used in many applications. We have also calculated NMR shielding constants with the strong-correlation local hybrid functional (scLH22t).80,81 The scLH22t functional was chosen because an accurate magnetizability of the ozone molecule is obtained with it. Further geometry optimizations and calculations of density matrices were performed for C3H3+ to investigate the convergence of the MIRC strength and orbital contributions to the MIRC strengths with respect to the basis-set size. The B3LYP functional and the def2-QZVP, cc-pVTZ, cc-pVQZ, cc-pV5Z*, and cc-pV6Z* basis sets were used.77,82–84 Basis functions with a larger angular momentum number than l = 4 were removed from the cc-pV5Z and cc-pV6Z basis set because the present implementation cannot handle h and higher angular momentum functions.

The MICD susceptibility tensor was calculated with the GIMIC program.85 The integration plane begins 8 bohr above the molecular plane and ends 8 bohr below it. In the other direction it begins at the symmetry axis and ends 15 bohr from the axis. The plane is divided into elements whose dimension is 0.05 × 0.05 bohr. We use 9 × 9 Gauss integration points in each element to ensure a high accuracy. The J(r) vector field was obtained by contracting the tensor with the direction of the external magnetic field. The external magnetic field was aligned along the main symmetry axis. The MIRC strengths were obtained by integrating the J(r) passing through a plane that cuts through the molecule.11–13 The molecular structures and orbitals are visualized with Chemcraft,86 and the J(r) pictures were made with ParaView.87 The plots were generated using Matplotlib, a Python library for data visualization.88

4 Results and discussion

4.1 Benzene

The orbital contributions to the MIRC strength calculated for the core, valence σ and π orbitals of benzene are given in Table 1. The total orbital contributions of all orbitals in each irreducible representation of the C6h point are divergence free in the limit of complete basis set when using the GIMIC method. The MICD susceptibility calculated with the GIMIC program is practically divergence free when using standard basis sets (def2-TZVP) as seen in Fig. 1B. The MIRC contributions of all orbitals in the irreducible representations are also given in Table 1.
Table 1 The average orbital contributions to the MIRC strength (〈Ii〉 in nA/T) and the total MIRC strengths of the core, σ, and π orbitals of benzene calculated at the B3LYP and scLH22t levels. The orbital contributions and the contributions of all orbitals in each irreducible representation of the C6h point group are compared to the orbital contributions of the isotropic shielding constant (σ(0) in ppm) previously calculated at the Hartree–Fock level in the center of the ring using the ipsocentric method.25 We also compare to orbital contributions to the σzz(0) value calculated at the DFT level using GIAOs89
Contribution
D6h C6h Type B3LYP scLH22t σ(0)25 σ zz(0)89
a The core contributions are not included.
1a1g 1ag Core 0.06 0.03 0.10
1e1u 1e1u Core 0.09 0.15 −0.64
1e2g 1e2g Core 1.19 1.22 0.73
1b2u 1bu Core 0.77 0.84 −0.06
2a1g 2ag σ 3.03 3.15 −0.40 14.4
2e1u 2e1u σ 4.88 5.26 0.09 11.4
2e2g 2e2g σ 4.35 4.56 8.94 8.8
3a1g 3ag σ 2.81 2.73 0.34 12.3
2b2u 2bu σ −0.51 −0.36 −2.25 −2.9
3e1u 3e1u σ −3.01 −3.39 −7.47 −14.0
1b1u 3bu σ 1.79 1.62 0.01 4.2
1a2u 1au π 4.12 4.14 2.92 12.8
3e2g 3e2g σ −15.16 −15.46 −11.34 −67.6
1e1g 1e1g π 7.53 7.59 18.90 23.4
Core 2.12 2.24 0.13 10.8
σ −1.82 −1.90 −12.06 −33.3
Core + valence σ 0.30 0.34 −11.94 −22.5
π 11.65 11.73 21.83 36.2
Ag 5.90 5.91 0.04 26.7a
Bu 2.05 2.10 −2.30 1.3a
E1u 1.96 2.02 −8.02 −2.6a
Au 4.12 4.14 2.92 12.8
E2g −9.62 −9.68 −1.67 −58.8a
E1g 7.53 7.59 18.90 23.4
Total 11.95 12.07 9.89 13.7


The MIRC strengths (in nA/T) calculated with GIMIC and the isotropic shielding constant in the center of the ring (σ(0) in ppm) calculated with the ipsocentric approach25 correlate suggesting that similar orbital contributions to the MICD susceptibility are obtained using GIMIC and the ipsocentric method. Orbital contributions to the MICD susceptibility and other magnetic properties of benzene have also been studied using the ipsocentric method,24,25 where they found that the HOMO contribution to J(r) of benzene is diatropic and mainly determines the MIRC strength.25 The GIMIC calculations yield a total π contribution of 11.65 nA/T, whereas the core + valence σ contribution is only 0.30 nA/T. These values are close to the ones of 11.7 nA/T and 1.1 nA/T, respectively, which were previously calculated at the Hartree–Fock level with the SYSMOIC program using the ipsocentric approach.48

The diatropic π contribution to the MIRC of benzene originates from the π orbitals in the Au and E1g irreducible representations of the C6h point group. The orbitals in the E2g irreducible representation sustain a strong paratropic MIRC of −9.62 nA/T, which is almost canceled by contributions of the valence σ and core orbitals in the other irreducible representations of the C6h point group. Since the separation into the core and valence σ contributions to the MIRC strength depends on the orientation and position of the integration plane, we compare average contributions that are obtained as described in Section 2.3. The average contribution of the core and valence σ orbitals are 2.12 nA/T and −1.82 nA/T, respectively. The contributions of the valence σ and core orbitals are tiny when the integration plane cuts the ring in the middle of a C–C bond as seen in Fig. 1C.

In Table 1, one can see that all occupied valence orbitals contribute significantly to the MIRC strength. The contribution of the core orbitals is diatropic, whereas the contributions of the valence σ orbitals are both positive and negative. The 〈Ii〉 values of the valence orbitals have in general the same sign as the corresponding orbital contributions to σ(0), since diatropic orbital contributions to the MIRC result in diamagnetic contribution to σ(0) and paratropic orbital contributions to the MIRC result in paramagnetic shielding contributions.49,50 Even though the orbital contributions to the MIRC strength and to σ(0) have the same trend, the relative size of the contributions significantly differs. The discrepancy is more severe when comparing the formally divergence-free contributions of all orbitals in each irreducible representation. The orbitals in the Ag irreducible representation sustain an MIRC of 5.90 nA/T at the B3LYP level, whereas the corresponding σ(0) value of 0.04 ppm is tiny.

For the Bu and E1u irreducible representations, the MIRC strength and the σ(0) contributions have opposite signs. The MIRC contribution of the orbitals in the E2g irreducible representation is strong and paratropic, whereas the σ(0) contribution is weakly paramagnetic. The π contribution to the MIRC strength of the only orbital in the Au irreducible representation is about 35% of the MIRC, which is of the same relative size as the Au contribution to σ(0). The MIRC contribution of the two π orbitals in the E1g irreducible representation is about 65% of the MIRC strength, whereas the corresponding contribution to σ(0) is a factor of two larger than the total σ(0) value. The three π orbitals (HOMO−4, HOMO−1 and HOMO) and their MIRC contributions as a function of the orientation angle of the integration plane are shown in Fig. 2. The MIRC strength of HOMO−4 with respect to the orientation of the integration plane is almost independent of the angle because HOMO−4 is the only orbital in the irreducible representation Au of the group of points C6h. The oscillating behavior of the contributions of the two π orbitals in E1g and that the sum of them is constant are also seen in Fig. 2. The small waves of the horizontal lines near the nuclei are due to a tiny J(r) leakage because finite basis sets have been used.11


image file: d5sc00627a-f2.tif
Fig. 2 The π orbitals of benzene (top), their MICD (middle) and the dependence of their contributions to the MIRC strength (in nA/T) as a function of the angle of the integration plane (bottom). The orbital pictures are made with ParaView.87

There is in general no clear correlation between the orbital contribution to the MIRC strength and the orbital contributions to the σzz(0) and σ(0) values of benzene. However, the MIRC strength of the π orbitals and the π-orbital contributions to the σzz(0) values correlate. In a study on 43 molecular rings, the best linear correlation was obtained for NICSπ,zz(z) (−σπ,zz(z)) values and MIRC strengths of the π orbitals when z is 1.25 Å.90

Even though the σ(0) contributions are calculated at the Hartree–Fock level using the ipsocentric approach,48 the comparison of the MIRC strengths and σ(0) values casts a shadow over the ability of NICS calculations to estimate orbital contributions to the MIRC strength. We discuss in the text only MIRC strengths calculated at the B3LYP level, since calculations on aromatic molecules using the scLH22t functional yield MIRC values that agree well with the ones obtained at the B3LYP level. See Table 1.

4.2 Borazine

The contributions to the MIRC strength of borazine of all orbitals of each irreducible representation of the C3h point group are given in Table 2, whereas the contributions of each occupied orbital are given in the ESI. The average MIRC contributions of the core, valence σ and π orbitals are also reported in Table 2. The orbital contribution to the MIRC strength of the HOMO is −0.47 nA/T. The HOMO is a doubly degenerate π orbital and the only orbital belonging to the E′′ irreducible representation. For benzene, the HOMO–LUMO transition is translationally allowed, resulting in a strong diatropic contribution to the MIRC of the HOMO, whereas the HOMO–LUMO transition of borazine is translationally and rotationally allowed,25 resulting in a weak MIRC contribution. The J(r) of the HOMO of borazine consists of local vortices circulating around the nitrogen atoms.91 The MIRC strength of the HOMO is independent of the orientation of the integration plane because it is the only occupied orbital belonging to the E′′ irreducible representation.
Table 2 The average contributions to the MIRC strength of borazine (〈Ii〉 in nA/T) of all orbitals in each irreducible representation of the C3h point group. The calculations were performed at the B3LYP and scLH22t levels. The contributions of the core, valence σ and π orbitals are also given. Orbital contributions are reported in the ESI
Contribution B3LYP scLH22t
A′ 7.14 7.06
E′ −6.94 −6.87
A′′ 3.52 3.50
E′′ −0.46 −0.63
Core 3.04 2.80
Valence σ −2.85 −2.61
Core + valence σ 0.20 0.19
π 3.06 2.87
Total 3.25 3.06


The MIRC contribution of HOMO−2 is 3.52 nA/T, which is the third π orbital. Borazine is weakly aromatic91 sustaining an MIRC of 3.25 nA/T, which is at the borderline between the MIRC strength of aromatic and nonaromatic molecules. Molecules with an absolute MIRC strength of <3 nA/T i.e., 25% of the benzene value can be considered nonaromatic. The average contribution of the core orbitals of 3.04 nA/T almost cancels the contribution of the valence σ orbitals of −2.85 nA/T. The contributions of the core + valence σ orbitals of 0.2 nA/T and of the π orbitals of 3.06 nA/T are close to the ones of 0.3 nA/T and 1.8 nA/T, respectively, which were obtained in SYSMOIC calculations at the Hartree–Fock level using the ipsocentric approach.48

4.3 1,5-Dibora−2,4-diazabenzene

1,5-Dibora-2,4-diazabenzene (C2B2N2H6) is isoelectronic with benzene and aromatic sustaining an MIRC strength of 7.67 nA/T, which is more than 60% of the one for benzene. Since it belongs to the C2v point group, Cs is its highest possible point group in the presence of an external magnetic field. The contributions to the MIRC strength of C2N2B2H6 of all orbitals in each irreducible representation of the Cs point group are given in Table 3, whereas the contributions of each occupied orbital are given in the ESI. The MICD susceptibility in the A′ (σ) and A′′ (π) irreducible representations is divergence free. The MIRC strength of the π orbitals is 7.59 nA/T. Thus, almost no MIRC is sustained in the σ orbitals. The average MIRC strength of the core orbitals of 2.67 nA/T is canceled by the contributions of the valence σ orbitals. The MIRC strengths of the core and valence σ orbitals depend on the orientation of the integration plane because the orbitals belong to the same irreducible representation of the Cs point group.
Table 3 The average contributions to the MIRC strength of C2B2N2H6 (〈Ii〉 in nA/T) of the orbitals in the two irreducible representations of the Cs point group. The calculations were performed at the B3LYP and scLH22t levels. The contributions of the core, valence σ and π orbitals are also given. Orbital contributions are reported in the ESI
Contribution B3LYP scLH22t
Core 2.67 2.52
Valence σ −2.59 −2.42
A′ (core + valence σ) 0.08 0.10
A′′ (π) 7.59 7.55
Total 7.67 7.65


4.4 Porphin

For porphin, the average contribution of the core orbitals of 11.22 nA/T is almost canceled by the contribution of the valence σ orbitals of −10.75 nA/T. The resulting contribution of the core + valence σ orbitals of 0.47 nA/T is much smaller than the average contribution of 27.2 nA/T of the π orbitals. The total MICD strength is 27.68 nA/T. The contributions to the MIRC strength of porphin of all orbitals in each irreducible representation of the C2h point group are given in Table 4. The contributions of the core, valence σ, and π orbitals are also reported.
Table 4 The average contributions to the MIRC strength (in nA/T) of porphin, tetraoxa-isophlorin, 1,4-cyclohexadiene and cyclobutadiene of the orbitals in each irreducible representation of the C2h point group. The calculations were performed at the B3LYP and scLH22t levels. The contributions of the core, valence σ and π orbitals are also given. Orbital contributions are reported in the ESI
Contribution Porphin Isophlorin 1,4-Cyclohexadiene Cyclobutadiene
B3LYP scLH22t B3LYP scLH22t B3LYP scLH22t B3LYP scLH22t
Ag 0.06 0.14 −0.26 −0.34 −4.08 −4.16 8.89 8.85
Bg 12.50 12.65 −76.80 −60.45 3.62 3.46 −19.88 −20.56
Au 14.70 14.90 13.62 12.95 −4.41 −4.13 3.84 3.85
Bu 0.41 0.48 −0.56 −0.59 4.23 4.26 −12.79 −12.91
Core 11.22 10.56 10.86 10.31 2.03 2.14 1.32 1.41
σ −10.75 −9.94 −11.69 −11.24 −1.88 −2.04 −5.22 −5.47
Core + valence σ 0.62 −0.83 −0.93 0.15 0.11 −3.90 −4.06
π 27.20 27.55 −63.17 −47.50 −0.79 −0.67 −16.05 −16.71
Total 27.68 28.17 −64.00 −48.44 −0.64 −0.56 −19.95 −20.76


About half of the π contribution to the MIRC strength (14.70 nA/T) originates from the orbitals in the Au irreducible representations of the C2h point group and the contribution of the Bg orbitals is 12.50 nA/T. The averaged contributions of the individual orbitals in the ESI show that nine π orbitals contribute more than 1 nA/T to the total MIRC. The individual MIRC contribution of these π orbitals is almost the same for seven of them (∼3.5 nA/T). The only exceptions are the contributions of HOMO of 2.20 nA/T and of HOMO−1 of 1.28 nA/T. The sum of their contributions is also ∼3.5 nA/T. The orbital contributions of four of the π orbitals are tiny, that is, in the range of ± 1nA/T. The contributions of HOMO−2 and HOMO−3 are −1.01 nA/T and −0.01 nA/T, respectively. A previous study of orbital contributions to the J(r) of porphin using the ipsocentric method suggested that four frontier π orbitals contribute to the MIRC strength,31 which does not agree with the results obtained in our study. The total MIRC strength is almost entirely determined by contributions of the π orbitals. However, many energetically low-lying π orbitals contribute to the total MIRC.

The diatropic contribution to the MIRC of the HOMO can be explained by the translationally allowed transition from HOMO (b1u) to LUMO (b3g), whereas the diatropic contributions of the energetically low-lying π (b1u) orbitals using selection rules may involve transitions to many virtual orbitals belonging to the same irreducible representation as the LUMO.

4.5 Nonaromatic molecule

The absence of continuous electron delocalization is generally the reason why molecular rings are nonaromatic. Here, we study 1,4-cyclohexadiene as an example of a nonaromatic molecular ring. The electron delocalization of 1,4-cyclohexadiene ring is disrupted by the out-of-plane hydrogen atoms at the two sp3-hybridized carbon atoms. The MIRC strengths of the σ and π orbitals of 1,4-cyclohexadiene can be separated because its molecular structure belongs to the D2h point group. In the presence of a magnetic field along the main symmetry axis, the core and valence σ orbitals belong to the Ag and Bu irreducible representations of the C2h point group. The sum of their contributions to the MIRC strength is 0.15 nA/T. An MIRC strength of the π orbitals of −0.79 nA/T is obtained as the sum of the contributions of the orbitals in the Au and Bg irreducible representations. Since 1,4-cyclohexadiene does not sustain any strong MIRC in the σ nor in the π orbitals it is nonaromatic, even though the transitions between occupied and unoccupied frontier orbitals are translationally allowed. The MIRC strengths of the orbitals in each irreducible representation are given in Table 4 where we also report the average contributions of the core, valence σ and π orbitals. Orbital contributions are reported in the ESI.

4.6 Antiaromatic molecules

Tetraoxa-isophlorin is a strongly antiaromatic porphyrinoid. Previous studies showed that it has a large magnetic transition moment from the electronic ground state to the lowest excited state resulting in a significant paramagnetic contribution to the magnetizability.92,93 Hexadehydro[12]annulene is also an antiaromatic molecule sustaining a strong paratropic MIRC.94–96 The average contributions to the MIRC strength of tetraoxa-isophlorin for all orbitals in the irreducible representations of the C2h point group are reported in Table 4. The MIRC strength of hexadehydro[12]annulene for all orbitals in the irreducible representations of the C3h point group are given in Table 5. Orbital contributions are reported in the ESI.
Table 5 The average contributions to the MIRC strength (in nA/T) of hexadehydro[12]annulene (HDH[12]A), C3H3+ and C3H6 of the orbitals in each irreducible representation of the C3h point group. The calculations were performed at the B3LYP and scLH22t levels. The contributions of the core, valence σ, and π orbitals are also given. Orbital contributions are reported in the ESI
HDH[12]A C3H3+ C3H6
B3LYP scLH22t B3LYP scLH22t B3LYP scLH22t
A′ −0.01 0.14 3.33 3.35 3.02 3.07
E′ −0.02 −0.26 3.14 3.04 6.42 6.42
A′′ −36.86 −31.05 3.90 3.92 2.84 2.81
E′′ 12.77 12.20 −2.36 −2.30
Core 11.18 9.26 1.07 1.30 0.85 0.79
σ −11.22 −9.40 5.40 5.09 8.59 8.70
Core + σ −0.04 −0.14 6.47 6.40 9.44 9.49
π −24.10 −18.86 3.90 3.92 0.49 0.50
Total −24.13 −18.99 10.37 10.32 9.92 9.99


The average contributions to the MIRC of the core + valence σ orbitals are nearly zero. The total MIRC strength originates mainly from the π orbitals, which belong to the Bg and Au irreducible representations for tetraoxa-isophlorin, and to the A′′ and E′′ irreducible representations for hexadehydro[12]annulene,.

The contributions of the individual orbitals show that the MIRC strength of nine π orbitals of tetraoxa-isophlorin is diatropic. Each of them contributes more than 1.5 nA/T to the total MIRC strength. Their total diatropic MIRC strength is 23.81 nA/T. The diatropic contribution of four other π orbitals is 2.22 nA/T. The total contribution of all π orbitals except one is 26.03 nA/T, which is almost as large as for porphin. The large paratropic contribution of −89.97 nA/T of the HOMO calculated at the B3LYP level makes tetraoxa-isophlorin antiaromatic. Calculations using the scLH22t functional yield a HOMO contribution of −73.54 nA/T showing the importance of significant long-ranged Hartree–Fock exchange in the functional to properly describe magnetic response of strongly antiaromatic molecules, whereas for aromatic molecules the MIRC strengths obtained with the two functionals are almost equal.

Five of the six π orbitals of hexadehydro[12]annulene sustain a total diatropic MIRC of 16.67 nA/T, whereas the MIRC contribution of the HOMO is −40.76 nA/T making the molecule antiaromatic. The strong paratropic MIRC contribution of the HOMO can be explained by using symmetry selection rules. The image file: d5sc00627a-t8.tif transition of hexadehydro[12]annulene is purely rotationally allowed as also for the HOMO (b2g)–LUMO (b3g) transition of tetraoxa-isophlorin.

The planar structure of cyclooctatetraene (COT), which belongs to the D4h point group, is reduced to the C4h point group in the presence of an external magnetic field along the main symmetry axis. The average contributions to the MIRC strength of COT of all orbitals in each irreducible representation are given in Table 6. The average contributions of the core, valence σ, and π orbitals are also given. Orbital contributions are reported in the ESI.

Table 6 The average contributions to the MIRC strength (〈Ii〉 in nA/T) of planar cyclooctatraene for all orbitals in the irreducible representations of the C4h point group and of the bent cyclooctatraene for all orbitals in the irreducible representations of the S4 point group are given. The calculations were performed at the B3LYP and scLH22t levels. Average contributions of the core, valence σ and π orbitals are also given. Orbital contributions are reported in the ESI
Contribution (C4h) B3LYP scLH22t
Ag 6.45 6.46
Bg −5.65 −5.83
Eg 7.48 7.30
Au 4.10 4.11
Bu −51.68 −46.47
Eu −1.11 −1.17
Core 4.06 3.74
σ −4.37 −4.28
Core + valence σ −0.30 −0.54
π −40.10 −35.06
Total −40.41 −35.60

Contribution (S4) B3LYP scLH22t
A −7.47 −6.37
B −2.02 −2.04
E 6.67 6.46
Total −2.82 −1.96


Planar COT is antiaromatic and characterized by a strong paratropic HOMO contribution of −51.68 nA/T to the MIRC strength, which can be explained by the rotationally allowed HOMO (b2u)–LUMO (b1u) transition.97 The three other π orbitals sustain diatropic contributions to the MIRC strength. The sum of their contributions is 11.58 nA/T, which is close to the MIRC strength of benzene. The core + valence σ contribution of only −0.30 nA/T can be ignored. Planar COT is strongly antiaromatic sustaining a paratropic MIRC of −40.41 nA/T at the B3LYP level. The MIRC strength of −35.60 nA/T calculated with the scLH22t functional is slightly smaller.

The MIRC contributions of the core + valence σ and π orbitals of planar cyclooctatetraene have been previously calculated at the Hartree–Fock level with the SYSMOIC program. The obtained core + valence σ contribution of −1.1 and the π contribution of −18.0 nA/T also show that the π contribution dominates making the molecule antiaromatic.48

Bent cyclooctatetraene belonging to the D2d point group is energetically more stable than the planar structure.98 The average MIRC contribution of the HOMO of −12.77 nA/T for the bent cyclooctatetraene structure shows that bending the molecule significantly reduces the HOMO contribution to the MIRC strength. The total MIRC strength of the bent COT is only −2.82 nA/T suggesting that it is nonaromatic. The paratropic MIRC contribution of the HOMO is the largest one. Many other orbitals contribute significantly. The average contribution of the core orbitals is 2.77 nA/T. Orbital MIRC strengths are reported in the ESI. The divergence-free contributions of all orbitals in each irreducible representation of the S4 point group are given in Table 6. The spatial distributions of the J(r) of the HOMO of the planar and bent cyclooctatetraene are shown in the ESI.

4.7 Strain

Molecular properties of planar ring-shaped molecules, whose ∠CCC angles significantly smaller than the standard ∠CCC angle of sp2 hybridized hydrocarbons, differ from those of ordinary molecular rings due to the ring strain.99 The strain leads to delocalization of the electrons in the σ orbitals and often to σ aromaticity. However, all strained molecular rings do not exhibit σ aromaticity.

We study orbital contributions to the J(r) of the planar cyclopropenium cation (C3H3+) and cyclobutadiene (C4H4) because their magnetic properties have significant contributions of the σ orbitals. Previous studies of the delocalization energy of the C3H3+ and C4H4 suggested that C3H3+ is σ + π aromatic and that C4H4 is σ + π antiaromatic.100,101 Calculations of the orbital contributions to the MIRC strength of C4H4 show that 20% of the MIRC strength is sustained by the core + valence σ orbitals. Since the core contribution is diatropic, the paratropic MIRC contribution of valence σ orbitals is even 25% of the total MIRC strength. The MIRC strength of all orbitals of C4H4 belonging to the irreducible representations of the C2h point group are given in Table 4. The total MIRC of the antiaromatic C4H4 is −19.95 nA/T at the B3LYP level and −20.76 nA/T at the scLH22t level. The HOMO (b2g) → LUMO (b3g) and the HOMO−1 (b2u) → LUMO+1 (b3u) transitions of C4H4 are rotationally allowed leading to σ + π antiaromaticity. Thus, a significant paratropic MIRC contribution of −12.25 nA/T arises from the HOMO−1 σ orbital due to the strain.

The core + valence σ contribution to the MIRC strength of C3H3+ is 6.47 nA/T, which is even larger than the π contribution of 3.90 nA/T. Since the average core contribution is 1.07 nA/T, the MIRC is dominated by the valence σ contribution. The MIRC contributions of C3H3+ for all orbitals of the irreducible representations of the C3h point group are given in Table 5. The total MIRC strength of 10.37 nA/T is almost as larger as the MIRC strength of benzene. All orbitals except the HOMO sustain a diatropic MIRC. The MIRC of the HOMO is unexpectedly paratropic with a MIRC strength of −3.49 nA/T because the HOMO (e′) to LUMO+1 (e′) transition is rotationally allowed.

Previous computational studies on cyclopropane54,56 (C3H6) suggested that it is not σ aromatic, whereas other MICD susceptibility studies on C3H6 yielded a strong MIRC in the σ orbitals and a weak MIRC contribution in the π orbitals.48,55 The present calculations of the orbital contributions to the MIRC of C3H6 confirm that the MIRC is mainly sustained in the σ orbitals as shown in Table 5. Even though the total MIRC of C3H3+ (10.37 nA/T) and of C3H6 (9.92 nA/T) are almost equal, the contribution of the core + valence σ orbitals to the MIRC strength of C3H6 of 9.44 nA/T dominates. It is an order of magnitude larger than the π-orbital contribution of 0.49 nA/T. The contribution of the HOMO to the MIRC strength of C3H6 is 1.57 nA/T. The HOMO (e′) → LUMO+2 (e′) transition is translationally allowed, whereas the HOMO (e′) → LUMO image file: d5sc00627a-t9.tif transition is both translationally and rotationally forbidden. The contributions to the MIRC strength of core + valence σ and π orbitals obtained in this study are also in good agreement with the previously reported ones of 9.7 nA/T and 0.5 nA/T, respectively, which were obtained using the SYSMOIC program.48

The extended Hückel rule for σ orbitals51 suggests that C4H4 is σ antiaromatic and C3H3+ is σ aromatic, since C4H4 has 4n core + valence σ electrons and C3H3+ has 4n + 2 core + valence σ electrons.

The angular dependence of the MIRC contributions of the core, valence σ, and π orbitals of C3H3+ and C4H4 is shown in Fig. 3. The core contribution is diatropic for all directions of the integration plane and vanishes when the plane passes through the middle of the C–C bond. The core contribution increases in the vicinity of the nuclei, where the valence σ contribution to the MIRC strength decreases. The sum of the two contributions is independent of the orientation angle of the integration plane. The core + valence σ contribution is tiny for most of the studied molecules, whereas for C3H3+ and C4H4 the sum of the two contributions is independent of the angle of the integration plane and it is of the same size as the π contribution to the MIRC strength.


image file: d5sc00627a-f3.tif
Fig. 3 The angular dependence of the contributions of the core, valence σ, and π orbitals of C3H3+ (top) and of cyclobutadiene (bottom). The MIRC strengths of C3H3+ were calculated using the def2-QZVP basis sets to show that the divergence is smaller when using larger basis sets. The positions of the nuclei are indicated by the vertical dashed lines.

4.8 Basis-set effects

The MICD susceptibility is formally divergence free. The charge conservation condition of the orbital contributions to J(r) is fulfilled when considering all orbitals in a given irreducible representation of the point group of the molecular structure in the presence of an external magnetic field. However, the condition is not fulfilled in every point in space when finite basis sets are used in the GIMIC calculation. Even though one can correct the divergence of the MICD susceptibility,45 it is generally not necessary because the divergence of the MICD susceptibilities obtained in GIMIC calculations is small and does not affect most interpretations. The MIRC strengths were calculated as a function of the orientation angle of the integration plane as discussed in Section 2.3. The difference between the maximum and minimum values of the total MIRC strength of C3H3+ and the orbital contributions to the MIRC strength of all orbitals in each irreducible representation of the C3h point group are reported in Table 7. The angular dependence of the MIRC strength decreases when using systematically larger basis sets from cc-pVTZ to cc-pV6Z, while the mean values are almost converged with the triple-ζ quality basis sets. The orbital contributions as a function of the orientation angle of the integration plane are not constant even with the largest basis set suggesting that larger basis sets than cc-pV6Z must be employed to obtain a MICD susceptibility that is divergence free in all points of space. The cc-pV6Z basis sets are large but they are probably not flexible enough near the nuclei.
Table 7 The mean (nA/T) and the range (Δ, nA/T) of the MIRC strength, calculated as the difference between the minimum and maximum values, are obtained from the angular dependence of the MIRC contributions of C3H3+ for all orbitals in each irreducible representation of the C3h point group. The MIRC strengths are calculated at the B3LYP level using the cc-pVTZ, cc-pVQZ, cc-pV5Z and cc-pV6Z basis sets
cc-pVTZ cc-pVQZ cc-pV5Z cc-pV6Z
Mean Δ Mean Δ Mean Δ Mean Δ
A′ 3.33 0.29 3.32 0.29 3.32 0.21 3.32 0.16
E′ 3.15 0.92 3.21 0.47 3.22 0.27 3.22 0.17
A′′ 3.90 0.10 3.91 0.05 3.91 0.03 3.91 0.01
Total 10.37 0.81 10.44 0.27 10.44 0.11 10.45 0.04


5 Summary and conclusions

We have developed and implemented a novel approach in the GIMIC program for calculating contributions of occupied orbitals to the magnetically induced current density (MICD) susceptibility tensor. The method has been applied to various aromatic, nonaromatic and antiaromatic monocyclic molecules. The orbital contributions to the MIRC strength of benzene, borazine, 1,4-cyclohexadiene, 1,5-dibora-2,4-diazabenzene, the cyclopropenium cation, cyclobutadiene, planar cyclooctatetraene, bent cyclooctatetraene, hexadehydro[12]annulene, porphin and tetraoxa-isophlorin have been thoroughly investigated. The strength of the magnetically induced ring current (MIRC) has been obtained by numerically integrating the MICD susceptibility passing through an integration plane that crosses the molecular ring using accurate integration grids. We have used various orientations of the integration plane to assess how well the charge-conservation condition of the MICD susceptibility is fulfilled. The present results agree with previously reported data that have been calculated for some of the studied molecules using alternative methods. The present work is more comprehensive presenting for the first time detailed investigations of the orbital contributions to the MICD susceptibility and the MIRC of the eleven studied molecules.

The orbital contributions to the MICD susceptibility are not divergence free implying that the orbital contributions to the MIRC strength vary around the molecular ring. Unique average orbital contributions to the MIRC strength were obtained by calculating the MIRC strength for various orientations of the integration plane. Average MIRC strengths were obtained by numerically integrating over the orientation angle around the whole ring. The total MIRC and the MIRC contributions of all orbitals belonging to a given irreducible representation of the molecular point group in the presence of an external magnetic field are formally divergence free. The charge leakage of the MICD susceptibility calculated using the GIMIC method is tiny and can in most applications be ignored. The leakage decreases when increasing the size of the basis set and vanishes in the complete basis-set limit.

For planar antiaromatic molecules, the MIRC contribution of the HOMO is strongly paratropic and much larger than the contributions of the other orbitals, which is not the case for aromatic molecules. The total MIRC of aromatic molecules consists of orbital contributions of many orbitals. The orbital contributions to the MIRC strengths of C4H4, C3H3+ and C3H6 showed that the ring strain leads to a strong MIRC strength in an energetically high-lying σ orbital. The MIRC of the σ orbital of C4H4 is paratropic and for C3H3+ and C3H6 it is diatropic.

The presented algorithm works well for molecules whose ground-state wave function is dominated by a single Slater determinant. The size of the molecules that can be studied is limited only by the calculations of the NMR shielding constants. Since molecular magnetic properties can be obtained with the GIMIC program by numerical integration of the MICD susceptibility multiplied by the vector potential of the magnetic property, the present work also opens new possibilities for investigating orbital contributions to molecular magnetic properties.

Data availability

The data supporting this article have been included as part of the ESI. Density functional theory calculations were performed with Turbomole version 7.8. The Turbomole webpage is https://www.turbomole.org/. GIMIC, version 2.0 can be freely downloaded from https://github.com/qmcurrents/gimic and https://zenodo.org/record/8180435. The new features of the GIMIC program are not yet available to everyone. ParaView can be downloaded from https://www.paraview.org/.

Author contributions

RTN implemented the methods, performed the calculations and prepared the pictures. DS proposed the project and developed the theory together with RTN. RTN and DS wrote the first version of the manuscript. All authors contributed to writing the article.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work has also been supported by the Academy of Finland through projects 340582 and 340583. The CSC – the Finnish IT Center for Science is acknowledged for computer time.

References

  1. R. Hoffmann, Am. Sci., 2015, 103, 18–22 Search PubMed .
  2. G. Merino, M. Solà, I. Fernández, C. Foroutan-Nejad, P. Lazzeretti, G. Frenking, H. L. Anderson, D. Sundholm, F. P. Cossío, M. A. Petrukhina, J. Wu, J. I. Wu and A. Restrepo, Chem. Sci., 2023, 14, 5569–5576 RSC .
  3. V. I. Minkin, Pure Appl. Chem., 1999, 71, 1919–1981 CAS .
  4. P. von Ragué Schleyer, C. Maerker, A. Dransfeld, H. Jiao and N. J. R. van Eikema Hommes, J. Am. Chem. Soc., 1996, 118, 6317–6318 Search PubMed .
  5. J. A. Pople, Mol. Phys., 1958, 1, 175–180 CrossRef CAS .
  6. R. McWeeny, Mol. Phys., 1958, 1, 311–321 CAS .
  7. T. A. Keith and R. F. W. Bader, Chem. Phys. Lett., 1993, 210, 223–231 CAS .
  8. P. Lazzeretti, M. Malagoli and R. Zanasi, Chem. Phys. Lett., 1994, 220, 299–304 Search PubMed .
  9. R. Zanasi, P. Lazzeretti, M. Malagoli and F. Piccinini, J. Chem. Phys., 1995, 102, 7150–7157 CAS .
  10. P. Lazzeretti, Prog. Nucl. Magn. Reson. Spectrosc., 2000, 36, 1–88 CAS .
  11. J. Jusélius, D. Sundholm and J. Gauss, J. Chem. Phys., 2004, 121, 3952–3963 Search PubMed .
  12. H. Fliegl, S. Taubert, O. Lehtonen and D. Sundholm, Phys. Chem. Chem. Phys., 2011, 13, 20500–20518 RSC .
  13. D. Sundholm, H. Fliegl and R. J. F. Berger, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2016, 6, 639–678 CAS .
  14. G. Monaco, F. F. Summa and R. Zanasi, J. Chem. Inf. Model., 2021, 61, 270–283 CrossRef CAS PubMed .
  15. M. Dimitrova and D. Sundholm, in Aromaticity: Modern Computational Methods and Applications, ed. I. Fernández López, Elsevier, 2021, pp. 155–194 Search PubMed .
  16. P. von Ragué Schleyer, M. Manoharan, Z.-X. Wang, B. Kiran, H. Jiao, R. Puchta and N. J. R. van Eikema Hommes, Org. Lett., 2001, 3, 2465–2468 CrossRef PubMed .
  17. T. Heine, P. v. R. Schleyer, C. Corminboeuf, G. Seifert, R. Reviakine and J. Weber, J. Phys. Chem. A, 2003, 107, 6470–6475 CrossRef CAS .
  18. T. Heine, R. Islas and G. Merino, J. Comput. Chem., 2007, 28, 302–309 CrossRef CAS PubMed .
  19. A. C. Tsipis, I. G. Depastas and C. A. Tsipis, Symmetry, 2010, 2, 284–319 CrossRef CAS .
  20. R. Islas, T. Heine and G. Merino, Acc. Chem. Res., 2012, 45, 215–228 CAS .
  21. N. D. Charistos, A. G. Papadopoulos, T. A. Nikopoulos, A. Muñoz-Castro and M. P. Sigalas, J. Comput. Chem., 2017, 38, 2594–2604 CAS .
  22. R. Gershoni-Poranne and A. Stanger, in Aromaticity: Modern Computational Methods and Applications, ed. I. Fernández López, Elsevier, 2021, ch. 4, pp. 99–154 Search PubMed .
  23. I. Pérez-Juste, M. Mandado and L. Carballeira, Chem. Phys. Lett., 2010, 491, 224–229 Search PubMed .
  24. E. Steiner and P. W. Fowler, J. Phys. Chem. A, 2001, 105, 9553–9562 CAS .
  25. E. Steiner and P. W. Fowler, Phys. Chem. Chem. Phys., 2004, 6, 261–272 CAS .
  26. E. Steiner and P. W. Fowler, Chem. Commun., 2001, 2220–2221 CAS .
  27. R. J. F. Berger and A. Viel, Z. Naturforsch., B:J. Chem. Sci., 2020, 75, 327–339 CAS .
  28. D. Sundholm, M. Dimitrova and R. J. Berger, Chem. Commun., 2021, 57, 12362–12378 CAS .
  29. H. Fliegl, D. Sundholm, S. Taubert, J. Jusélius and W. Klopper, J. Phys. Chem. A, 2009, 113, 8668–8676 CAS .
  30. A. Soncini and P. W. Fowler, Chem. Phys. Lett., 2004, 396, 174–181 CAS .
  31. E. Steiner and P. W. Fowler, ChemPhysChem, 2002, 3, 114–116 CAS .
  32. F. London, J. Phys. Radium, 1937, 8, 397–409 CAS .
  33. R. Ditchfield, Mol. Phys., 1974, 27, 789–807 CAS .
  34. J. R. Cheeseman, G. W. Trucks, T. A. Keith and M. J. Frisch, J. Chem. Phys., 1996, 104, 5497–5509 CAS .
  35. T. Helgaker, M. Jaszuński and K. Ruud, Chem. Rev., 1999, 99, 293–352 CrossRef CAS PubMed .
  36. L. N. Wirz, M. Dimitrova, H. Fliegl and D. Sundholm, J. Phys. Chem. Lett., 2018, 9, 1627–1632 CrossRef CAS PubMed .
  37. H. Fliegl, D. Sundholm and F. Pichierri, Phys. Chem. Chem. Phys., 2011, 13, 20659–20665 RSC .
  38. K. Reiter, F. Weigend, L. N. Wirz, M. Dimitrova and D. Sundholm, J. Phys. Chem. C, 2019, 123, 15354–15365 CrossRef CAS .
  39. R. R. Valiev, T. Kurten, L. I. Valiulina, S. Y. Ketkov, V. N. Cherepanov, M. Dimitrova and D. Sundholm, Phys. Chem. Chem. Phys., 2022, 24, 1666–1674 RSC .
  40. D. Blasco and D. Sundholm, Dalton Trans., 2024, 53, 10150–10158 RSC .
  41. A. Rabe, Q. Wang and D. Sundholm, Dalton Trans., 2024, 53, 10938–10946 RSC .
  42. Y. C. Lin, J. Jusélius, D. Sundholm and J. Gauss, J. Chem. Phys., 2005, 122, 214308 CrossRef PubMed .
  43. M. Orozco-Ic, L. Soriano-Agueda, D. Sundholm, E. Matito and G. Merino, Chem. Sci., 2024, 15, 12906–12921 RSC .
  44. R. J. F. Berger, M. Dimitrova, R. T. Nasibullin, R. R. Valiev and D. Sundholm, Phys. Chem. Chem. Phys., 2022, 24, 624–628 RSC .
  45. G. Monaco, F. F. Summa, R. Zanasi and R. J. F. Berger, J. Chem. Phys., 2024, 161, 194105 CrossRef CAS PubMed .
  46. J. A. Bohmann, F. Weinhold and T. C. Farrar, J. Chem. Phys., 1997, 107, 1173–1184 CrossRef CAS .
  47. M. Orozco-Ic, N. D. Charistos, A. Muñoz-Castro, R. Islas, D. Sundholm and G. Merino, Phys. Chem. Chem. Phys., 2022, 24, 12158–12166 RSC .
  48. G. Monaco, R. Zanasi, S. Pelloni and P. Lazzeretti, J. Chem. Theory Comput., 2010, 6, 3343–3351 CrossRef CAS PubMed .
  49. R. K. Jinger, H. Fliegl, R. Bast, M. Dimitrova, S. Lehtola and D. Sundholm, J. Phys. Chem. A, 2021, 125, 1778–1786 CAS .
  50. H. Fliegl, M. Dimitrova, R. J. F. Berger and D. Sundholm, Chemistry, 2021, 3, 1005–1021 CAS .
  51. M. J. S. Dewar, Bull. Soc. Chim. Belg., 1979, 88, 957–967 CAS .
  52. M. J. S. Dewar, J. Am. Chem. Soc., 1984, 106, 669–682 CAS .
  53. D. Moran, M. Manoharan, T. Heine and P. V. R. Schleyer, Org. Lett., 2003, 5, 23–26 CAS .
  54. S. Pelloni, P. Lazzeretti and R. Zanasi, J. Phys. Chem. A, 2007, 111, 8163–8169 CrossRef CAS .
  55. P. W. Fowler, J. Baker and M. Lillington, Theor. Chem. Acc., 2007, 118, 123–127 Search PubMed .
  56. W. Wu, B. Ma, J. I-Chia Wu, P. von Ragué Schleyer and Y. Mo, Chem. - Eur. J., 2009, 15, 9730–9736 CrossRef CAS PubMed .
  57. R. Carion, B. Champagne, G. Monaco, R. Zanasi, S. Pelloni and P. Lazzeretti, J. Chem. Theory Comput., 2010, 6, 2002–2018 CrossRef CAS PubMed .
  58. Y. Hao, J. Wu and J. Zhu, Chem. - Eur. J., 2015, 21, 18805–18810 Search PubMed .
  59. M. Orozco-Ic and D. Sundholm, Phys. Chem. Chem. Phys., 2023, 25, 12777–12782 RSC .
  60. P. Pulay, in Advances in Chemical Physics: Ab Initio Methods in Quantum Chemistry Part 2, ed. P. K. Lawley, John Wiley & Sons Ltd, Singapore, 1987, pp. 241–286 Search PubMed .
  61. M. Häser, R. Ahlrichs, H. P. Baron, P. Weis and H. Horn, Theor. Chim. Acta, 1992, 83, 455–470 Search PubMed .
  62. K. Wolinski, J. F. Hinton and P. Pulay, J. Am. Chem. Soc., 1990, 112, 8251–8260 CrossRef CAS .
  63. J. Gauss, Chem. Phys. Lett., 1992, 191, 614–620 CrossRef CAS .
  64. J. Olsen, K. L. Bak, K. Ruud, T. Helgaker and P. Jørgensen, Theor. Chim. Acta, 1995, 90, 421–439 CrossRef CAS .
  65. J. Gerratt and I. M. Mills, J. Chem. Phys., 1968, 49, 1730–1739 CrossRef CAS .
  66. A. Pausch, M. Gebele and W. Klopper, J. Chem. Phys., 2021, 155, 201101 CAS .
  67. S. Pelloni and P. Lazzeretti, Int. J. Quantum Chem., 2011, 111, 356–367 CrossRef CAS .
  68. A. J. Ceulemans, Group theory applied to chemistry, Springer, 2nd edn, 2024 Search PubMed .
  69. H. Fliegl, J. Jusélius and D. Sundholm, J. Phys. Chem. A, 2016, 120, 5658–5664 CrossRef CAS .
  70. Q. Wang, J. Pyykkö, M. Dimitrova, S. Taubert and D. Sundholm, Phys. Chem. Chem. Phys., 2023, 25, 12469–12478 RSC .
  71. R. Ahlrichs, M. Bär, H. S. Plitt and H. Schnöckel, Chem. Phys. Lett., 1989, 161, 179 CrossRef CAS .
  72. S. G. Balasubramani, G. P. Chen, S. Coriani, M. Diedenhofen, M. S. Frank, Y. J. Franzke, F. Furche, R. Grotjahn, M. E. Harding, C. Hättig, A. Hellweg, B. Helmich-Paris, C. Holzer, U. Huniar, M. Kaupp, A. Marefat Khah, S. Karbalaei Khani, T. Müller, F. Mack, B. D. Nguyen, S. M. Parker, E. Perlt, D. Rappoport, K. Reiter, S. Roy, M. Rückert, G. Schmitz, M. Sierka, E. Tapavicza, D. P. Tew, C. van Wüllen, V. K. Voora, F. Weigend, A. Wodyński and J. M. Yu, J. Chem. Phys., 2020, 152, 184107 Search PubMed .
  73. Y. J. Franzke, C. Holzer, J. H. Andersen, T. Begušić, F. Bruder, S. Coriani, F. Della Sala, E. Fabiano, D. A. Fedotov, S. Fürst, S. Gillhuber, R. Grotjahn, M. Kaupp, M. Kehry, M. Krstić, F. Mack, S. Majumdar, B. D. Nguyen, S. M. Parker, F. Pauly, A. Pausch, E. Perlt, G. S. Phun, A. Rajabi, D. Rappoport, B. Samal, T. Schrader, M. Sharma, E. Tapavicza, R. S. Treß, V. Voora, A. Wodyński, J. M. Yu, B. Zerulla, F. Furche, C. Hättig, M. Sierka, D. P. Tew and F. Weigend, J. Chem. Theory Comput., 2023, 19, 6859–6890 CrossRef CAS PubMed .
  74. A. D. Becke, J. Chem. Phys., 1993, 98, 5648–5652 Search PubMed .
  75. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B, 1988, 37, 785–789 CAS .
  76. E. Caldeweyher, C. Bannwarth and S. Grimme, J. Chem. Phys., 2017, 147, 034112 Search PubMed .
  77. F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305 RSC .
  78. K. Reiter, F. Mack and F. Weigend, J. Chem. Theory Comput., 2018, 14, 191–197 Search PubMed .
  79. Y. J. Franzke, J. Chem. Theory Comput., 2023, 19, 2010–2028 CAS .
  80. C. J. Schattenberg, A. Wodyński, H. Åström, D. Sundholm, M. Kaupp and S. Lehtola, J. Phys. Chem. A, 2023, 127, 10896–10907 Search PubMed .
  81. C. J. Schattenberg and M. Kaupp, J. Phys. Chem. A, 2024, 128, 2253–2271 Search PubMed .
  82. J. Dunning and H. Thom, J. Chem. Phys., 1989, 90, 1007–1023 Search PubMed .
  83. D. E. Woon, J. Dunning and H. Thom, J. Chem. Phys., 1994, 100, 2975–2988 Search PubMed .
  84. A. Wilson, T. van Mourik and T. H. Dunning Jr, J. Mol. Struct., 1997, 388, 339–349 CrossRef .
  85. GIMIC, version 2.0, a current density program, can be freely downloaded from https://github.com/qmcurrents/gimic and https://zenodo.org/record/8180435 (assessed 1.1.2025), 2019.
  86. Chemcraft, Chemcraft - graphical software for visualization of quantum chemistry computations, https://www.chemcraftprog.com, assessed: 23.01.2025.
  87. J. Ahrens, B. Geveci and C. Law, ParaView: An End-User Tool for Large Data Visualization, Visualization Handbook, Elsevier, 2005, ISBN-13: 978-0123875822, see also: http://www.paraview.org (assessed 18.1.2025) Search PubMed .
  88. J. D. Hunter, Comput. Sci. Eng., 2007, 9, 90–95 Search PubMed .
  89. N. D. Charistos, A. G. Papadopoulos and M. P. Sigalas, J. Phys. Chem. A, 2014, 118, 1113–1122 CrossRef CAS PubMed .
  90. S. Radenković and S. Ðorđević, Phys. Chem. Chem. Phys., 2021, 23, 11240–11250 RSC .
  91. R. Báez-Grez and R. Pino-Rios, RSC Adv., 2022, 12, 7906–7910 RSC .
  92. R. R. Valiev, H. Fliegl and D. Sundholm, Phys. Chem. Chem. Phys., 2017, 19, 25979–25988 RSC .
  93. R. R. Valiev, G. V. Baryshnikov, R. T. Nasibullin, D. Sundholm and H. Ågren, J. Phys. Chem. C, 2020, 124, 21027–21035 CrossRef CAS .
  94. J. A. Pople and K. G. Untch, J. Am. Chem. Soc., 1966, 88, 4811–4815 CrossRef CAS .
  95. J. Jusélius and D. Sundholm, Phys. Chem. Chem. Phys., 2001, 3, 2433–2437 RSC .
  96. J. Jusélius and D. Sundholm, Phys. Chem. Chem. Phys., 2008, 10, 6630–6634 RSC .
  97. P. Fowler and E. Steiner, Chem. Phys. Lett., 2002, 364, 259–266 CrossRef CAS .
  98. R. W. Havenith, P. W. Fowler and L. W. Jenneskens, Org. Lett., 2006, 8, 1255–1258 CrossRef CAS PubMed .
  99. R. Huisgen, Angew. Chem., Int. Ed., 1986, 25, 297–311 CrossRef .
  100. D. Cremer and J. Gauss, J. Am. Chem. Soc., 1986, 108, 7476–7477 Search PubMed .
  101. G. Hohlneicher, L. Packschies and J. Weber, Phys. Chem. Chem. Phys., 2007, 9, 2517–2530 RSC .

Footnote

Electronic supplementary information (ESI) available: Visualization of magnetically induced current densities, average orbital contributions to the magnetically induced ring (MIRC) currents, the angular dependence of the MIRC strength of the core, valence σ and π orbitals, the angular dependence of the sum of the MIRC strengths of all orbitals in each irreducible representation, pictures of the molecular structures of the studied molecules, and the Cartesian coordinates of the studied molecules (zipped). See DOI: https://doi.org/10.1039/d5sc00627a

This journal is © The Royal Society of Chemistry 2025
Click here to see how this site uses Cookies. View our privacy policy here.