Excitonic effects in third harmonic generation: the case of carbon nanotubes and nanoribbons

TL;DR

This paper investigates how excitonic effects influence third harmonic generation in carbon nanotubes and nanoribbons, revealing that many-body interactions significantly alter the nonlinear optical response by reducing peak intensities and redistributing spectral weight.

Contribution

It introduces a new ab-initio methodology to analyze excitonic effects on nonlinear susceptibilities in low-dimensional nanostructures.

Findings

Excitonic effects significantly modify third-harmonic generation spectra.

Many-body effects reduce main peak intensity in nonlinear response.

Spectral weight is redistributed among excitonic resonances.

Abstract

Linear and nonlinear optical properties of low dimensional nanostructures have attracted a large interest in the scientific community as tools to probe the strong confinement of the electrons and for possible applications in optoelectronic devices. In particular it has been shown that the linear optical response of carbon nanotubes [Science 308, 838 (2005)] and graphene nanoribbons [Nat. Comm. 5, 4253 (2014)] is dominated by bounded electron-hole pairs, the excitons. The role of excitons in linear response has been widely studied, but still little is known on their effect on nonlinear susceptibilities. Using a recently developed methodology [Phys. Rev. B 88, 235113 (2013)] based on well-established ab-initio many-body perturbation theory approaches, we find that quasiparticle shifts and excitonic effects significantly modify the third-harmonic generation in carbon nanotubes and graphene nanoribbons. For both systems the net effect of many-body effects is to reduce the intensity of the main peak in the independent particle spectrum and redistribute the spectral weight among several excitonic resonances.

Equations

$${\epsilon }_{n\mathbf{k}}^{\text{QP}}={\epsilon }_{n\mathbf{k}}+Z_{n\mathbf{k}}\Delta {\Sigma }_{n\mathbf{k}}({\epsilon }_{n\mathbf{k}}).$$

$$Z_{n\mathbf{k}}=[1-\partial \Delta {\Sigma }_{n\mathbf{k}}(\omega )/\partial \omega {|}_{\omega ={\epsilon }_{n\mathbf{k}}}{]}^{-1},$$

$$\Delta {\Sigma }_{n\mathbf{k}}\equiv ⟨n\mathbf{k}|\Delta \hat{\Sigma }|n\mathbf{k}⟩,$$

$$\Delta \hat{\Sigma }=\hat{\Sigma }-{\hat{V}}^{\text{xc}},$$

$$i\hslash \frac{d}{dt}|v_{m\mathbf{k}}⟩$$

$$P_{\parallel }=-\frac{ef|\mathbf{a}|}{2\pi {\Omega }_c}\text{mbox}Imlog\prod _{\mathbf{k}}^{N_{\mathbf{k}}-1}\text{mbox}detS\left(\mathbf{k},\mathbf{k}+\mathbf{q}\right),$$

$$\mathbf{P}={\chi }^{(1)}\mathcal{E}+{\chi }^{(2)}{\mathcal{E}}^2+{\chi }^{(3)}{\mathcal{E}}^3+...$$

$$H_{\mathbf{k}}^{\text{IP}}\equiv H_{\mathbf{k}}^{\text{KS}},$$

$$H_{\mathbf{k}}^{\text{QP}}\equiv H_{\mathbf{k}}^{\text{KS}}+\Delta H_{\mathbf{k}}.$$

$$H_{\mathbf{k}}^{GW+\text{BSE}}\equiv H_{\mathbf{k}}^{\text{KS}}+\Delta H_{\mathbf{k}}+V_h(\mathbf{r})[\Delta \rho ]+{\Sigma }_{\text{SEX}}[\Delta \gamma ],$$

$$\Delta \rho \equiv \rho (\mathbf{r};t)-\rho (\mathbf{r};t=0)$$

$$\Delta \gamma \equiv \gamma (\mathbf{r},{\mathbf{r}}^{\prime };t)-\gamma (\mathbf{r},{\mathbf{r}}^{\prime };t=0)$$

