We perform structural relaxation and compute initial electronic wavefunctions and energies, using density functional theory (DFT) (Kohn and Sham, 1965) and the Perdew-Burke-Ernzerhof (PBE) (Perdew et al., 1996) exchange-correlation functional, as implemented in the Quantum Espresso package (Giannozzi et al., 2017) (see full computational details in the SI). The DFT Kohn-Sham eigenvalues and eigenfunctions are then taken as the first guess for the MBPT calculations. We use the GW approximation to compute the quasiparticle energies and bandstructure (Hybertsen and Louie, 1986; Deslippe et al., 2012; Hedin, 1965), where we calculate self-energy corrections via Σ = iGW, for G the single-particle Green’s function and W the screened Coulomb interaction, W G G ′ ( q ; 0 ) = ε G G ′ − 1 ( q ; 0 ) v ( q + G ′ ) , where G, G′ are reciprocal lattice vectors, v is the bare Coulomb interaction, and ɛ _GG′(q; 0) is the dielectric function of the system for interaction wavevector q and zero frequency, evaluated via (Hybertsen and Louie, 1986; Deslippe et al., 2012): ε G G ′ ( q ; 0 ) = δ G G ′ − v ( q + G ) χ G G ′ 0 ( q ; 0 ) , (1) with χ G G ′ 0 ( q ; 0 ) the non-interacting electronic polarizability calculated using the random-phase approximation.

Optical excitations and excitonic properties are computed through the BSE formalism (Rohlfing and Louie, 1998; Rohlfing and Louie, 2000) E c k + Q − E v k A v c k Q S + ∑ v ′ c ′ k ′ ⟨ v k ; c k + Q | K e h | v ′ k ′ ; c ′ k ′ + Q ⟩ A v ′ c ′ k ′ Q S = Ω Q S A v c k Q S (2) for a hole state v k ⟩ and an electron state c k + Q ⟩ , where k is the crystal momentum and Q is the exciton center-of-mass momentum. S indexes the exciton state at momentum Q. A v c k Q S is the amplitude of the free electron-hole pair. E _{c k+Q} and E _{v k} are the quasiparticle energies calculated within the GW approximation; Ω Q S is the excitation energy; and K ^eh is the electron-hole interaction kernel. Exciton dispersion is obtained by solving the BSE at different exciton momenta Q following the methodology developed in Refs. (Gatti and Sottile, 2013; Qiu et al., 2015). The GW-BSE exciton wavefunction can be represented as Ψ S , Q ⟩ = ∑ v c k A v c k Q S c k + Q ⟩ v k ⟩ .

We compute the dielectric screening at the C8-PDI@Au interface using two approaches. The first is a direct calculation using G₀W₀ for the interface, which is highly computationally challenging and is made possible due to recent advances both in large-scale computing capabilities and associated code development (Del Ben et al., 2019). As a comparison, we also compute the dielectric function of the interface using a recently developed substrate screening GW approach (Liu et al., 2019). In this approach, the non-interacting polarizability of the interface is approximated by the summation of the separately calculated polarizabilies of the substrate (Au) and adsorbate (C8-PDI), i.e., χ t o t , G G ′ 0 ( q ; 0 ) ≈ χ m o l , G G ′ 0 ( q ; 0 ) + χ s u b , G G ′ 0 ( q ; 0 ) , (3) where χ mol 0 and χ sub 0 are the polarizabilities associated with the standalone molecular layer and metal substrates, respectively. This approximation holds for systems with weak hybridization between the adsorbate and the substrate (Liu et al., 2019; Xuan et al., 2019; Adeniran et al., 2020). The goal is to assess the validity of Eq. 3 against a direct GW calculation of the interface, to determine the nature of the surface effect on quasiparticle and exciton properties of the C8-PDI.

As shown in Figure 1B, we mainly study four different structures: (i) “bulk”, a bulk crystal in its experimentally resolved triclinic unit cell and relaxed atomic coordinates; (ii) “com. bulk”, a modified bulk structure in an orthorhombic lattice, so that it is commensurate with the lattice of the Au substrate; (iii) “com. layer”, a layered orthorhombic structure with large vacuum, taken as one molecular layer along the c axis from (ii); and (iv) “interface”, the C8-PDI monolayer from (iii) adsorbed on four layers of Au (111) surface, with each layer consisting of 2 × 3 Au atoms and an adsorption height of 3.18 Å. Such thickness of the simulated Au surface was previously found to sufficiently capture wavefunction hybridization and dielectric screening effects (Liu et al., 2017; Refaely-Abramson et al., 2019). We apply periodic boundary conditions to all structures, and fully relax the internal coordinates for each one of the structures. For (i) and (ii), the calculations were carried out using the converged parameters of a (8 × 4 × 4) k-point grid and 600 bands in the summation to compute the dielectric matrix, as well as a kinetic energy cutoff of 80 Ry and a dielectric cutoff of 10 Ry. The calculation of the absorption spectrum was carried out via an interpolation to a finer (12 × 6 × 6) k-point grid. For (iii) and (iv), the calculations used a k-point grid of (8 × 4 × 1) and 1,000 bands in the summation to compute the dielectric matrix, as well as the same cutoffs as (i) and (ii). (i)-(iii) used a semiconductor screening treatment of the q → 0 limit and (iv) used a metallic screening, as implemented in the BerkeleyGW package (Deslippe et al., 2012). Furthermore, (iii) and (iv) used a slab truncation for the Coulomb interaction.

