High-harmonic generation (HHG) in the gas phase^[1–3] has long been studied as a coherent soft X-ray radiation source with an extremely short pulse width^[4–6]. Its ability to probe the molecular structure^[7,8] and ultrafast electron dynamics^[9–11] has also been proved. Recently, HHG from solid state materials has attracted much attention^[12,13], as it provides an opportunity to probe electron dynamics in solid states and a potential to produce attosecond pulses with an enhanced yield^[14–24].

Because of the complexity of solid materials, HHG shows a variety of new phenomena and corresponding mechanisms. Especially, even-order harmonics are generated by breaking the inversion symmetry of the system^{[16,17,25,26]}. From a single crystal of gallium selenide (GaSe), even-order harmonics are produced with a comparable yield as the odd orders, which are attributed to the interference between multiple interband transitions^[16]. It has also been proved that the transition dipole phase (TDP) is a necessity for modeling HHG, with which even-order harmonics can be produced by only considering a direct two-band transition^[25,27]. Recently, the Berry connection is found crucial to calculate a physical HHG spectrum that includes even-order harmonics^[28,29]. It is noted that all available descriptions of even-order harmonic generation rely on the length gauge, and multiple parameters increase the complexity to interpret experimental observations.

As an electromagnetic gauge to describe the light field under dipole approximation, the velocity gauge preserves the translational symmetry of the crystal and decouples electron dynamics at each crystal momentum (k) point^[30], which can be advantageous to study laser–solid interaction. In this work, we study the generation of even-order harmonics by calculating the high-harmonic spectrum from an inversion-asymmetric one-dimensional periodic potential. Under the picture of velocity gauge, we proved that the even-order harmonics only originate from the multiple-step transition, while the odd orders are mainly from the direct one. The lower valence band is found crucial to the generation of even-order harmonics, which differs from the inversion-symmetric case. We find the relative intensity of even versus odd harmonics increases over pump laser intensity, which results from the growing hole population in the lower valence band.

We model the inversion-asymmetric crystal structure as a one-dimensional chain constituted by three kinds of atoms equally spaced [Fig. 1(a)]. The Coulomb potential of each atom is treated as a localized potential well^[31,32]. The total periodic potential is given byV(x)=∑R∑i=A,B,CVi[1+tanh(x−xi,R)][1−tanh(x−xi,R)],(1)in which VA=0.6, VB=0.7, and VC=0.8 in atomic units (a.u.). Atomic centers are separated by 1/3 of the lattice constant a=xi,R−xi,R−1=12a.u. Note that three different atoms in each lattice are the simplest asymmetric case for equally spaced atoms. We calculate the Bloch functions as eigenfunctions of the field-free Hamiltonian, and the energies of the first five bands (εn,k, n=1,2,3,4,5) are shown in Fig. 1(b).

Under dipole approximation, the multiple-electron dynamics of the system is calculated using the semiconductor Bloch equation (SBE) in the velocity gauge with a phenomenological relaxation rate^[29]:∂∂tPn,n′(k,t)=−i(ϵ˜n,k−ϵ˜n′,k−i/T2)Pn,n′(k,t)+iA(t)[ρn(k,t)−ρn′(k,t)]pn,n′(k)−iA(t)∑m[Pm,n′(k,t)pn,m(k)−Pn,m(k,t)pm,n′(k)],(2)∂∂tρn(k,t)=−ρn(k,t)/T1−2A(t)Im[∑mPn,m(k,t)pm,n(k)],(3)where pn,n′(k)=−i〈ϕn,k|∇k|ϕn′,k〉 is the transition matrix element calculated from Bloch functions, Pn,n′(k,t) is the microscopic polarization, ρn(k,t) is the electron population, and ε˜n,k is the modified band energy, ε˜n,k=εn,k+A(t)pn,n(k).