Full text

I Introduction

Carbon nanotubes (CNTs) and graphene nanoribbons (GNRs) have remarkable electronic and optical properties due to their one-dimensional structure that combines solid-state characteristics with molecular dimensions. In these nanostructures light absorption produces strongly correlated electron-hole states in the form of excitons. Evidences of excitons have emerged from experimental studies of optical spectra and excited-state dynamics.Wang et al. (2005); Denk et al. (2014) The key role of excitons in the interpretation of the optical absorption of these materials has been confirmed by ab-initio computational studiesPrezzi et al. (2008); Spataru et al. (2004) based on many-body perturbation theory (MBPT).Strinati (1988); Aulbur et al. (1999) In fact the formation of strongly bounded excitons in GNRs has been theoretical predictedPrezzi et al. (2008) before the experimental measurement.Denk et al. (2014) In recent years, also the nonlinear optical response of these low dimensional structures has attracted a large attention, both from a fundamental and an applicative point of view. In particular due to their strong nonlinear response, one dimensional nanostructures find application as nanoantennas and optical switches.Chen et al. (2002); Tatsuura et al. (2003)

Experimentally, the absolute measure of nonlinear optical responses of these nanostructures is, however, not straightforward.111See for example Y. R. Shen “The Principles of Nonlinear Optics”, Wiley Interscience (2003) The first non-zero non-linear response function in carbon nanoribbon and non-chiral nanotubes is the third-order susceptibility.

In CNTs, only few measurements on third-order nonlinear susceptibility have been reported so far and most of them on the $\chi^{(3)}(-\omega;\omega,-\omega,\omega)$Maeda et al. (2005a); Liu et al. (1999); Ichida et al. (2009); Lauret et al. (2004), responsible for the intensity dependence of the refractive index. The one study on the third-harmonic-generation, $\chi^{(3)}(-3\omega;\omega,\omega,\omega)$,Stanciu et al. (2002); Nemilentsau et al. (2006a) explored the nonperturbative regime. Regarding the GNRs, to our knowledge there are no available experimental measures of the THG though recently measurements of the THG in graphene have been obtainedKumar et al. (2013); Hong et al. (2013); Säynätjoki et al. (2013). Interestingly these measurements found that graphene’s THG is of the order of $10^{-15}-10^{-16}m^{2}/V^{2}$ ($10^{-7}-10^{-8}$ esu), thus comparable to resonant THG in bulk materials.

Theoretically, state-of-the-art calculations of nonlinear optical responses of periodic systems neglect excitonic effects that are deemed essential to describe optical properties in these carbon nanostructures. For CNTs theoretical studies on THG—which mostly focused on radius and chirality dependence—have been performed at the independent particle approximationMargulis and Sizikova (1998); Zhang et al. (2006); Rezania and Daneshfar (2012). The only work which includes many-body effects to our knowledge is from Lacivita and coworkers Lacivita et al. (2016) that computed the static second hyperpolarizability of CNTs within the coupled Kohn-Sham (KS) equations formalism.222Note that none of the mentioned studied included temperature effects, except Ref. Rezania and Daneshfar, 2012 which includes temperature effects by means of Green’s function theory. For GNRs the nonlinear properties have been addressed in few studiesWang and Andersen (2016); Zhang et al. (2014); Karamanis and Pouchan (2013); Ahmadia and Asgarib (2011) that focused mainly on the possibility of enhancing the nonlinear response by engineering the ribbon edges.Zhang et al. (2014); Karamanis and Pouchan (2013) The works in Refs.Zhang et al., 2014; Karamanis and Pouchan, 2013 are based on time-dependent density functional theory with an hybrid approximation for the exchange-correlation functional, which approximately accounts for excitonic effects. Those studies however have been performed on finite length GNRs.