We begin with investigating the electronic and excitonic properties of the bulk C8-PDI system, using the experimentally-resolved crystal structure (Briseno et al., 2007), shown in Figure 2A (see full structural details in the SI). Figure 2B shows the computed GW quasiparticle bandstructure for this system (solid, red line). The three frontier electronic bands shown, valence (v₁), conduction (c₁), and second conduction (c₂), dominate the low-lying excitonic spectra we discuss in this work. The computed GW quasiparticle gap is 2.6 eV, with an expected large self-energy correction on top of the DFT (PBE) gap of 1.1 eV. The quasiparticle bands are well isolated, as typical in molecular crystals due to the molecular-like nature of the material (Refaely-Abramson et al., 2013; Kronik and Neaton, 2016); still, we note a significant band dispersion of ∼0.7 eV for the valence band between k-points X and Γ and 0.3 eV for the conduction bands between C and X, reflecting non-negligible intermolecular electronic coupling.

To explore the effect of intermolecular orientation, we computed different bulk structures through variation of relative intermolecular distance and the crystallographic parameters of the unit cell (see SI for full details). Our motivation to vary these parameters is twofold: first, to explore the effect of molecular packing on the quasiparticles and excitons; second, to look into the effect of changes in unit cell vectors, which are needed to achieve proper commensurateness between the unit cell of the C8-PDI crystal and that of the Au surface. In order to maintain minimal strain within the Au surface, we only vary the PDI cell parameters. For the case of a single C8-PDI molecule per cell, this results in a modified cell (“com. bulk” in Figure 1B), in which the intermolecular distance grows compared to the original structure. Here we explore the limiting case of an orthorhombic molecular crystal unit cell, although a monoclinic cell leads to very similar results (see SI). By allowing more molecules per cell, commensurate layers can be achieved with smaller atomic modifications, however, performing GW-BSE calculations on these large supercells in an accurate manner is highly computationally challenging. The associated structures and the computed GW bandstructures for each of these cells in its bulk form are given in Supplementary Table S1; Supplementary Figure S1.

Since the electronic and excitonic states stem from the molecular building blocks, the effect of such unit-cell modification is not trivial and depends on the level of the electronic and excitonic wavefunction localization. Upon structural reorganization from the “bulk” structure to the “com. bulk” cell that is commensurate with the Au surface, the intermolecular distance increases from 3.4 Å to 3.6 Å, and the cell angle γ (Figure 2A) changes from 82^◦ to 90^◦. As a result, the computed GW quasiparticle gap increases by 0.4 eV and the band dispersion decreases by 0.35 eV, as shown in Figure 2B. Two other bulk structures with intermediate variations reveal a gradual modification between these values (see SI).

The strong structural dependence of the quasiparticle band dispersion reveals a somewhat surprisingly large effect of the molecular packing and orientation on the quasiparticle and exciton picture. It stems from intermolecular interactions and screening, which vary significantly with intermolecular distances and relative orientation. Figure 2C shows the computed BSE absorption spectra of the “bulk” structure, with a low-energy singlet peak at 2.05 eV, in good agreement with the experimental value of ∼2.2 eV (Le et al., 2018; Felter et al., 2019). The absorption spectra along the three main polarization directions are shown in different colors, with the main optically-active dipole transitions along the a ̂ axis, namely through the π − π stacking direction, and the least active direction along the c axis. The low-lying triplet exciton is found at 1.4 eV, shown with a black line. As a comparison, the absorption computed for the commensurate bulk structure is shown as well (dashed black line). As shown in Table 1, comparing the “com. bulk” structure with the “bulk” structure, the lowest singlet excitation energy increases by 0.19 eV; however, the lowest triplet excitation energy increases by only 0.06 eV.

