Diabatic Mechanisms of Higher-Order Harmonic Generation in Solid-State Materials under High-Intensity Electric Fields

We theoretically investigate mechanisms of higher-order harmonic generation in solid-state materials under a high-intensity ac electric field. A new theoretical framework presented in this Letter holds the legitimacy of Bloch theorem even under the influence of the high-intensity electric field and provides an exact treatment of the diabatic processes of Bloch electrons. Utilizing this framework, we first discovered that the diabatic processes, namely, ac Zener tunneling and semimetallization of semiconductors, are key factors for nonperturbative mechanisms of HHG. These mechanisms are classified by the field intensity and could be understood by an extended simple man model based on an analogy between tunnel ionization in gaseous media and Zener tunneling in semiconductors. These conclusions would stimulate the universal understanding of HHG mechanisms in both atomic and solid cases.

spectra [4,5], whose definition is a threshold between a constant intensity region (plateau region) and strong decay region, while, in atomic cases, that is determined by the relation with ionization energy of a single atom. Thus, HHG in solid-state materials is expected to have a different mechanism from that in atomic cases, and a clear understanding of such mechanism will provide a possibility for opening new research fields of high-intensity optical technology.
For the clear understanding of mechanisms of HHG in solid-state materials, we should develop an exact treatment on dynamics of the collective electrons under a high-intensity electric field. In solid-state materials where crystalline structure ensures a periodicity of atomic potential, Bloch's theorem enables us to reduce the collective behavior of electrons to single quasiparticle, which is called Bloch electron [9]. Based on this scheme, we expect an exact investigation of dynamics of Bloch electrons under the high-intensity AC electric field will unravel mechanisms of HHG. The highintensity AC electric field causes a concurrence of excitation and transport processes of Bloch electrons, and therefore, simultaneous treatment of these processes is necessary whose considerations enable discussions on diabatic processes of Bloch electrons such as Zener tunneling [10]. However, in spite of the simplicity of this problem, previous studies [5,[11][12][13][14][15][16][17] have not provided reasonable treatments yet because of their invalid assumptions in constructing frameworks which would be summarized as follows.
In constructing the frameworks, we should take care of legitimacy of Bloch's theorem under an influence of the high-intensity electric field. If the electric field is introduced in the form of a nonperiodic scalar potential = − ( ) • , the periodicity of atomic potential would be damaged and justification of Bloch's theorem is no longer expected to be adequate. Here, is the electron charge, the external electric field, and the position of the electrons. In this case, Bloch functions cannot be regarded as a complete basis set and concepts of band structures and Bloch electrons become invalid [18,19]. These invalidities would greatly affected theoretical treatments of diabatic processes, as referred in the previous work [18,19], and their influence on perspective of HHG becomes emphasized especially when employing the high-intensity AC electric field. In spite of these failures, some previous works [5,11,12], where transport processes are introduced by the scalar potential, constructed their theoretical frameworks based on assumptions of the completeness of Bloch functions as well as the concepts of band structures and Bloch electrons. For correct perspective of HHG in Bloch electron picture, it is essential to hold the periodicity of atomic potential even under an influence of the high-intensity electric field.
In the simultaneous treatments of excitation and transport processes of Bloch electrons, we also should take care of justification of classical kinetic equation ℏ = − , where ℏ and is reduced Planck constant and Bloch wavevector. This equation is usually assumed in the treatments of transport processes of Bloch electrons and can be employed in theoretical frameworks by simple replacement of Bloch wavevector with − ( /ℏ ) , where the electric field is described as = − (1/ ) ( ) . Here, and are the velocity of light and a vector potential. For the simultaneous treatments of excitation and transport processes, the previous works [13][14][15][16] utilized this replacement in order to introduce transport processes of Bloch electrons, while excitation processes is involved in usual dipole transition. However, as indicated by several researchers [20,21], this replacement is true only for single-band processes and should be violated by excitation processes such as Rabi flopping [22]. In the present situation, Bloch electron states are always described as a superposition of ground and excited states, and consequently, we should conclude that the previous works could not provide an accurate perspective of HHG due to an inaccurate expression of nonperturbative interaction of Bloch electrons. Moreover, the physical interpretation of HHG mechanism based on Bloch oscillations [4,5,17], whose concept is true only for single-band processes, is no longer expected to be adequate. For correct perspective, the theoretical frameworks should be constructed without assuming the conventional replacement of Bloch wavevector.
In this paper, avoiding invalid assumptions employed in the previous works [5,[11][12][13][14][15][16][17], we clarify non-perturbative mechanisms of HHG in solid-state materials. To treat the competition between excitation and transport processes of Bloch electrons accurately, our theoretical framework started from the Hamiltonian = (1/2 6 )( − / ) 8 + Σ ; − ; , where Bloch's theorem has not been applied yet. Here, 6 is the electron mass, the momentum of bare electrons, and − ; the periodic core potential of atoms located at ; . In this model, we employed a homogeneous electric field in the form of a vector potential rather than a scalar one, as has been proposed by Krieger and Iafrate in the case of DC electric field [19]. An advantage of this treatment is that the Hamiltonian maintains its periodicity even including an influence of the high-intensity electric field. This ensures the validity of Bloch's theorem under the high-intensity electric field, and consequently, Bloch functions are expected to be a complete basis set and concepts of band structures and Bloch electrons become adequate. Legitimacy of Bloch's theorem including the high-intensity electric field enables us to expect that there will be temporally changed band structures linked to an AC electric field, whereas Bloch wavevector will not change and be a kinetic constant [21]. We should take care of the Hamiltonian not always leading to the classical kinetic equation for Bloch wavevector ℏ = − , especially when treating excitation and transport processes together. Based on this scheme, we firstly discovered diabatic processes, namely, AC Zener tunneling and semimetallization of semiconductors, determine characteristics of HHG in solid-state materials. These considerations would provide new understandings of the fundamental mechanisms of HHG in semiconductors.
In constructing our theoretical framework, we suppose two-dimensional semiconductors which are referred in the recent experiment [5], while conclusions presented in this paper have no dependence on the dimensionality. Starting from this assumption and focusing only on conduction and valence bands, we could introduce a Hamiltonian described by = 6 + ? [23], where Here, R = ℏ 8 8 /2 R ( = , ℎ) are the kinetic energies of electrons and holes, B is the band-gap energy, (ℎ ) and C (ℎ C ) are annihilation and creation operators of electrons (holes), is the argument of Bloch wavevector , and Ω I = Ω I6 exp − − 6 8 / 8 cos 6 is the Rabi frequency [28] where the intensity of the electric field is renormalized in Ω I6 . Throughout this paper, we will fix the parameters of the incident electric field as t 6 = 12π/ω 6 and τ = 4π/ω 6 . The first term on the right hand side in Eq. (2) indicates the intraband transition in which Bloch wavevector is a kinetic constant and it can be renormalized in the single-particle energy R , where R = R + B /2 + ℏΩ I cos . These modifications caused by an AC electric field mean a temporal variation of the band structure, which is consistent with the above discussion. The variations can be derived from temporal changes of diagonal matrix elements originating from the light-matter interaction that is usually ignored in nonlinear optics [22,28]. Moreover, comparing with Landau-Zener model [29], we can interpret the Hamiltonian to include diabatic processes such as above-threshold ionization and Zener tunneling [10]. The second term of the right side in Eq. (2) indicates a dipole transition causing multiphoton absorption and Rabi flopping, which are nonlinear optical phenomena [22,28]. The factors sin and cos are derived from the form factor reflecting microscopic information on the arrangement of atoms in crystals.
Considering the above-mentioned Hamiltonian, we can derive time evolution equations for populations R = C and polarization = ℎ F C with Bloch wavevector as Here, c and g are the transverse and longitudinal relaxation constants, and throughout this study, they will be fixed to c = 0.1 6 and g = 0.01 6 , respectively. These values are estimated by the recent experimental situation [5] where the incident frequency are several tens of terahertz, while the dephasing and the energy relaxation time are regarded as a few tens of femtoseconds and a few picoseconds, respectively. The numerical solutions of these equations give the time evolutions of distributions of the carrier densities and polarization in two-dimensional space. The distributions arise from a dipole transition under the anisotropic modification of the band structure characterized by the factor ℏΩ I cos , and consequently, they show anisotropic behavior in space linked to the direction of the electric field [23]. Therefore, we expect that the intraband processes lead to carrier transport (wave packet dynamics), and the time evolution of the current can be calculated using the definition Accordingly, we can derive higher-order harmonic spectra from the definition = 8 , where is the Fourier transform of , and discuss the dependence of higher-order harmonic spectra on the intensity of the incident electric field corresponding to the Rabi frequency Ω I6 .
The numerical results show characteristics of HHG are changed depending on the Rabi frequency. In this study, we always assume a situation where ω 6 ≪ B /ℏ. At the beginning, in multiphoton absorption regime where a perturbative treatment of light-matter interaction is assured (FIG. 1(a)), we could identify the well-known characteristics in HHG spectra (FIG. 2(a)), i.e., the conventional relation of nonlinear optics n ∝ n 8 ∝ 6 8n [22] and strongly generation of only odd-order harmonics which is explained by HHG emitted in each half-cycle period of the incident electric field [1]. With increasing the Rabi frequency, unconventional non-perturbative mechanisms of HHG should be emphasized, as shown in FIG. 1(b) and FIG. 1(c), where excitation processes are dominated by Zener tunneling and semimetallization processes, respectively. To understand these regimes easily, in these figures, we assume uniform modifications of band structures described by k R = k R + B /2 + ℏΩ I , while the anisotropic factor cos only causes dynamics of the wave packet in space. The threshold of the field intensity between multiphoton absorption regime and AC Zener regime can be found as the boundary when a perturbative theory becomes inadequate, which is derived as Ω I6 = 0.5ω 6 [30], while that between AC Zener regime and semimetal regime can be estimated as the boundary when the sum of the renormalized single-particle energy is zero: q6 A + q6 D = 0 , i.e., B /2ℏ = Ω r6 . In the following, we will discuss characteristics and mechanisms of HHG in each regime. FIG. 2(b) shows the higher-order harmonic spectra in AC Zener regime, where the band-gap energies are B = 10ℏ 6 (blue line) and B = 5ℏ 6 (red line) in the case of Ω I6 = 2 6 . We find that the spectra are divided into plateau and decay regions, and the plateau region becomes broader with increasing band-gap energy. This mechanism can be understood by an analogy between Zener tunneling and tunnel ionization processes in gaseous media [1]. The threshold between the plateau and decay regions is called cutoff energy, which corresponds to maximum energy of electron-hole pair accelerated by electric field after the pair excitation. We will consider a simple man model [1,31] in semiconductors, as shown in FIG. 3(a), and estimate the cutoff energy s as s = 2× B /2 + ℏΩ r6 + 3.2 w = B + 2×3.2 w , where B = B − 2ℏΩ r6 is renormalized band-gap energy modulated by the effect of intraband transition and w is a ponderomotive energy which represents the quiver energy of a free particle averaged over one cycle. Here, the factor 2 is derived from the two kinds of particles (electrons and holes). In addition, considering the correspondence of AC Stark shift and ponderomotive shift, we can introduce the relation w = (1/4)ℏΩ r6 [30], and then the cutoff energy can be derived as s ≈ B + 1.6ℏΩ r6 . This equation reveals that the cutoff energy shifts toward higher-energy side, in short, the plateau region becomes broader, as the bang-gap energy becomes larger. Assuming an experimental situation where an AC electric field of 30 THz is applied to GaSe and the perturbative treatments are becoming inadequate, we can expect the Rabi-frequency and band-gap energy to be Ω I6 ≈ 0.5ω 6 and B ≈ 16ℏ 6 . Accordingly, we can conclude that the cutoff energy is roughly determined by the band-gap energy, which is consistent with recent experimental results [5].
FIG . 2(c) shows the higher-order harmonic spectra in semimetal regime, where the band-gap energies are B = 10ℏ 6 (blue line) and B = 5ℏ 6 (red line) in the case of Ω I6 = 8 6 . Different from multiphoton absorption and AC Zener regimes, the higher-order harmonics show a characteristic of white spectrum which masks the odd-order harmonics, and moreover, the cutoff energy no longer depends on the band-gap energy. The generation of the white higher-order harmonics means the emission law in each half cycle is broken. This half-cycle emission law can also be interpreted as temporal interference of HHG waves emitted dependently in a half cycle of the incident electric field [32,33]. On the basis of this idea, the breaking of the emission law is due to the overlap between the conduction and valence bands by the semimetallization for a long interval, which disturbs the interference that generates only odd-order harmonic spectra. To explain the independences of the cutoff energy from the band-gap energy, we extend the simple man model, as shown in FIG. 3(b). In this model, carriers are generated when the band gap becomes zero, after which they are accelerated and finally recombine. From such processes, we can derive the cutoff energy as s = 2× ℏΩ r6 + 3.2 w ≈ 3.6ℏΩ I6 , and this equation reveals that the cutoff energy no longer depends on the bandgap energy. These considerations are based on simple two-band models and their justification is assured by a condition that the band-gap energy is much smaller than energy differences between the first and the second conduction bands. Therefore, our results would be identified in experiments with several materials such as GaAs and InN.
FIGs. 4 (a) and (b) show the dependence of the cutoff frequency s on the Rabi frequency, as estimated from the numerical results (purple solid circles) and the simple man model (blue and red lines) in the cases of (a) B = 5ℏ 6 and (b) B = 10ℏ 6 . In the numerical calculation, the cutoff frequency s are determined as the threshold energy between the plateau and decay regions, where the plateau region is defined from a condition that the harmonic intensities become less than 10% compared with the previous and later orders. We found a crossover between AC Zener and semimetal regimes from the dependence of the cutoff energy on the Rabi frequency.
Finally, we will clarify dephasing effect on HHG spectra whose importance is discussed in the previous works [13,14]. The increasing dephasing effects enable us to presume disappearance of temporal interference of HHG waves, each of which is emitted dependently in a half cycle. Therefore, when employing the dephasing value of c ≈ 6 , we expect HHG spectra to be characterized by the single-cycle excitation processes, and consequently, the spectra widths become broader [34]. We identified this conjecture by performing our numerical calculations. This characteristic becomes important when considering HHG in semimetal regime. In this regime, the increasing dephasing effect changes noisy contimuumlike spectra into clean odd-ordered ones. We emphasize the cutoff laws derived in this paper are not influenced by the dephasing effect as long as considering a range of 0 ≤ c ≤ 6 , which covers recent experimental situations [4][5][6][7][8].
In this paper, based on an exact treatment of diabatic processes of Bloch electrons, we firstly discovered non-perturbative mechanisms of HHG in solid-state materials. In our theoretical framework, employing an external electric field in the form of a vector potential rather than the scalar one, justification of Bloch's theorem is ensured even under an influence of the high-intensity electric field. In this scheme, temporally changed band structures linked to an AC electric field characterize diabatic processes of Bloch electrons. Utilizing this framework, we revealed these diabatic processes, namely, AC Zener tunneling and semimetallization of semiconductors, determine properties of HHG spectra whose cutoff laws are in the form of s = B + 1.6ℏΩ r6 and s = 3.6ℏΩ I6 , respectively. The cutoff laws could be understood based on analogies between Zener tunneling and tunnel ionization processes, as well as semimetallization of semiconductors and over-the-barrier ionization processes [30,35], where band-gap energy and Bloch electrons correspond to the ionization energy and bare electrons, respectively. The analogies presented in this paper would provide possibilities for universal understandings of HHG mechanisms in both atomic and solid cases, and consequently, propose availability of the same method for coherent control such as a polarization gating [2,36].

Derivation of the Hamiltonian
To derive the Hamiltonian employed in our theory, we start from the well-known formula = (1/2 6 )( − / ) 8 + Σ ; − ; , where 6 is the electron mass, the electron charge, the velocity of light, the vector potential of incident electric fields, the momentum of the bare electrons, and − ; the periodic core potential of atoms located at ; . Ignoring the quasi- is the field operator of electrons. To consider the most basic structure of solid-state materials, we assume a two-dimensional covalent crystal in which two atoms A and B exist in a unit cell keeping the space-inversion symmetry. This assumption is equivalent to only focusing on the conduction and valence bands in a semiconductor. The incident light is assumed to be linearly polarized electric fields described as = .
Using the tight-binding model only considering nearestneighbor hopping of electrons and following the transformation In the following formulation, we will ignore the dependence of the Rabi frequency, which is usually permitted in semiconductor physics [2]. In this Hamiltonian, the diagonal matrix elements in Eqs. (5) and (6) become zero because of the space-inversion symmetry in a semiconductor and no carrier doping. The transformation to the band-structure picture can be performed by diagonalization of the single-particle part 6 through the use of the electron-hole picture defined as = 1/ 2 + ;‹ ( ) and ℎ F C = 1/ 2 − + ;‹ ( ) . The Hamiltonian in the electron-hole picture are derived as = 6 + ? , where 6 = Σ C + ℎ F C ℎ F and ? = ℏΩ I Σ cos Ž( ) the annihilation operator of electrons (holes). In this expression, we have ignored the constant energy term that gives the total-energy shift. To obtain a more specific description, we should include microscopic information on the arrangement of the atoms in the crystal. To consider the most simplified case, we assume ≈ (ℏ 8 8 /2 R + B /2) and Ž( ) = . Thus, we can derive a Hamiltonian in the form of The similar derivation of the Hamiltonian in the electron-hole picture was performed in graphene systems, where the hexagonal lattice structure is the essential property [3,4].

Some characteristics of the Hamiltonian
We should take care of the differences between the derived Hamiltonian and the conventional one employed in semiconductor physics. The conventional Hamiltonian in semiconductor physics usually takes the form [5], Here, we can see that the main difference of Eqs. (4) and (5) from Eqs. (6) and (7) is the first term in Eq.
(2) that is scaled by the Rabi frequency Ω I . This term indicates the intraband transition of Bloch electrons where the Bloch wavevector is a kinetic constant. In semiconductor physics, this process is usually omitted because of the selection rule of the dipole transition, where the band structure is assumed to be fixed even under a strong AC electric field (the rigid-band assumption) [5]. In contrast, our theory naturally includes this term which is renormalized in the single-particle energy R as R = R + B /2 + ℏΩ I cos , where the band structure varies with the strength of the AC electric field. Comparing these expressions and focusing on the characteristics of the first term in Eq.
(2) that is scaled by the Rabi frequency Ω I , we can conclude that the rigid-band assumption is only true for a weak electric field, in short, ℏΩ I ≪ B . Therefore, in treating the carrier dynamics under high-intensity electric fields, the modification of the band structure becomes important and leads to a violation of the selection rule, as has been pointed out in the discussions of gaseous media [6]. By comparing with the Landau-Zener model [7], we can also interpret the first term in Eq. (2) describes the diabatic processes such as the above-threshold ionization and Zener tunneling [1,8]. Other differences between the derived and conventional Hamiltonians are the factors sin and cos which are derived from the form factor defined by = Σ ; ; • Š . These factors reflect microscopic information on the arrangement of atoms in the crystal, which are never included in the conventional Hamiltonian employed in semiconductor physics. In particular, the factor cos , which characterizes the single-particle energy, R = R + B /2 + ℏΩ I cos , indicates an anisotropic deformation of the band structure, and the dipole transition under such anisotropic modifications leads to the aspect of transport phenomena.

Numerical results of kinetic equations
The above Hamiltonian enables us to derive time evolution equations for the populations R = C and polarization = ℎ F C with the Bloch wavevector in the form of ℏ R = 2Im Ω R sin † − g R .
Here, Ω I = Ω I6 exp − − 6 8 / 8 cos 6 is the Rabi frequency in which the intensity of the electric field is renormalized in Ω I6 , R = R + B /2 + ℏΩ I cos is the single-particle energy modulated by the AC electric field, and c and g are the transverse and longitudinal relaxation constants. Throughout this study, the parameters of the incident electric field will be fixed to t 6 = 12π/ω 6 and τ = 4π/ω 6 (see the top panel of FIG. 1), and the relaxation constants to c = 0.1ℏ 6 and g = 0.01ℏ 6 , respectively. As discussed in the main paper, the temporal variation of the single-particle energy indicates the anisotropic modification of the band structure, and the dipole transition under such modifications leads to the anisotropic distributions of the carrier densities and polarization in two-dimensional space. FIGs. 1(a1) and (b1) show the densities A = D and (a2) and (b2) show the distributions of the polarizations 2Im , where B = 5ℏ 6 and Ω I6 = 2 6 in the case of t = 12π/ω 6 and t = 13π/ω 6 . These figures mean that the anisotropic distributions in space are linked to the direction of the electric field and that they arise from a dipole transition under the anisotropic modification of the band structure characterized by the factor ℏΩ I cos . Therefore, we expect that the intraband processes lead to carrier transport (wave packet dynamics), and the time evolution of the current can be calculated using the definition = − ? / = 1 − A − k D cos − 2Im sin .