Here we present a theoretical/computational combined study of THG in CNTs and GNRs (Fig. 1) that includes excitonic effects. Specifically, we use an ab-initio approach based on many-body perturbation theory (MBPT) which has been shown to provide accurate results for both linear and nonlinear optical properties of periodic systems.Strinati (1988); Aulbur et al. (1999); Attaccalite et al. (2015); Grüning and Attaccalite (2014)

With respect to the works mentioned above which include correlation beyond the independent particle approximation, our study addresses infinite periodic nanostructures (at difference with Ref.Zhang et al., 2014; Karamanis and Pouchan, 2013) and frequency dependent third order response (at difference with Ref. Lacivita et al., 2016). Furthermore our MBPT approach does not depend on any semi-empirical parameter at difference with hybrid density-functional based approaches used for previous works. In the latter, the semi-empirical mixing parameter eventually determines both the fundamental band gap and the amount of dielectric screening.

The purpose of this work is threefold: first, to provide an accurate theoretical estimate of the THG in a small semiconducting CNT and GNR; second and more importantly—as we can switch on and off excitonic and many-body effects in the (effective) Hamiltonian of the electron system—to evaluate how those effects affect the THG of the material; finally, to provide a benchmark to assess the reliability and accuracy of calculations at the independent-particle level which neglects these effects. The paper is organized as follows: in section II we summarize the computational methods employed in the calculation of the electronic structure and nonlinear response; in section III we present results for both the CNTs and GNRs. We consider here only centrosymmetric systems for which the THG is the first non-negligible nonlinear response.333A centrosymmetric system has an inversion center and thus the second-harmonic generation is zero within the electric dipole approximation. When multipole corrections are considered the second-harmonic generation is nonzero, but still very small.

II Computational Methods

The nanotubes and nanoribbons atomic structures have been generated from the ideal graphene with bonding length of $1.421$Å. Subsequently atomic positions have been optimized by means of Density Functional Theory (DFT), using the local density approximation for the exchange correlation functional.Perdew and Zunger (1981); Perdew and Wang (1992)

All DFT calculations have been performed with the Quantum Espresso codeGiannozzi et al. (2009) where the wave functions are expanded in plane waves with a cutoff of $60$ Ry and the effect of core electrons is simulated by norm-conserving pseudopotentials.Troullier and Martins (1991) We used a $1\times 1\times 22$ ${\bf k}$-points grid to converge the density in the CNT and GNR (both structures are oriented along the $z$ axis). Valence and conduction orbitals that enter in the Green’s function theory are obtained from the diagonalization of the KS eigensystem. The KS eigensolutions $\{\varepsilon_{n{\bf k}};|n{\bf k}\rangle\}$ correspond to the energies and Bloch wave functions (with ${\bf k}$ the crystal wave vector and $n$ the band index) of the independent particle system that reproduces the electronic density of the system under study.

In order to simulate isolated nanotubes and nanoribbons we used a supercell approach with a tube-tube distance larger than the nanotube diameter and a distance between the ribbons of 16 a.u. in the perpendicular direction and larger than the ribbon size in the plane.

II.1 Quasiparticle band structure

Starting from the KS eigenvalues and eigenfunctions, we obtained the quasiparticle (QP) band structure by means of MBPT within the so-called $GW$ approximation. Specifically, we use non-self consistent $GW$, often denoted as $G_{0}W_{0}$, in which the screened Coulomb potential, $W$, and the Green’s function, $G$, are built from the KS eigensolutions $\{\varepsilon_{n{\bf k}};|n{\bf k}\rangle\}$ and the quasiparticle energies are obtained from:

[TABLE]

In Eq. 1

[TABLE]

is the renormalization factor and

[TABLE]

where

[TABLE]