We further examine the effect of intermolecular coupling on the exciton nature by investigating the exciton dispersion, where the center-of-mass momentum Q represents the momentum difference between the hole and the excited electron, which is taken into account in the BSE, Eq. 2; (Qiu et al., 2015). Figure 2D shows the computed exciton bandstructure along the optically-active a ̂ direction (other crystal directions are shown in the SI). The exciton band shape varies along the reciprocal space due to indirect transitions occurring between the Γ and X points. At small exciton momentum Q, the lowest singlet exciton, S ₁, shows a parabolic behavior, with effective mass of 0.35m _e (for m _e the electron mass). Higher singlet excitons are also shown, revealing varying levels of localization. As expected, the low-lying triplet state has a larger effective mass of 2.7m _e , supporting its higher degree of real-space localization, as typical in organic molecular crystals (Rangel et al., 2016; Refaely-Abramson et al., 2017).

Figure 2E shows the exciton coefficients ( A v c k S in Eq. 2) for the Q = 0 transition as a function of quasiparticle momentum in the (K _x  , K _y ) crystal plane. The lowest singlet state shows high localization in momentum space at the C point ([0.5,0.5,0] in reciprocal space), suggesting spatial delocalization. On the contrary, for the lowest triplet state, the exciton coefficients are nearly uniform within the Brillouin zone, pointing to spatial localization. In the commensurate bulk structure, the same singlet excitons are more spread in reciprocal space, reflecting increased real-space localization. The triplet state is localized in both structures and experiences a smaller change due to the structural modification, compared to the singlet state. The sizable differences between singlet and triplet localization is coupled to the short-range exchange in C8-PDI, dominating intermolecular electron-hole coupling (Rohlfing and Louie, 2000; Rangel et al., 2016). Upon changes in intermolecular packing, induced exciton localization leads to enhanced exchange interactions. Our results thus reveal strong state- and structure-dependence of the exciton nature in bulk C8-PDI, as an outcome of the intermolecular interactions dominating the electron-hole coupling.

Next, we study the effect of crystal layering and surface adsorption on the excitation energies. As shown in Table 1, for the commensurate structure, the quasiparticle gap of a freestanding C8-PDI layer increases by ∼0.4 eV, while both the singlet and triplet exciton energies decrease by ∼0.1 eV, compared to the commensurate bulk structure. This results in a significant increase of the exciton binding energy, an expected result due to the strong dimensionality effect on the dielectric screening, as we further elaborate in the discussion section below. From a computational point of view, constructing a monolayer is motivated by its subsequent adsorption on Au substrate within a computationally tractable periodic cell. In the following we directly compare the charged and neutral excitations of the freestanding layer with the adsorbed one, to gain direct insight into the Au screening effect on exciton binding.

To capture the dielectric screening at the C8-PDI@Au interface, we compare two approaches: an explicit GW calculation of the entire interface, with the simulation cell presented in Figure 1; and a substrate screening GW approach (Liu et al., 2019), in which the polarizabilities of the two parts of the interface are computed separately and then combined in the interface cell, as we discuss in the methods section. Our motivation is to verify the validity of the substrate screening approximation for the C8-PDI@Au interface. If this is true, we can conclude that the effect of the Au substrate merely provides a dielectric media that renormalizes the quasiparticle and excitation energies within the C8-PDI molecular layer, rather than altering the nature of the quasiparticle orbitals and excitons of the C8-PDI via orbital hybridization.

Figure 3A shows the computed GW quasiparticle projected density of states (pDOS) of the full interface at high symmetry k-points onto C8-PDI (purple) and Au (orange) atomic orbitals. Dashed black lines represent the pDOS associated with the freestanding C8-PDI layer, which we align with the pDOS onto C8-PDI at the resonance corresponding to the highest occupied molecular orbital (at about −0.8 eV in Figure 3A). The quasiparticle gap associated with the molecular levels is reduced by 0.46 eV compared to the monolayer, due to the Au surface. The associated valence and conduction electronic charge distributions for the PDI-localized interface wavefunctions are well separated from the Au states, demonstrating the negligible interface hybridization in this system. The change in quasiparticle energies is a direct outcome of interface screening. Figure 3B compares the computed dielectric function ɛ ₀₀ (Eq. 1) of the commensurate bulk, layered, and interface systems for the head elements, namely G = G′ = 0. At large interaction distance q, corresponding to short-range interactions in real space, the dielectric function reaches the expected limit of ɛ = 1. At small q, however, the dielectric function increases significantly due to the Au screening, while remaining close to unity for the freestanding layer due to its reduced dimensionality.

The computed GW quasiparticle gap resulting from the substrate screening approach is 2.95 eV, which is in very good agreement with the direct calculation (2.94 eV). Additionally, in Supplementary Figure S5, we compare the diagonal elements of the non-interacting polarizability of the interface ( χ tot 0 ) computed using the two approaches, and find that they agree very well. We can infer from these results that there is negligible orbital hybridization at the C8-PDI@Au interface, and that the Au substrate simply provides a dielectric environment that effectively screens the molecular crystal. We expect that the substrate screening approach, based on the additivity of the non-interacting polarizability of the interface (Liu et al., 2019; Xuan et al., 2019), is applicable to other systems without significant orbital hybridization or covalent bonding.