In our calculation, the laser pulse has an electric field E(t)=ELcos4(t/τ)cos(2ωLt) with a wavelength of λL=2πc/ωL=4μm and pulse duration of τ=63fs. The relaxation times T1 and T2 are set to be 7 fs and 1.1 fs, respectively^[16]. For initial conditions, we set ρ1(k,0)=ρ2(k,0)=1, while other ρn(k,0) and Pn,n′(k,0) are all zero [see the Fermi level shown in Fig. 1(b)]^[31,33,34].

After solving Eqs. (2) and (3), we calculate the electric current byj(t)∝∑k,n,n′Pn,n′*(k,t)pn,n′(k)+∑k,nρn(k,t)[pn,n(k)+A(t)].(4)

The high-harmonic spectrum is given byS(ω)∝ω2|j(ω)|2.(5)

The high-harmonic spectrum at the pump laser field amplitude of EL=2.2V/nm is shown in Fig. 2(a). Five bands in total have been involved in the calculation. The convergence of results in this spectrum range has been checked for up to 10 bands. The spectrum contains both odd and even orders. While being generally weaker than adjacent odd orders, a plateau of even-order harmonics, ranging from the 6th to 14th order, is clearly shown in Fig. 2(a). The abrupt increase of harmonic intensity between the 4th and 6th order coincides with the first Van Hove singularity point shown by the joint density of states (JDoS) [Fig. 2(b)]. For orders higher than the 14th, the intensity of even-order harmonics begins to decrease with increasing order. Such pattern complies with the interband mechanism^[20,35].

To determine the transition path that leads to the generation of even-order harmonics, we turn off the sequential transition by setting pn,n′=0 for all |n−n′|>1 in the calculation. As a result, one can see the even-order harmonics disappear completely in the spectrum (Fig. 3), which indicates that they originate from the multiple-band mechanism. The spectrum of odd harmonics also slightly changes in Fig. 3, which implies minor contributions from sequential transitions as well. The calculation is done under several laser field strengths (0.8–2.4 V/nm), and the result stays the same.

It is noted that a direct two-band transition can produce even-order harmonics in the length gauge, in which it is crucial that the TDP be included^[25,27]. In our calculation, since we have obtained the Bloch functions, the phase of the transition dipole is included. But, the spectrum produced by forbidding all sequential transitions in Fig. 3 does not include even orders. This is because in the velocity gauge, the two-band transition with the transition matrix element phase does not contribute to the even-order harmonics.

To prove the above statement analytically, we reduce Eq. (2) into a two-band case:$$\begin{gathered} ∂∂tPc,v(k,t)=−i(ε˜c,k−ε˜v,k−i/T2)Pc,v(k,t) \\ +iA(t)[1−ρ(k,t)]pc,v(k). \end{gathered}$$(6)

Note that in velocity-gauge SBE [Eqs. (2), (3), and (6)] the system evolution at each k is completely decoupled. Following Keldysh’s approach^[3,36,37], Eq. (6) can be written asPc,v(k,t)=i∫−∞t[1−ρ(k,t)]A(τ)pc,v(k)dτeiSe−t−τT2,(7)where S is the dynamical phase, S=−∫τtε˜n,kdt′.

Applying Eq. (4), the corresponding current isj(t)∝i∑kpc,v*(k)∫−∞t[1−ρ(k,t)]A(τ)pc,v(k)dτeiSe−t−τT2+c.c.(8)

In Eq. (8), the complex phase of pc,v inside and outside the integral sign cancels out and does not affect the result current. Following a similar routine, one can prove that the current generated by band population ρn(k,t) is not affected by the phase of pc,v either. Since, in such one-dimensional two-band systems, the asymmetry of the material is encoded in the phase of pc,v (which is linked to TDP by a phase shift of π/2)^[25], the irrelevance of the phase of pc,v means the asymmetry of the material is omitted in the high-harmonic current, and thus no even-order harmonics can be generated.

We want to clarify that our conclusion does not conflict that of Ref. [25], which finds that the TDP dominates the even-order harmonic generation process in ZnO, because the electromagnetic gauges we use are different. Under electromagnetic-gauge transformation, the corresponding electron transition path is known to be different^[38].