is the difference between $\Sigma=GW$, the $GW$ self-energy, and $V^{\text{xc}}$, the exchange-correlation functional used in the KS calculation.Aulbur et al. (1999) The screened Coulomb potential $W$ has been evaluated within the random-phase approximation. In the $GW$ approach we used the Godby-Needs plasmon-pole model to approximate the dynamical behavior of $W$,Rojas et al. (1995) while in the Bethe-Salpeter equation (BSE) framework described in the next subsection we use the static approximation.Strinati (1988) In the $GW$ (and BSE) calculations we used a truncated Coulomb potential to reduce the interaction between the periodic replica.Rozzi et al. (2006) The Green’s function and the self-energy that appear in the $GW$ calculations are expanded in the basis of the KS eigensolutions: we used a 40 ${\bf k}$-points and 320 bands for the CNT and 60 ${\bf k}$-points and 200 bands for the GNR.

II.2 Linear and nonlinear response functions

Linear and nonlinear optical properties are obtained by means of a real-time implementation of the BSE.Strinati (1988); Schmidt et al. (2003) We used an effective Hamiltonian that includes electron-hole interaction through a screened exchange interaction and the coupling between electrons and the external field is described by means of modern theory of polarization.Souza et al. (2004); Attaccalite and Grüning (2013) This formulation allows us to correctly describe response functions beyond the linear one.Attaccalite and Grüning (2013)

Specifically we solve a set of coupled one-particle effective time-dependent Schrödinger equations:

[TABLE]

where $|v_{m{\bf k}}\rangle$ is the periodic part of the Bloch states that determine the system polarization Souza et al. (2004) as discussed below. In the r.h.s. of Eq. 2, $H^{\text{sys}}_{\mathbf{k}}$ is the system Hamiltonian—which is discussed later in this Section; the second term, $\boldsymbol{\cal E}\cdot\tilde{\partial}_{\bf k}$, describes the coupling with the external field $\boldsymbol{\cal E}$ in the dipole approximation. As we imposed Born-von Kármán periodic boundary conditions, the coupling takes the form of a ${\bf k}$-derivative operator $\tilde{\partial}_{\bf k}$. The tilde indicates that the operator is “gauge covariant” and guarantees that the solutions of Eq. 2 are invariant under unitary rotations among occupied states at ${\bf k}$ (see Ref. Souza et al., 2004 for a discussion on this point).

From $|v_{m{\bf k}}\rangle$, the time-dependent polarization of the system $P_{\parallel}$ along the lattice vector $\mathbf{a}$ is calculated as

[TABLE]

where $S({\bf k},{\bf k}+\mathbf{q})$ is the overlap matrix between the valence states $|v_{n{\bf k}}\rangle$ and $|v_{m{\bf k}+{\bf q}}\rangle$. Furthermore, $\Omega_{c}$ is the unit cell volume, $f$ is the spin degeneracy, $N_{{\bf k}}$ is the number of ${\bf k}$ points along the polarization direction, and $\mathbf{q}=2\pi/(N_{{\bf k}}{\mathbf{a}})$. Finally, the third harmonic coefficient is extracted from the expansion of the polarization in the laser field $\boldsymbol{\cal E}$ power series:

[TABLE]

as detailed in Ref. Attaccalite and Grüning, 2013.

Notice that our approach to calculate nonlinear susceptibilities does not work for metallic systems or systems with a very small gap, since it uses the polarization [Eq.3] as fundamental quantity. For metallic systems other approaches based on the electron current-density should be used instead.Springborg and Kirtman (2008)

In Eq. 2 the model Hamiltonian chosen for $H^{\text{sys}}_{\mathbf{k}}$ determines the level of approximation in the description of correlation effects in the linear and nonlinear spectra. In this work we use different models for the system Hamiltonian: (i) the independent-particle (IP) model,

[TABLE]

where $H^{\text{KS}}_{\bf k}$ is the unperturbed KS Hamiltonian and we consider the KS system as a system of independent particles. The QP model (ii),

[TABLE]

where a scissor operator $\Delta H_{\bf k}$, estimated from the MBPT [Eq. 1], has been applied to the KS eigenvalues. Finally (iii) the full $GW$+BSE model reads,

[TABLE]

where

[TABLE]

