ORIGINAL RESEARCH article
Sec. Theoretical and Computational Chemistry
Volume 9 - 2021 | https://doi.org/10.3389/fchem.2021.746162
The GW/BSE Method in Magnetic Fields
- 1Institute of Theoretical Solid State Physics, Karlsruhe Institute of Technology (KIT), Karlsruhe, Germany
- 2Institute of Physical Chemistry, Karlsruhe Institute of Technology (KIT), Karlsruhe, Germany
- 3Institute of Nanotechnology, Karlsruhe Institute of Technology (KIT), Eggenstein-Leopoldshafen, Germany
The GW approximation and the Bethe–Salpeter equation have been implemented into the Turbomole program package for computations of molecular systems in a strong, finite magnetic field. Complex-valued London orbitals are used as basis functions to ensure gauge-invariant computational results. The implementation has been benchmarked against triplet excitation energies of 36 small to medium-sized molecules against reference values obtained at the approximate coupled-cluster level (CC2 approximation). Finally, a spectacular change of colour from orange to green of the tetracene molecule is induced by applying magnetic fields between 0 and 9,000 T perpendicular to the molecular plane.
The description of excited states of molecules in strong magnetic fields poses a major challenge for quantum chemical methods. (Delos et al., 1983; Turbiner and López Vieyra, 2004; Hampe and Stopkowicz, 2017; Stopkowicz, 2018; Hampe and Stopkowicz, 2019; Wibowo et al., 2021). On the one hand, it is well known that introducing magnetic fields also introduces a gauge-dependence when standard, real-valued Gaussian-type basis functions are used. As a solution, as proposed by London, a complex phase factor countering the gauge-dependence of the magnetic field, can be used. (London, 1937; Helgaker and Jørgensen, 1991; Ruud et al., 1993; Tellgren et al., 2008). This in turn leads to complex-valued basis functions, which significantly increase the cost of subsequent calculations. On the other hand, many “work-horse” methods used to describe excited states as linear-response (LR) time-dependent density functional theory (TD-DFT) cannot be straightforwardly adapted to include arbitrary magnetic fields due to instabilities occurring in the respective non-collinear exchange-correlation (XC) kernel. The instabilities in the XC kernel are related to the same instabilities that also plague other non-collinear TD-DFT kernels in, for example, relativistic two-component TD-DFT. (Gao et al., 2005; Egidi et al., 2017; Komorovsky et al., 2019). While solutions to these problems have been proposed, they inevitably lead to XC kernels that do not exhibit full rotational invariance. (Egidi et al., 2017; Komorovsky et al., 2019). Contrary to TD-DFT, coupled-cluster methods are not plagued by any instabilities, but suffer from their steep cost, which increases exponentially with their accuracy. Furthermore, the complex gauge-independent London atomic orbitals lead to another steep increase in the computational complexity, effectively preventing calculations on systems with more than a few electrons. (Hampe and Stopkowicz, 2017; Hampe et al., 2020). Even though the computational limitations are severe, the investigation of molecular properties in strong external magnetic fields has become an increasingly popular topic within the field of quantum chemistry in recent years. Several field-dependent properties including non-linear effects on the electronic structure of small molecules, (Tellgren et al., 2008; Tellgren et al., 2009; Lange et al., 2012; Stopkowicz et al., 2015), molecular geometries, (Tellgren et al., 2012; Irons et al., 2021), spin-phase transitions (Sun et al., 2019a) and excited state properties (Sun et al., 2019b; Sen et al., 2019; Stetina et al., 2019; Wibowo et al., 2021) have been explored using quantum-chemical methods at different levels of theory.
Since the largest magnetic field currently created on Earth exhibits a field strength of about 100 T, (Sims et al., 2008), there is hardly any need to treat strong magnetic fields in more than a perturbative manner from an experimental point of view. Still, scientific curiosity has for a long time been a strong motor to investigate also situations which are (currently) not directly accessible. Given the lack of experimental data, highly accurate quantum-chemical methods are desirable in order to explore molecular properties in the field regime of
With the GW/Bethe–Salpeter equation (BSE) method, a suitable way of calculating properties from Kohn-Sham (KS) reference states has emerged within the last few years. (Bruneval et al., 2015; Jacquemin et al., 2015; Leng et al., 2016; Holzer and Klopper, 2017; Krause and Klopper, 2017; Gui et al., 2018; Blase et al., 2020; Kehry et al., 2020). It has seen great success, exhibiting a more favourable behaviour than TD-DFT on many occasions. While both TD-DFT and the GW/BSE method start from the same Kohn-Sham reference, GW/BSE fully accounts for charge-transfer and Rydberg excitations due to its correct asymptotic long-range behaviour. (Sagmeister and Ambrosch-Draxl, 2009; Blase and Attaccalite, 2011; Blase et al., 2011; Blase et al., 2018). Furthermore, the description of core excitations is significantly improved within the GW/BSE method. (Olovsson et al., 2009; Vinson et al., 2011; Kehry et al., 2020). The accuracy of the GW/BSE method is an improvement over TD-DFT. Therefore, adapting the GW/BSE method to be applicable to arbitrary molecules in arbitrary magnetic fields is worthwhile. It allows for an investigation of the effects of strong magnetic fields in sizable molecular systems while still retaining a certain robustness with respect to accuracy.
Within this paper we therefore aim at describing a fully consistent formulation and implementation of the GW/BSE method for the description of optical spectra of sizable molecules within strong magnetic fields. In the following chapters, the general formulas for the G0W0 and the eigenvalue self-consistent GW (evGW) methods as well as the BSE in strong magnetic fields are outlined. The resulting implementation is able to describe excited states of molecules of significant size. As such strong external magnetic fields are not accessible in experimental setups, a set of benchmark values obtained from truncated coupled cluster theory is provided for 36 small to medium-sized molecules. Finally, we demonstrate the capabilities of the GW/BSE equation in strong magnetic fields by predicting the colour change of tetracene in a strong uniform magnetic field.
2.1 GW Approach in Magnetic Fields Using London Atomic Orbitals
GW quasiparticle (QP) energies form the basis for calculating excitation energies from the Bethe–Salpeter equation. The principal theory to obtain GW QP energies in a magnetic field has been outlined in Ref. (Holzer et al., 2019). for atoms and complex-valued spinors. For molecules, to retain full gauge-invariance, instead of real Gaussian-type atomic orbitals, complex London-type atomic orbitals (LAOs) have to be used. These are obtained as a direct product of a Gaussian-type orbital ϕμ(r) and a complex phase factor:
The complex phase factor is used in order to cancel the dependency of all observable properties on the gauge origin O which naturally arises from the choice of a Coulomb gauge (∇ ⋅A = 0) for a magnetic vector potential (∇ ×A = B). In a two-component (2c) framework, complex spinors can be constructed as a linear combination of LAOs:
Non-collinear spin densities are well represented in this 2c spinor framework. Therefore, uniform and non-uniform magnetic fields can be included in this way. (Sen et al., 2019). More generally, within the notation used in this paper, any arbitrary non-collinear spin density can be employed. Furthermore, the complex phase-factor including LAOs are strictly needed to ensure gauge-independence for GW quasiparticle energy evaluations of multi-atomic systems, as well as for consecutive calculations of excitation energies using the Bethe–Salpeter equation. As the magnetic field is represented by a one-electron operator within 2c Kohn–Sham equations, the according information is fully absorbed into the complex spinors expanded in LAOs. Therefore, all quantities occurring in the BSE in a magnetic field must generally be assumed to be complex, unless further symmetries can be exploited.
where m denotes an excited state, and where the space-spin-coordinate x ≡ (r, σ) includes both space and spin coordinates. The charge fluctuation can be expressed using molecular spinors as
In the GW approximation, the correlation self-energy is obtained from the expression
where G is the one-electron Green’s function
As usual, to avoid instabilities and to make Eq. 6 integratable, a small positive number δ is added to the denominator. ɛp is the eigenvalue of the pth spinor that solves the Kohn-Sham equation for the underlying density functional approximation. The Fermi-level chemical potential μ is chosen to lie between the energy levels of the lowest unoccupied and highest occupied spinors, and Wc is the correlation contribution to the linearly screened potential,
Evaluating the integral on the right-hand side of Eq. 5 yields
where η = 2δ. We thus obtain the following working equation for the real-valued correlation contribution to the quasiparticle energy:
The two-electron integrals (pq|ρm) are computed as
The exchange self-energy is
where ΣX is the exchange self-energy and Vxc is the exchange-correlation potential of the underlying density functional theory. To obtain eigenvalue self-consistent quasiparticle energies (i.e., evGW quasiparticle energies), Eq. 13 is evaluated repeatedly with Z = 1 until the obtained eigenvalues are converged.
It was found that DIIS (direct inversion in the iterative subspace) (Pulay, 1980) procedures can speed up this process considerably. Usually less than ten consecutive evaluations of Eq. 15 are then needed to obtain converged evGW quasiparticle energies.
Finally, we note that also the analytic continuation (AC) and contour deformation (CD) GW variants described in Ref. (Holzer and Klopper, 2019) can be adapted to LAOs in the same manner. However, unlike the previous formulas derived for the analytic GW variant in an magnetic field, our current AC-GW and CD-GW variants are approximate in the sense that they ignore the lack of time-reversal (Kramers) symmetry. While we expect our AC-GW and CD-GW variants to be well behaved in a system with a vanishing spin expectation value (⟨S2⟩ ≈ 0), more research on these methods has to be performed in cases of non-vanishing ⟨S2⟩.
2.2 The Bethe–Salpeter Equation in a Magnetic Field
Starting from the gauge-invariant quasiparticle energies described in the previous section, the gauge-invariant excitation energies can be obtained from the Bethe–Salpeter equation also making use of LAOs. The BSE can be expressed in terms of complex spinors as
The orbital rotation matrices A and B are defined as
where ɛi is the quasiparticle energy of the ith Kohn–Sham eigenstate from a preceding GW computation, via,bj is a Coulomb integral over complex spinors,
and Wpq,rs is the static screened potential from the BSE. Properties such as for example oscillator strengths, excited state dipole moments, or nuclear forces can be obtained in a straightforward manner from the solutions of the eigenvalue problem of the BSE, again expressed in a basis of LAOs. Using complex LAO-based spinors, the static (i.e., ω = 0) screened potential W, which is given in its real-space expression in Eq. 7, takes the form
where ϵ is the dielectric function. The non-interacting response-function χ0 is diagonal and real if the quasiparticle energies are real, even if complex LAO-based spinors are used:
From the response function, and using the resolution-of-the-identity (RI) approximation
the screened potential Wpq,rs can be evaluated as (Krause and Klopper, 2017)
For the 3-index intermediate
the auxiliary functions ϕP are chosen to be real as ordinary Gaussian-type atomic orbitals without losing gauge invariance of the results obtained from computations in a magnetic field. This considerably simplifies the inner part of Eq. 22, representing the response function in the auxiliary subspace,
which is symmetric and real for the special case of the static BSE even in a (uniform or nonuniform) magnetic field. Finally, the efficient evaluation of the 3-index integrals (ϕQ|φpφq) has been described in Ref. (Pausch and Klopper, 2020). Therefore, the evaluation of the BSE in magnetic fields can proceed in a straightforward manner, making it an invaluable tool to assess excited stats of molecules in magnetic fields at roughly the same cost as required for linear-response Hartree–Fock computations, with the advantage of being significantly more accurate.
3 Computational Details
All implementation work in this work have been carried out in the framework of the Turbomole (Ahlrichs et al., 1989; Furche et al., 2014; Balasubramani et al., 2020) program package. Consequently, all calculations have also been done using Turbomole.
To assess the correctness of the implementation and theory of the GW/BSE method in magnetic fields, the first excited triplet states of 36 small molecules have been evaluated using evGW/BSE in a magnetic field of 1,000 T in z-direction with respect to the coordinates supplied in the Supporting Information (SI). The molecules used in this evaluation are acetaldehyde, acetylene, CCl2, CClF, CF2, cyanoacetylene, cyanoformaldehyde, cyanogen, diacetylene, difluorodiazirine, formaldehyde, formic acid, formyl chloride, formyl fluoride, glyoxal, H2C3, HCN, HCP, HNO, HPO, HPO, HPS, HSiF, isocyanogen, nitrosamine, nitrosylcyanide, phosgene, propynal, pyrazine, selenoformaldehyde, SiCl2, silylidene, tetrazine, thioformaldehye, thioformylchloride, thionylcarbonylfluoride, and thiophosgene. The set of molecules is taken from Ref. (Suellen et al., 2019). Starting from the geometries provided in Ref. (Suellen et al., 2019), the geometries have been re-optimized using RI-MP2 (resolution-of-identity Møller–Plesset perturbation theory to second order) in the corresponding magnetic field of 1,000 T in z-direction using numerical gradients. At the RI-MP2-optimized geometries, the three first excited states have been evaluated using evGW/BSE and PBE0, (Perdew et al., 1997; Adamo and Barone, 1999), LC-ωPBE, (Vydrov and Scuseria, 2006), BHLYP (also known as BH&HLYP), (Becke, 1993), and CAM-B3LYP (Yanai et al., 2004) as underlying functionals. In all calculations the energy and norm of the difference density matrix were converged to 10−8 hartree and 10−7, respectively. The def2-TZVP basis set (Weigend and Ahlrichs, 2005) was used throughout, in conjunction with the resolution-of-identity (RI) approximation for the Hartree and exchange terms with the corresponding auxiliary fitting basis sets for the Kohn–Sham ground-state (Weigend, 2006; Weigend, 2008) as well as the appropriate auxiliary fitting basis sets for the RI-MP2 and GW/BSE calculations. (Hättig, 2005).
For further comparison, the corresponding excited state energies have also been determined using the approximate coupled-cluster RI-CC2 method, (Hättig and Weigend, 2000), which has been adapted to calculations in finite magnetic fields in the course of the present work. It is closely related to the equation-of-motion coupled-cluster singles-and-doubles (EOM-CCSD) method in magnetic fields that has been described by Hampe and Stopkowicz, (Hampe and Stopkowicz, 2017), and to the two-component RI-CC2 implementation of Krause and Klopper. (Krause and Klopper, 2015). Compared to EOM-CCSD, RI-CC2 is computationally significantly less involved. This allows for the assessment of the larger molecules in the test set in the applied magnetic field.
For tetracene, we performed calculations on the evGW (10)/BSE@DFT level using the contour deformation (CD) technique. (Holzer and Klopper, 2019). CD-evGW (10) denotes that only the highest 10 occupied and lowest 10 unoccupied spinor energies have been corrected using CD-evGW, while the remaining spinor energies are shifted (“scissoring”) accordingly. Testing the self-energy obtained from the CD-GW variant reveals that indeed for the systems and magnetic field strengths investigated in this paper, CD-GW exhibits errors of the order of 1 meV or less, making it perfectly feasible for (not too strong) magnetic fields. As reference density functionals for the evGW(10)/BSE@DFT calculations, we used PBE0, B3LYP, BHLYP and CAM-B3LYP. Also the tetracene calculations were carried out in the def2-TZVP basis.
4 Results and Discussion
4.1 Test Set of Small Molecules
For the 36 molecules tested, in a field of 1,000 T, all ground states retain their closed-shell character, yielding no spin polarization. Therefore, the spacing in-between the three triplet states (T−1,0,1), which are non-degenerate in the magnetic field, are solely determined by the Zeeman effect. The T−1,1 components of the triplet are found exactly at
Figure 1 compares the RI-CC2 and evGW/BSE excitation energies of the full set. It exhibits a near-linear shift between the two methods, with the difference getting more pronounced for excited states with higher energy. Furthermore, RI-CC2 consistently yields blue-shifted excitation energies when compared to evGW/BSE. This is in accordance with the finding of Suellen et al., who also found CC2 to yield too high excitation energies on average for this test set, though in the field-free case. In contrast, evGW/BSE was found to yield too low excitation energies on average, especially for triplet excited states. (Gui et al., 2018). Too low triplet excited states are well known phenomenon for the GW/BSE method. (Rangel et al., 2017; Jacquemin et al., 2017; Holzer and Klopper, 2018). While Jacquemin et al. proposed to use the Tamm-Dancoff approximation (TDA) to improve this, (Jacquemin et al., 2017), two of us proposed adding the correlation kernel of the underlying density functional approximation to improve triplet excitations, yielding the correlation-kernel augmented BSE (cBSE) approach. (Holzer and Klopper, 2018). While a linear-response time-dependent DFT implementation in magnetic fields would be needed to apply the cBSE method, using the TDA is straightforwardly obtained by setting B = 0 in Eq. 16. For the tested molecules, these findings can be partly confirmed, with the TDA leading to a significantly blue shift of especially the lower lying excited states, improving the agreement between evGW/TDA-BSE and RI-CC2 (Figure 2). The improvement of the GW/BSE method when using the TDA is generally also observed in field-free cases. (Rangel et al., 2017).
FIGURE 1. Comparison of excitation energies from the RI-CC2 and evGW/BSE@DFT methods. evGW/BSE calculations have been performed using either PBE0, LC-ωPBE, BHLYP or CAM-B3LYP spinors. All calculations used the def2-TZVP basis set. All values in eV.
FIGURE 2. Comparison of excitation energies from the RI-CC2 and evGW/TDA-BSE@DFT methods. evGW/BSE calculations have been performed using either PBE0 or LC-ωPBE spinors. All calculations used the def2-TZVP basis set. All values in eV.
A closer inspection of Figure 1 reveals that for the evGW/BSE method, two distinct groups, one with smaller deviations from RI-CC2 and one with larger deviations, are found. The “high-error” group is composed of the molecules nitrosamine (1.02 eV), HCP (3.61 eV), diacetylene (4.20 eV), cyanoacetylene (4.55 eV), cyanogen (5.03 eV), isocyanogen (5.42 eV), acetylene (5.65 eV), and HCN (6.58 eV). The values in parenthesis are the T0 excitation energies of the corresponding RI-CC2 references. Except for the low-energy excited state of nitrosamine, all these molecules feature triple bonds in their respective ground states. Furthermore we find instabilities for the molecules HNO (0.74 eV) and nitrosylcyanide (0.72 eV) with rather low lying triplet excited states. This is further hinting at triple bonds and nitrosyl groups being described with sub-par quality within the evGW/BSE methods. For the remaining molecules significantly smaller errors are found.
Employing the TDA removes the instabilities encountered in the evGW/BSE calculations for the molecules HNO and nitrosylcyanide, in both cases yielding excitation energies that are lower by
TABLE 1. Mean average error (MAE), mean signed error (MSE), standard deviation (SD), and maximum error (MAX) of evGW/BSE@DFT (“BSE”) and evGW/TDA-BSE@DFT (“TDA”) excitation energies with respect to CC2 excitation energies. All values in eV.
As displayed also in Figure 1 and Figure 2, Table 1 reveals that BHLYP, which incorporates a relatively large amount of Hartree–Fock exchange, performs best for the investigated molecules. The range-separated hybrids LC-ωPBE and CAM-B3LYP yield comparable results, and generally perform better than PBE0 but worse than BHLYP. This is in line with observations for field-free cases, indicating that (at least for moderate field strength) conclusions drawn from field-free benchmarks are still applicable. (Holzer et al., 2021).
As shown in Figure 2, the class of molecules with triple bonds or nitrosyl groups exhibits a significantly reduced error within the TDA for all investigated functionals. Triplet excitation energies from the latter class of molecules are now in line with all other molecules. We therefore expect the TDA to be especially valuable for molecules with triple bonds or nitrosyl groups. Still, regarding the TDA, there are some caveats left. While some of the improvements can indeed be related to error compensation, where the blue-shift of the TDA counteracts the general red-shift of the evGW/BSE method with respect to CC2 excitation energies, this can not fully explain the strong reduction of the error regarding the class of molecules with triple bonds or nitrosyl groups, which indicates that also the correlation from the BSE is sometimes insufficient to describe triplet excitations sufficiently well. Given the overall increase in accuracy from the TDA, it may be advisable to even use it by default in magnetic fields until the cBSE method becomes available. (Holzer and Klopper, 2018). However, it shall be noted that the usage of CC2 as reference method is not the best possible but a pragmatic choice for this test set. While its accuracy is comparable or even slightly better than that of EOM-CCSD, (Suellen et al., 2019), more refined methods as CC3 or EOM-CCSDT would be needed to obtain true reference values with errors significantly below 0.1 eV. Given the immense computational cost of the latter two methods, only results for a single diatomic molecule, namely CH+, have been reported for EOM-CCSDT so far in a finite magnetic field. (Hampe et al., 2020). RI-CC2 as computational efficient method is therefore a suitable compromise, providing robust values. However, as shown in Ref. (Suellen et al., 2019), CC2 has a tendency to deliver too high excitation energies when compared to experimental and CC3 excitation energies. In contrast, evGW/BSE tends to underestimate excitation energies as shown in Ref. (Gui et al., 2018), especially for triplet excited states. This has to be taken into account when comparing the CC2 and evGW/BSE methods. Concerning the reference state, evGW is able to even out the differences between the underlying functionals completely. The difference between excitation energies obtained from either evGW/BSE@PBE0 or evGW/BSE@LC-ωPBE is statistically insignificant. The presented results suggest that the performance of the evGW/BSE method in magnetic fields is similar to its performance in field-free situations, yielding good to excellent excitation energies, at a considerably reduced effort when compared to coupled-cluster methods.
4.2 Optical Properties of Tetracene in a Magnetic Field
As established in the last section, the evGW/BSE method quite accurately predicts molecular excitation energies in the presence of an external magnetic field. Our implementation into the Turbomole package thus appears to be an efficient yet reliable method of predicting excited state properties of sizable molecules in strong magnetic fields. As a consequence, real-world properties such as the absorption and emission spectra, and therefore also the colour, of a substance can now be obtained by simulating the vertical excitations and related oscillator strengths of a molecule under such extreme conditions. In this section, the effects of a strong external magnetic field on the excited states of tetracene are studied in detail, exemplifying the effects such extreme environments can have on chemical substances.
The optical spectrum of tetracene in the absence of an external magnetic field is mainly composed of three bands. The p-band (peak at 455–477 nm) corresponds to the HOMO → LUMO transition (B2u symmetry). The α- and β-bands correspond to the HOMO → LUMO+1 and HOMO −1 → LUMO transitions (both B3u symmetry), respectively, and show peaks at 373–393 nm (α-band) as well as 272–275 nm (β-band). While the two bands in the visible (p) and near-UV (α) region of the spectrum exhibit relatively small oscillator strengths, the β excitation is associated with an oscillator strength several orders of magnitude larger than that of the two other transitions. (Guidez and Aikens, 2013; Sony and Shukla, 2007; Klevens and Platt, 1949; Biermann and Schmidt, 1980; Bree and Lyons, 1960).
In order to investigate the optical properties of tetracene for the field-free case, we first optimized the geometry at the PBE0, B3LYP, BHLYP and CAM-B3LYP levels, respectively. Using these structures, subsequent CD-evGW/BSE@DFT calculations were carried out. The resulting wavelengths for the excitations are presented in Table 2. The different reference functionals provide similar values for the α- and β-excitations, slightly overestimating the energies of both excitations. The energies of the p-excitation, however, vastly differ with respect to the reference functional, ranging from 465 nm at the CD-evGW/BSE@BHLYP level to 614 nm at the CD-evGW/BSE@PBE0 level. Further investigations reveal that this is almost exclusively an effect of the geometry and not the method itself as using the CAM-B3LYP reference geometry yields wavelengths between 468 and 487 nm for the p-excitation in all cases. This is to be expected as the impact of the reference functional is usually not very large for evGW/BSE calculations. (Holzer et al., 2021).
TABLE 2. Wavelengths of p-, α- and β-excitations in nm as calculated at the CD-evGW/BSE@DFT level employing different reference functionals. In order to highlight the dependence on the geometry, all calculations were performed both on the geometry as optimized using the reference functional as well as on the geometry as optimized using the range-separated hybrid (RSH) functional CAM-B3LYP. Experimental values are also given.
In order to gain a better insight into the electronic processes behind the optical spectrum, the transition densities are examined for the p-, α- and β-excitations. In accordance with the literature, (Lim et al., 2004; Guidez and Aikens, 2013), it is found that while the p-band corresponds to a transition polarized along the short axis of the molecule, the α- and β-excitations can be described by transitions polarized along the long axis of the molecule. The transition densities of the three excitations are depicted in Figure 3.
FIGURE 3. Transition densities of the p-, α- and β-excitations as generated by a two-component CD-evGW/BSE@CAM-B3LYP calculation. The plots were generated with VMD using an isovalue of 0.0001
Having gained a general overview of the electronic excitations that mainly constitute the optical spectrum of tetracene, it is now possible to extend these findings in order to understand the effect a strong external magnetic field may have on such a system. Furthermore, it is now possible to use these findings in order to predict the colour shift of the tetracene molecule under the influence of such a strong magnetic field.
This investigation contains three steps: Firstly, it is necessary to generate the optical spectrum of tetracene in the absence of a magnetic field. This was done using the CD-evGW (10)/BSE@DFT methods employing BHLYP and CAM-B3LYP as reference functionals as they most accurately describe the electronic excitations in the zero-field compared to the experimental values. The peaks are broadened using a damped sum-over-states formalism which translates to a Lorentzian line shape with full width at half maximum of 0.15 eV for all excitations. (Norman et al., 2004; Barron, 2004; Fernandez-Corbaton et al., 2020). From this information, the RGB colour code of the substance can be computed. By calculating the integrals over the entire visible part of the spectrum and arbitrarily setting it to 1 for the zero-field case, the relative intensities of the colour may also be calculated.
Secondly, the immediate influence of an external magnetic field has to be assessed. Applying the external field perpendicular to the molecular plane lowers the symmetry of the system. The point group of tetracene in such an external field is C2h instead of D2h. The excitations of the p-, α- and β-bands are all of Bu symmetry. A further investigation reveals that the subsequent two excitations (here denoted as γ and δ) are also of Bu symmetry.
Thirdly, by slowly raising the magnetic field strength in steps of 1,000 T and generating the optical spectrum at each field strength as previously described, it is possible to track how the excitations are influenced by the external magnetic field. The resulting UV/Vis spectra and the energies of the five lowest excitations are plotted in Figures 4A,B. While the excitation energies of the α–δ-excitations are only slightly shifted between 0 T and 9,000 T, the p-band is strongly red-shifted. At the same time, the oscillator strength of the p-excitation decreases with an increasing field strength. The resulting UV/Vis spectrum is mostly dominated by both the location and intensity of the p-band while only a small section of the violet and blue part of the spectrum between 400 and 500 nm is caused by the other excitations. Only the five lowest excitations (p, α–δ) are depicted in Figure 4B, but the UV/Vis spectrum was generated by the 75 lowest excitations. However, none of the other excitations contribute significantly to the visible part of the spectrum.
FIGURE 4. (A) UV/Vis spectra of tetracene as predicted at different magnetic field strengths between 0 T and 9,000 T. Solid lines denote calculations using the CD-evGW (10)/BSE@BHLYP method while dashed lines denote CD-evGW (10)/BSE@CAM-B3LYP calculations. (B) Wavelengths of the relevant lowest vertical excitations of tetracene at different magnetic field strengths between 0 T and 9,000 T. The p-excitation, which is predominantly responsible for the colour of tetracene, is most affected by the external field. Solid lines denote calculations using the CD-evGW (10)/BSE@BHLYP method while dashed lines denote CD-evGW (10)/BSE@CAM-B3LYP calculations. (C) Colour of tetracene as predicted at the CD-evGW (10)/BSE@CAM-B3LYP level of theory at different magnetic field strengths between 0 T and 9,000 T. To obtain the depicted colours, the vertical excitations of the optical spectrum were broadened by 0.15 eV and converted into a RGB colour code while the intensity was scaled relative to the zero-field by integrating over the visible region of the spectrum.
The resulting predicted colours of tetracene at various magnetic field strengths are shown in Figure 4C. As the p-band is red-shifted, the main absorption band moves from the blue part of the spectrum towards the red. The orange-red colour of tetracene therefore shifts towards the colour blue. Finally, between 5,000 T and 9,000 T, the contribution of the p-band does not contribute significantly to the spectrum anymore. Since the α–δ-bands are still active and relatively unchanged in their location, part of the blue light is still absorbed, resulting in a turquoise to green colour which becomes less and less intense as the magnetic field strength is increased.
It is worth noting that while the peak positions of the α–δ-excitations exhibit only a minor dependence on the applied field, their respective oscillator strengths change significantly as the β-excitation becomes less important. Subsequently, the γ- and δ-excitations become dominant at different magnetic field strengths. Furthermore, certain additional transitions slowly start to arise as they are no longer symmetry forbidden due to a lowering of the point group symmetry in the magnetic field. At a magnetic field strength of approximately 6,000 T, specifically, the γ- and δ-excitations are very close energetically to the usually predominant β-excitation, leading to resonance phenomena such as a splitting into multiplets.
Finally, the strong influence of the external field on the excitation energy of the p-band can best be understood by examining the changes to the electronic structure of tetracene. In the magnetic field, the energy of the HOMO increases and the energy of the LUMO decreases. Thus, the HOMO-LUMO gap decreases significantly. As the p-excitation corresponds to the HOMO → LUMO transition, the resulting excitation energy is subsequently lowered.
The transition density of the p-excitation at a magnetic field strength of 8,000 T is depicted in Figure 5. It exhibits the effects the magnetic field has on this most important transition, showing a slightly more delocalized nature of this transition in the magnetic field compared to the zero-field case. In order to ensure gauge-origin invariance, the transition density plot was generated employing London atomic orbitals.
FIGURE 5. Transition density of the p-excitation at a magnetic field strength of 8,000 T as generated by a two-component CD-evGW/BSE@CAM-B3LYP calculation. The plot was generated with VMD using an isovalue of 0.0001
In this paper, we have presented a gauge-invariant formulation of the GW/BSE method for excited states in strong magnetic fields. The resulting implementation was benchmarked against reference values obtained from approximate coupled-cluster (CC2) theory. The obtained results indicate that the GW/BSE method provides a similar accuracy in strong magnetic fields as in the field-free case. The known issue of an underestimation of the excitation energy of triplet excited states is also present in magnetic fields. Like in the field-free case, it is shown that the Tamm-Dancoff approximation is able to remove a significant amount of this underestimation, improving the overall accuracy when compared to coupled-cluster values. The remaining error is nearly linear, making it easy to be accounted for.
Furthermore, using the tetracene molecule as showcase example, it was demonstrated that the GW/BSE method is able to tackle systems far beyond the possibilities of any prior ansatz that has been used to describe excited states in strong magnetic fields. For the tetracene molecule, we analyzed the shift of the main absorption peaks in magnetic fields ranging from 0 to 9,000 T. It was found that some excited state energies are more affected than others, leading to prominent changes in the spectrum. Ultimately, the colour of tetracene was estimated from the calculated spectra in the assessed magnetic fields. Starting from the bright orange colour of tetracene, we predict the compound to exhibit a blue colour at 5,000, which is converted towards a green colour at 9,000 T. While the dependence of the excited states on the external magnetic field are interesting on their own, the example of tetracene also outlines the fascinating world that even moderately strong magnetic fields could open for the broad field of photochemistry.
To summarize, the GW/BSE method has proven once more that it has become a formidable member of the toolbox of quantum chemistry.
Data Availability Statement
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.
CH, AP, and WK have formulated the theory and derived the equations. CH and AP have implemented the formalism in a basis set of London atomic orbitals. CH has performed the calculations on the test set of 36 molecules and AP has performed the calculations on tetracene. CH, AP and WK have written parts of the manuscript.
Conflict of Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
CH gratefully acknowledges the Volkswagen Stiftung for financial support. AP gratefully acknowledges financial support by Fonds der chemischen Industrie and Studienstiftung des deutschen Volkes.
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fchem.2021.746162/full#supplementary-material
Ahlrichs, R., Bär, M., Häser, M., Horn, H., and Kölmel, C. (1989). Electronic Structure Calculations on Workstation Computers: The Program System Turbomole. Chem. Phys. Lett. 162, 165–169. doi:10.1016/0009-2614(89)85118-8
Balasubramani, S. G., Chen, G. P., Coriani, S., Diedenhofen, M., Frank, M. S., Franzke, Y. J., et al. (2020). TURBOMOLE: Modular Program Suite for Ab Initio Quantum-Chemical and Condensed-Matter Simulations. J. Chem. Phys. 152, 184107. doi:10.1063/5.0004635
Blase, X., and Attaccalite, C. (2011). Charge-transfer Excitations in Molecular Donor-Acceptor Complexes Within the Many-Body Bethe-Salpeter Approach. Appl. Phys. Lett. 99, 171909. doi:10.1063/1.3655352
Blase, X., Attaccalite, C., and Olevano, V. (2011). First-PrinciplesGWcalculations for Fullerenes, Porphyrins, Phtalocyanine, and Other Molecules of Interest for Organic Photovoltaic Applications. Phys. Rev. B. 83, 115103. doi:10.1103/physrevb.83.115103
Bree, A., and Lyons, L. E. (1960). 1002 The Intensity of Ultraviolet-Light Absorption by Monocrystals. Part IV. Absorption by Naphthacene of Plane-Polarized Light. J. Chem. Soc. 1960, 5206–5212. doi:10.1039/JR9600005206
Bruneval, F., Hamed, S. M., and Neaton, J. B. (2015). A Systematic Benchmark of the Ab Initio Bethe-Salpeter Equation Approach for Low-Lying Optical Excitations of Small Organic Molecules. J. Chem. Phys. 142, 244101. doi:10.1063/1.4922489
Egidi, F., Sun, S., Goings, J. J., Scalmani, G., Frisch, M. J., and Li, X. (2017). Two-Component Noncollinear Time-Dependent Spin Density Functional Theory for Excited State Calculations. J. Chem. Theor. Comput. 13, 2591–2603. doi:10.1021/acs.jctc.7b00104
Fernandez-Corbaton, I., Beutel, D., Rockstuhl, C., Pausch, A., and Klopper, W. (2020). Computation of Electromagnetic Properties of Molecular Ensembles. ChemPhysChem. 21, 878–887. doi:10.1002/cphc.202000072
Gao, J., Zou, W., Liu, W., Xiao, Y., Peng, D., Song, B., et al. (2005). Time-Dependent Four-Component Relativistic Density-Functional Theory for Excitation Energies. II. The Exchange-Correlation Kernel. J. Chem. Phys. 123, 054102. doi:10.1063/1.1940609
Gui, X., Holzer, C., and Klopper, W. (2018). Accuracy Assessment of GW Starting Points for Calculating Molecular Excitation Energies Using the Bethe-Salpeter Formalism. J. Chem. Theor. Comput. 14, 2127–2136. doi:10.1021/acs.jctc.8b00014
Hampe, F., Gross, N., and Stopkowicz, S. (2020). Full Triples Contribution in Coupled-Cluster and Equation-Of-Motion Coupled-Cluster Methods for Atoms and Molecules in strong Magnetic fields. Phys. Chem. Chem. Phys. 22, 23522–23529. doi:10.1039/d0cp04169f
Hampe, F., and Stopkowicz, S. (2019). Transition-Dipole Moments for Electronic Excitations in Strong Magnetic Fields Using Equation-Of-Motion and Linear Response Coupled-Cluster Theory. J. Chem. Theor. Comput. 15, 4036–4043. doi:10.1021/acs.jctc.9b00242
Hättig, C. (2005). Optimization of Auxiliary Basis Sets for RI-MP2 and RI-CC2 Calculations: Core-Valence and Quintuple-ζ Basis Sets for H to Ar and QZVPP Basis Sets for Li to Kr. Phys. Chem. Chem. Phys. 7, 59–66. doi:10.1039/b415208e
Hedin, L. (1991). Properties of Electron Self-Energies and Their Role in Electron Spectroscopies. Nucl. Instr. Methods Phys. Res. Section A: Acc. Spectrometers, Detectors Associated Equipment. 308, 169–177. doi:10.1016/0168-9002(91)90619-2
Holzer, C., Franzke, Y. J., and Kehry, M. (2021). Assessing the Accuracy of Local Hybrid Density Functional Approximations for Molecular Response Properties. J. Chem. Theor. Comput. 17, 2928–2947. doi:10.1021/acs.jctc.1c00203
Holzer, C., and Klopper, W. (2017). Communication: Symmetry-Adapted Perturbation Theory With Intermolecular Induction and Dispersion Energies From the Bethe-Salpeter Equation. J. Chem. Phys. 147, 181101. doi:10.1063/1.5007929
Holzer, C., and Klopper, W. (2019). Ionized, Electron-Attached, and Excited States of Molecular Systems With Spin-Orbit Coupling: Two-Component GW and Bethe-Salpeter Implementations. J. Chem. Phys. 150, 204116. doi:10.1063/1.5094244
Holzer, C., Teale, A. M., Hampe, F., Stopkowicz, S., Helgaker, T., and Klopper, W. (2019). GW Quasiparticle Energies of Atoms in Strong Magnetic Fields. J. Chem. Phys. 150, 214112. doi:10.1063/1.5093396
Kehry, M., Franzke, Y. J., Holzer, C., and Klopper, W. (2020). Quasirelativistic Two-Component Core Excitations and Polarisabilities From a Damped-Response Formulation of the Bethe-Salpeter Equation. Mol. Phys. 118, e1755064. doi:10.1080/00268976.2020.1755064
Komorovsky, S., Cherry, P. J., and Repisky, M. (2019). Four-Component Relativistic Time-Dependent Density-Functional Theory Using a Stable Noncollinear DFT Ansatz Applicable to Both Closed- and Open-Shell Systems. J. Chem. Phys. 151, 184111. doi:10.1063/1.5121713
Krause, K., Harding, M. E., and Klopper, W. (2015). Coupled-Cluster Reference Values for the GW27 and GW100 Test Sets for the Assessment of GW Methods. Mol. Phys. 113, 1952–1960. doi:10.1080/00268976.2015.1025113
Krause, K., and Klopper, W. (2015). Description of Spin-Orbit Coupling in Excited States With Two-Component Methods Based on Approximate Coupled-Cluster Theory. J. Chem. Phys. 142, 104109. doi:10.1063/1.4908536
Lim, S.-H., Bjorklund, T. G., Spano, F. C., and Bardeen, C. J. (2004). Exciton Delocalization and Superradiance in Tetracene Thin Films and Nanoaggregates. Phys. Rev. Lett. 92, 107402. doi:10.1103/physrevlett.92.107402
Norman, P., Ruud, K., and Helgaker, T. (2004). Density-Functional Theory Calculations of Optical Rotatory Dispersion in the Nonresonant and Resonant Frequency Regions. J. Chem. Phys. 120, 5027–5035. doi:10.1063/1.1647515
Olovsson, W., Tanaka, I., Puschnig, P., and Ambrosch-Draxl, C. (2009). Near-Edge Structures From First Principles All-Electron Bethe-Salpeter Equation Calculations. J. Phys. Condens. Matter. 21, 104205. doi:10.1088/0953-8984/21/10/104205
Rangel, T., Hamed, S. M., Bruneval, F., and Neaton, J. B. (2017). An Assessment of Low-Lying Excitation Energies and Triplet Instabilities of Organic Molecules With an Ab Initio Bethe-Salpeter Equation Approach and the Tamm-Dancoff Approximation. J. Chem. Phys. 146, 194108. doi:10.1063/1.4983126
Sagmeister, S., and Ambrosch-Draxl, C. (2009). Time-Dependent Density Functional Theory Versus Bethe-Salpeter Equation: an All-Electron Study. Phys. Chem. Chem. Phys. 11, 4451–4457. doi:10.1039/b903676h
Sims, J. R., Rickel, D. G., Swenson, C. A., Schillig, J. B., Ellis, G. W., and Ammerman, C. N. (2008). Assembly, Commissioning and Operation of the NHMFL 100 Tesla Multi-Pulse Magnet System. IEEE Trans. Appl. Supercond. 18, 587–591. doi:10.1109/tasc.2008.922541
Sony, P., and Shukla, A. (2007). Large-scale Correlated Calculations of Linear Optical Absorption and Low-Lying Excited Statesof Polyacenes: Pariser-Parr-Pople Hamiltonian. Phys. Rev. B. 75, 155208. doi:10.1103/physrevb.75.155208
Stetina, T. F., Sun, S., Williams‐Young, D. B., and Li, X. (2019). Modeling Magneto‐Photoabsorption Using Time‐Dependent Complex Generalized Hartree‐Fock. ChemPhotoChem. 3, 739–746. doi:10.1002/cptc.201900161
Stopkowicz, S., Gauss, J., Lange, K. K., Tellgren, E. I., and Helgaker, T. (2015). Coupled-Cluster Theory for Atoms and Molecules in strong Magnetic Fields. J. Chem. Phys. 143, 074110. doi:10.1063/1.4928056
Suellen, C., Freitas, R. G., Loos, P.-F., and Jacquemin, D. (2019). Cross-Comparisons Between Experiment, TD-DFT, CC, and ADC for Transition Energies. J. Chem. Theor. Comput. 15, 4581–4590. doi:10.1021/acs.jctc.9b00446
Sun, S., Williams-Young, D. B., Stetina, T. F., and Li, X. (2019a). Generalized Hartree-Fock With Nonperturbative Treatment of Strong Magnetic Fields: Application to Molecular Spin Phase Transitions. J. Chem. Theor. Comput. 15, 348–356. doi:10.1021/acs.jctc.8b01140
Sun, S., Williams-Young, D., and Li, X. (2019b). An Ab Initio Linear Response Method for Computing Magnetic Circular Dichroism Spectra With Nonperturbative Treatment of Magnetic Field. J. Chem. Theor. Comput. 15, 3162–3169. doi:10.1021/acs.jctc.9b00095
Tellgren, E. I., Reine, S. S., and Helgaker, T. (2012). Analytical GIAO and Hybrid-Basis Integral Derivatives: Application to Geometry Optimization of Molecules in Strong Magnetic fields. Phys. Chem. Chem. Phys. 14, 9492–9499. doi:10.1039/c2cp40965h
Weigend, F., and Ahlrichs, R. (2005). Balanced Basis Sets of Split Valence, Triple Zeta Valence and Quadruple Zeta Valence Quality for H to Rn: Design and Assessment of Accuracy. Phys. Chem. Chem. Phys. 7, 3297–3305. doi:10.1039/b508541a
Wibowo, M., Irons, T. J. P., and Teale, A. M. (2021). Modeling Ultrafast Electron Dynamics in Strong Magnetic Fields Using Real-Time Time-Dependent Electronic Structure Methods. J. Chem. Theor. Comput. 17, 2137–2165. doi:10.1021/acs.jctc.0c01269
Keywords: GW, Bethe-Salpeter, excitation energy, magnetic field, density functional theory
Citation: Holzer C, Pausch A and Klopper W (2021) The GW/BSE Method in Magnetic Fields. Front. Chem. 9:746162. doi: 10.3389/fchem.2021.746162
Received: 23 July 2021; Accepted: 19 October 2021;
Published: 25 November 2021.
Edited by:Linn Leppert, University of Twente, Netherlands
Reviewed by:Xinguo Ren, Institute of Physics (CAS), China
Antonios Alvertis, Lawrence Berkeley National Laboratory, United States
Copyright © 2021 Holzer, Pausch and Klopper. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Wim Klopper, email@example.com