We believe this conclusion could help us understand the physical interface where a C8-PDI monolayer in its pristine bulk lattice (rather than the artificial commensurate lattice as we did in this work) is adsorbed on Au (111). From a computational perspective, such a physical interface is incommensurate and is prohibitively expensive to calculate. However, based on the conclusion achieved in this section, we can infer what would happen if we were modelling the incommensurate physical interface: all the properties associated to the exciton wavefunctions of the C8-PDI bulk crystal as reported in Figure 2 will be qualitatively unchanged in the interface, thanks to the negligible orbital hybridization. The exciton binding energies will be renormalized due to the dielectric screening of the Au substrate.

Our results demonstrate the effect of intermolecular and surface screening on both quasiparticle and exciton properties via a step-by-step structural variation. Specifically, we look closely at the exciton nature as a function of molecular crystal packing and dielectric environment, as summarized in Figure 4. For the case of Frenkel-like molecular excitons, the exciton nature should stay roughly unchanged upon changes in crystal packing, as the environmental effect will be mainly manifested through an effective dielectric constant (Sharifzadeh et al., 2012; Cudazzo et al., 2013). Nonetheless, we find that the bulk C8-PDI excitons do have dispersion and are crystal-momentum dependent, reflecting non-negligible excitonic coupling beyond the Frenkel picture. The dispersion is more enhanced in singlet states compared to triplet states due to the exchange interaction, as observed before in related organic molecular crystals (Rangel et al., 2016; Refaely-Abramson et al., 2017; Rangel et al., 2018). On top of this effect, the dielectric screening induces further variations: upon crystal layering, the quasiparticle gap increases, and upon surface adsorption it renormalizes, while the exciton energies remain largely unchanged. This is a direct outcome of the non-local screening that dominates the quasiparticle energies, compared to local screening that dominates the electron-hole binding in the low-lying excitons. These dielectric effects are captured in the GW calculation through an explicit evaluation of the dielectric function (Hybertsen and Louie, 1986; Deslippe et al., 2012). As shown in Figure 3B, the dielectric functions of the various structures are similar at the short-range interaction regime and differ greatly at the long-range interaction regime.

Table 1 indicate the strong sensitivity of the charged quasiparticle excitation energies (as reflected in gaps) to the dielectric environment and the weak sensitivity of the neutral excitation energy (as reflected in optical transitions in the case of strong electron-hole coupling) to the environment (Kronik et al., 2012; Refaely-Abramson et al., 2015). The change in the optical gaps and the exciton binding energies between the C8-PDI monolayer and the C8-PDI@Au interface points to a significant decrease in the exciton binding energy upon formation of the interface with Au. This effect is expected, as the enhanced screening has a significant influence on the quasiparticle gap (Neaton et al., 2006; Egger et al., 2015; Liu et al., 2017). On the other hand, simple metal surfaces tend to have a much weaker influence on the optical excitations localized within the adsorbate (Spataru, 2013; Deilmann and Thygesen, 2019).

Importantly, the computational approach we employed here is not limited to cases of weak hybridization at the interface, and can be hence further used to investigate other types of PDI interfaces where the exciton separation mechanisms are expected to involve significant charge-transfer components. For example, triplet exciton transfer across C8-PDI and SiO₂ interfaces was recently suggested to be strongly coupled to changes in intermolecular interactions due to surface hybridization (Cotton et al., 2020). In addition, few recent studies explored the role of surface passivation in modifying and controlling the efficiency of exciton transfer from acene molecular crystals adsorbed on silicon substrates (Einzinger et al., 2019; Daiber et al., 2020). The relation between interface bonding, charge and exciton localization far from and at the junction, and the resulting energy-transfer efficiency at PDI-based junctions is intriguing, and the results we present here for the case of weak metal-organic coupling can be thought of as a computational test case at the weakly interacting regime. Our results thus demonstrate that the GW-BSE approach offers a reliable tool to explore interface effects on excitonic properties at organic-inorganic interfaces with various levels of interface hybridization.

To conclude, we studied the effect of structural modifications on quasiparticle and exciton nature in the C8-PDI molecular crystal and its interface with Au. We explored the excitonic properties of the bulk system in detail, and investigated the effect of crystal packing and interface dielectric screening, building a structural modification route from the bulk structure to a heterogeneous interface. Our results demonstrate that while the quasiparticle band gap undergoes significant variations upon the structural modification and surface adsorption, the optical gap is much less affected by them, leading to strong structural sensitivity of the exciton binding energies. Our methods allow us to quantify this effect, and relate it to the specifics of local and non-local structural modifications.