is the variation of the electronic density and

[TABLE]

is the variation of the density matrix induced by the external field $\boldsymbol{\cal E}$.

In Eq. (7) the term $V_{h}(\mathbf{r})[\Delta\rho]$ is the time-dependent Hartree potentialAttaccalite et al. (2011) and is responsible for the local-field effectsAdler (1962) originating from system inhomogeneities. In the same equation, $\Sigma_{\text{SEX}}$ is the screened-exchange self-energy that accounts for the electron-hole interaction,Strinati (1988) and is given by the convolution between the screened interaction $W$ and $\Delta\gamma$. In the limit of small perturbation Eq. 7 reproduces the optical absorption calculated with the standard $GW$ + BSE approachStrinati (1988), as shown both analytically and numerically in Ref. Attaccalite et al., 2011.

We calculate the $|\chi^{(3)}_{zzzz}(-3\omega;\omega;\omega;\omega)|$, i.e. the magnitude of the third order nonlinear susceptibility at a frequency of $3\omega$ along the $z$ axis (the orientation axis of nanotubes and the ribbons) when a monochromatic electric field of frequency $\omega$ is applied along $z$. This quantity is obtained by integrating Eq. (2) numerically for a time interval of $120$ fs using the numerical approach described in Ref. Souza et al., 2004 (originally taken from Ref. Koonin and C., 2008) with a time step of $\Delta t=0.01$ fs, which guarantees for numerically stable and sufficiently accurate simulations. $|\chi^{(3)}_{zzzz}|$ is finally extracted from the total polarization after 100 fs as described in Ref. Attaccalite and Grüning, 2013. A dephasing term with a time $\tau=8.78$ fs is introduced in Eq. (2) in order to simulate a finite broadening of about $0.15$ eV.Attaccalite and Grüning (2013) Finally note that our calculations provide the volume third-order nonlinear susceptibility, but because of the use of the supercell, two of the dimensions are not physical. We then rescale by the effective physical dimensions of the systems: for both systems we use $0.335$ nm as effective thickness and for the GNRs we consider as effective width the width of the ribbon plus $0.2$ nm.

III Results

III.1 Carbon nanotube

A single-walled CNT is formed by rolling a sheet of graphene into a cylinder along an $(m,n)$ lattice vector in the graphene plane. These two indexes determine the diameter and chirality, which are key parameters of a nanotube. Depending on the chirality (the chiral angle between hexagons and the main CNT axis), nanotubes can be either metals or semiconductors, with band gaps that may vary between few meV to the eV, even if they have nearly identical diameters.Smalley et al. (2003) We consider here semiconducting zig-zag CNTs only.

We analyze the effect of correlation in THG within the $GW$+BSE approach as detailed in Sec. II.2. Due to the computational cost we limit the calculation to the (10,0) CNT. We expect the analysis to be valid for larger zig-zag CNTs.

In Fig. 3 we report the THG magnitude $|\chi^{(3)}_{zzzz}(\omega)|$ at different levels of approximation [panels (a) and (c)] and the corresponding optical absorption [panels (b) and (d)] so to identify the resonances in the THG spectrum from the comparison.

At the IP level [panel (b) of Fig. 3] the absorption spectrum is dominated by a large peak at about $0.9$ eV corresponding to the transition between the highest valence and lowest conduction bands at $\Gamma$ [see band structure, panel (a) of Fig. 2]. Transitions between second-highest valence and the lowest conductions bands give rise to the shoulder around $1.3$ eV, that is visible also in the THG spectra. Peaks at higher energies are originated from transitions from the third-highest valence bands.