As the generation of even-order harmonics results from sequential transitions, we further specify the role of the individual band. In the inversion-symmetric system with a similar band structure^[33], it has been found that the lower valence band plays a negligible role in HHG. We test the applicability of the finding in our system by forbidding all transitions involving band 1 and compare the spectrum with the real one (Fig. 4). The result shows that when band 1 is excluded it does not affect much of the odd harmonics, which agrees with the previous finding. Meanwhile, nearly most of the even-order harmonics disappeared, which distinctively differs from the symmetric case. We have checked the result under several laser field strengths (0.8–2.4 V/nm). The key role of band 1 in even-order harmonic generation can be understood, as it includes the most bounded states of the system, and thus is sensitive to the asymmetry of the periodic potential.

We then calculate the laser field dependence of the yield of even-order harmonics on the odd orders. The result is shown in Fig. 5(a). At the field strength of EL=1.0V/nm, the even-order harmonics are weaker than adjacent odd harmonics by two to three orders of magnitude. As the laser field increases to EL=2.4V/nm, the yield of the 12th and 14th becomes of the same order of magnitude as adjacent odd orders. To show the dependence more clearly, we calculate the ratio between the even-order harmonics and their neighboring two odd orders, I2s/(I2s+1I2s−1)0.5. As shown in Fig. 5(b), the calculated ratio has an increasing trend and minor oscillation.

Our analysis above shows that the even-order harmonics are generated by sequential transitions involving band 1, while the odd orders are generated mostly from direct transitions. Next, we will show how the two processes depend differently on the increase of the laser field. A schematic diagram is shown in Fig. 6. For the odd harmonics, the harmonic intensity is determined by the probability of direct transition m→n, Om→n, while for the even-order harmonics, by the probability of sequential transition 1→m→n, O1→m→n. Note that the only difference of the two transitions is the involvement of band 1. Under low laser amplitude, band 1 is hardly involved with being tightly bounded. Under higher laser amplitude, the involvement of band 1 is increased, which gives rise to the ratio between the probabilities of the two transitions O1→m→n/Om→n. This explains the increasing trend of relative yield of even-order harmonics. Since the yield of high harmonics also depends on dynamical phase^[37], which changes as well with pump intensity, the relative intensity of even-order harmonics as a function of laser field is not monotonous, as is shown in Fig. 5(b). But, as we average the ratio of different even orders, the increasing trend is clear.

We check the above explanation by comparing the hole population in band 1, 1−ρ1, with averaged even/odd ratios. As shown in Fig. 7, the two curves roughly coincide with one another, which supports the validity of our explanation.

As the phenomenological relaxation rate in Eqs. (2) and (3) can suppress some electron trajectories for HHG^[39], the numerical results are quantitatively dependent on the relaxation rate used in our calculation. But, as even and odd harmonics mainly come from different electron transition paths, one can imagine that their yield should depend differently on laser parameters. Though our analysis is based on velocity gauge, considering gauge freedom, analysis including the TDP and Berry connection in length gauge^[25,28,29] should lead to the same result. In real systems, the second harmonics generated from the material front surface can also drive the electron and lead to the generation of even-order harmonics^[26]. This might be the reason that in some experimental results even and odd harmonics have similar laser intensity dependence^[16,40]. It has been shown that with reflection geometry such nonlinear propagation effects can be suppressed^[41].

In conclusion, we calculated the high-harmonics spectra from an inversion-asymmetric one-dimensional periodic potential. The result is understood in velocity gauge, where the even-order harmonics contribute solely to sequential transitions, and the odd harmonics mostly to direct transitions. The vital role of the lower valence band in even-order harmonic generation is found in our result. As the involvement of the lower valence band requires high laser intensity, the yield of even-order harmonics increases more rapidly with the laser field than odd-order harmonics. The difference between the laser intensity dependence of even and adjacent odd harmonics should be independent of gauge and can be examined by experiments free from propagation effects. Our work offers an alternative perspective to study the interplay between the intense laser field and inversion-asymmetric solids. As even-order harmonics are linked only to sequential transition processes, our work implies a way to resolve the electron motion using high-harmonic spectroscopy.