Turning to the THG, at the level of the IP model [panel (a) of Fig. 3] we observe at about $0.3$ eV and $0.4$ eV three photon-resonances corresponding to the main peaks in the absorption spectrum. Notice that in the THG the intensity of the peaks is reversed respect to the linear optics. At $0.6$ eV we found a smaller three photon-resonance peak that stems from the transition from the third-highest valence band. One-photon resonance peaks are barely visible on this scale. Our results are consistent with those obtained by Nemilentsau et al.Nemilentsau et al. (2006b) and Xu et al. Xu and Xiong (2004). With respect to Nemilentsau and coworkers our calculated linear and nonlinear response functions have a richer structure as we use a full ab-initio band structure rather than the two-bands model employed in their work.

Inclusion of QP corrections rigidly blue-shifts the absorption spectrum by approximately $1.6$ eV (not shown). When one includes as well the electron-hole interaction [$GW$+BSE model in Fig. 3 panel (d)] the spectrum is red-shifted with respect to the QP model by $1.5$ eV (excitonic binding energy). As a result the resonance peak is blue-shifted by about $0.1$ eV with respect to the IP level. More remarkably the spectrum is dominated by a strong excitonic peak at $1.1$ eV which almost doubles in intensity the van Hove singularity in the IP spectrum.

In the THG spectrum the inclusion of QP corrections blue-shifts it by about $0.5$ eV and reduces the spectral intensity by almost one order of magnitude due to the sum rules constraints.Saarinen (2002) Similarly to what is observed for the absorption spectrum the inclusion of electron-hole interaction [panel (c) of Fig. 3] red-shifts the spectrum and enhances the spectral intensity with respect to the QP model. The spectral enhancement is however not uniform: the intensity of the main peak is doubled—though it is still about 1/3 of the main peak intensity in the IP model. On the other hand peaks at higher frequencies aquire a comparable intensity as the main peak. The one-photon resonance with the main peak remains very weak.

Notice that the $GW$+BSE calculations include also the so-called local field effects, namely the response of the time-dependent Hartree term to the external field. These corrections are large for inhomegenous systems, such as isolated molecules or localized orbitals and are exactly zero for an homogeneous electron gas. In the linear response of one-dimensional systems (e.g. CNTs) local field effects are negligible along the periodic direction and large for the perpendicular directions. We found that in the periodic direction the local fields effects are negligible for the THG as well.

Finally, at the IP level only we study how the THG varies with the size of the CNT [Fig.4 panel (a)]. Consistently with what found previouslyMargulis and Sizikova (1998); Zhang et al. (2006); Xu and Xiong (2004) the THG increases superlinearly with the CNT radius. The resonant peaks shift to lower energy as a consequence of the gap shrinking for larger CNTs. The two peaks however do not shift of the same amount and they become more apart in larger CNTs. We found a factor of about 1.4 between the peaks of the (13,0) and (16,0) CNTs in agreement with Refs. Margulis and Sizikova (1998); Xu and Xiong (2004) which predicted a factor of about 1.5. The peaks magnitude is in agreement with that of Margulis and co-workers Margulis and Sizikova (1998) while other works have reported larger magnitude especially Xu and coworkers whose reported THG at resonances is larger by almost two orders of magnitude. Regarding the spectral shape, it is clear that the two-band model employed in previous works is not sufficient and at least the second-highest valence band needs to be included.

III.2 Armchair nanoribbon

Armchair nanoribbons are divided in three distinct families depending on the ribbon width, namely $N=3p$, $N=3p+1$, and $N=3p+2$, N being the number of dimer lines along the width, with $p$ a positive integer. Within each family, the fundamental band gap decreases with increasing ribbon width.Son et al. (2006) We consider here the semiconducting $N=3p$, $N=3p+1$ families only.

We analyze the many-body effects on the THG [panels (a) and (c) of Fig.5] and on the absorption spectrum [panels (b) and (d) of the same figure] for the armchair GNR with $N=9$ (9-AGNR). Within the IP level of approximation [Eq.(5)] the optical spectrum shows the characteristic 1D van Hove singularity at about $0.8$ eV [panels (b) of Fig.5] resulting from transitions at the $\Gamma$ point [see band structure, panel (c) of Fig. 2] and two shoulders at about $1.0$ eV and $1.2$ eV. When the quasiparticle corrections and the electron-hole interaction are turned on [$GW$+BSE model in Eq. (7)] the absorption spectrum below $2$ eV [panels (d) of Fig.5] presents a single excitonic peak well below the onset of the continuum [which is shifted of $1.2$ eV by quasiparticle corrections, see vertical lines in panels (b) and (d)] with an exceedingly large excitonic binding energy for semiconducting materials ($1$ eV). The main peak is enhanced by about 40% by the electron-hole interaction.Prezzi et al. (2011) To sum up, the inclusion of the quasiparticle corrections and electron-hole interaction modifies the absorption line-shape narrowing the main peak. These effects are indeed known to be important for a qualitative and quantitative predictions of the optical spectra.Prezzi et al. (2008)

The effect of correlation on the THG intensity is even more dramatic and not entirely predictable from what observed for absorption. Panel (a) of Fig. 5 shows the third-harmonic intensity obtained from the IP-model. The peaks observed in the $\chi^{(3)}_{zzzz}$ correspond to three photon-resonances of the van-Hove singularities observed in the absorption and in fact the spectral shape roughly resembles that of the absorption [panel (b) of Fig 5]. When we apply QP correction to the KS band structure [QP model in Eq. (6)], the spectrum (same panel) is shifted to higher energies by about $0.35$ eV which is roughly 1/3 of the quasiparticle correction to the fundamental band gap (see also Ref. Prezzi et al., 2011). Since by construction $\chi^{(3)}$ contains the product of terms with poles at different energies the $GW$ correction effect is not a simple energy shift. The spectral shape is modified too. More strikingly there is a substantial reduction of the peak intensity (about 75%) as one can expect from sum rule constraints.Saarinen (2002)

When excitonic effects are turned on ($GW$+BSE model, panel (c) of Fig 5), the spectra is red-shifted with respect to the QP model. As observed for the absorption, the cancellation between the QP corrections and the exciton binding energy is partial and the main peak in the spectrum is slightly blue-shifted ($0.05$ eV) with respect to the IP model. More importantly the spectral weight is redistributed among few excitonic three-photon resonances and the spectral shape is distinctly different from the one obtained from the IP model. The main peak is in this case significantly broadened.

In analogy to the CNTs, we study how the THG varies as a function of the ribbon size (width) for $N=3p+1$, with $p=2,3,4$ [Fig. 4b panel (b)] and $N=3p$ with $p=3,4,5$ [Fig. 4b panel (c)]. Similarly to what is found for the CNTs the intensity increases with the size. This is a consequence of the bandgap decreasing when the width increases consistently with the findings in Ref. Ahmadia and Asgarib, 2011. We also observe that the THG of the $N=3p+1$ family, which has the largest bandgaps, is smaller by about one-two orders of magnitude than the THG of the $N=3p$ family, which is related to the suppresion of the intrabands terms in the THG. As well in the $N=3p+1$ the increase of the THG with the size is more pronounced since the bandgap decreases rapidly with $p$. For both families the increase is not uniform, but is larger for the first peak. In particular for the $N=3p$ family the first peak increases by a factor $2$ passing from $N=9$ to $N=12$ and by a factor $1.5$ passing from $N=12$ to $N=15$. In both cases it is significantly blue-shifted. The shoulder at $0.4$ eV in the 9-AGNR instead is only slightly blue-shifted and increased in intensity in the larger GNRs.

Finally our calculations predict a THG of the order of $10^{-7}esu$ at resonances for small GNRs of the same order of magnitude as the THG in graphene.Kumar et al. (2013); Hong et al. (2013); Säynätjoki et al. (2013)

IV Discussion and Conclusions

Using state-of-the-art ab-initio MBPT approaches ($GW$+BSE) we have studied the THG of two paradigmatic one-dimensional semiconducting nanostructures, a carbon nanotube and a graphene nanoribbon. By comparing the results from the QP and $GW$+BSE model Hamiltonians with the simple IP model Hamiltonian we were able to single out the effect of the electron-electron and electron-hole interaction in the THG of these two systems. For both systems the inclusion of many-body effects modifies significantly the THG: first, the intensity of the main peak is reduced by about a factor 3; second, additional structures have a significant spectral weight and as a result the spectrum has a much richer structure and covers a much larger range of frequencies than its counterpart at the IP level. These results indicate that it is important to include many-body effects for a qualitative and quantitative description of nonlinear optical properties of nanostructures.

Our results agree with what observed by Lacivita and coworkers in a recent study on static hyperpolarizabilities of carbon nanotubes using a coupled perturbed KS approach. Lacivita et al. (2016) Their calculations show a gradual suppression of the longitudinal hyperpolarizabilities of the carbon nanotube when increasing the non-local exchange contribution (Hartree-Fock) to the exchange-correlation potential from [math] to $100$%. In fact the inclusion of non-local exchange in the DFT functional opens the fundamental band gap—similarly to the QP corrections in our work—which reduces the third-harmonic intensity. Lacivita and coworkers further observed that the reduction due to the band gap opening is only partially compensated by the inclusion of electron-hole interaction in the response function, in agreement with our findings. Note that in experiments it is difficult, if not impossible to isolate the effects contributing to the THG. For example one could in principle compare the $\chi^{(3)}_{zzzz}$ measured in an isolated CNT with that in CNTs bundles—increasing the dielectric screening and thus ‘suppressing’ the many-body effects. However other effects, as for instance the variation of the phase relaxation time with the environment, would make the interpretation of these experiments not straighforward.Maeda et al. (2005b)

Remarkably, including QP corrections and electron-hole interaction has a very different effect on nanostructures optical absorption spectra: QP corrections usually correspond to a rigid blue-shift of the spectra, while electron-hole interaction produces a red-shift and a substantial enhancement of the intensity of the first peak, which is usually more intense than the corresponding van Hove singularity at the IP level. Calculations of the second-harmonic generationAttaccalite et al. (2015); Pedersen and Cornean (2007) showed that electron-hole interaction enhances by a 20% to 200% the intensity of the main spectral features when compared with spectra calculated at IP level. In the case of SHG in fact addition of QP corrections shifts the spectrum and reduces the overall spectral intensity (as in the case of THG), however this reduction is overcompensated by the enhancement from the electron-hole interaction.

We argue that the different behavior of the THG, SHG and absorption follows from sum rule constraints.Saarinen (2002) On the other hand dimensionality could also play an important role since in one-dimensional systems both QP corrections and exciton binding energy are particularly large due to geometrical confinement and poor screening. At this regard it would be of interest to systematically study the THG for systems of different dimensionality. At present this study is hindered by the computational cost of solving the time-dependent Bethe-Salpeter equationAttaccalite et al. (2011) for systems with a large number of ${\bf k}$-grid points, as it is the case of bulk semiconductors. Other approaches, based for example on the generalization of time-dependent density functional theoryGrüning et al. (2016); Luppi et al. (2010) or extension of coupled Kohn-Sham equations to dynamical electric fieldsKirtman et al. (2015) could make this kind of studies affordable, once accuracy issues have been addressed.

Finally, we have also studied the dependence of the THG on the size of the nanostructures (i.e. the radius of the CNT and the width of the GNR). Consistently with previous studies we have found that THG increases by increasing the radius and width of the nanostructures mainly because of the shrinking of the bandgap. Our results show as well that the two-band model that has been used in previous studies does not capture the details of the THG spectrum and more bands close to the Fermi energy need to be taken into account.

V Acknowledgments

Computing time has been provided by the national GENGI-IDRIS supercomputing centers at Orsay under contract $n^{o}$ i2012096655. CA and MG thank the EUspec COST Action for Short Term Scientific Missions that allowed to initiate and carry out this project. CA acknowledges support from the “Fondation Aix-Marseille Université”.