text:
# Theory of coherent acoustic phonons in InGaN/GaN multi-quantum wells

## I Introduction

In recent years, experiments have shown that optical excitation of electron-hole pairs in semiconductors by ultrafast lasers can coherently excite longitudinal optical phonon modes in semiconductors.  In uniform bulk semiconductors, since the laser wavelength is much larger than the lattice spacing, the photo-generated carriers are typically excited by the optical pump over spatial areas that are much larger than the lattice unit cell. As a result, the excited carrier populations are generated in a macroscopic state and the carrier density matrix has only a $`q0`$ Fourier component. Coupling of the photoexcited carriers to the phonons leads only to coherent optical phonon modes with $`q0`$. Since the frequency of the $`q0`$ acoustic phonon is zero, coherent acoustic phonons are not excited in bulk semiconductors.

In semiconductor superlattices, even though the laser pump has a wavelength large compared to the lattice spacing, the pump can preferentially generate electron-hole pairs in the wells. The result is to create photoexcited carrier distributions that have the periodicity of the superlattice. Since the density matrix of the photo-excited carrier populations now has a $`q0`$ Fourier component, the photo-excited carriers can not only couple to the optical phonon modes, but they can also generate coherent acoustic phonon modes with a nonzero frequency and wavevector $`q2\pi /L`$ where $`L`$ is the superlattice period. In superlattices, the coherent phonon oscillation of zone folded acoustic phonons has been observed in $`AlAs/GaAs`$ superlattices. . However, the reflection modulation, observed to be on the order of $`\mathrm{\Delta }R/R10^510^6`$, is very small.

Recently, C.-K. Sun et. al.  reported studies of coherent acoustic phonon oscillations in wurtzite (0001) $`InGaN/GaN`$ multi-quantum well samples with strain induced piezoelectric fields. Owing to the strong piezoelectric fields at the interfaces, huge coherent acoustic phonon oscillations were observed. The oscillations were strong enough to be seen in the transmision (rather than the usual reflectivity) with $`\mathrm{\Delta }T/T10^210^3`$. The oscillation frequency, in the THz range, corresponding to the LA phonon frequency with $`q2\pi /L`$, varied between samples in accordance with their different superlattice periods, $`L`$.

In this paper, we formulate a microscopic model for the generation of coherent acoustic phonons in strained wurtzite superlattices via ultrafast laser photo-excitation of real carriers. Whereas in bulk systems the microscopic theory of coherent LO phonons can be mapped onto a forced oscillator model , we show that coherent LA phonon generation in superlattices, under appropriate conditions, can be mapped onto a loaded string model which is readily solved for the lattice displacement. Since acoustic phonons are almost the same in the well as in the barrier, to lowest order we can treat the string as being uniform . The forcing term on the string, however, is not uniform since photoexcitation of carriers occurs only in the wells.

Our paper thus provides justification for using a simple, uniform string model with a nonuniform forcing term, rather than a more complicated microscopic theory. In addition, we provide a microscopic expression for the forcing term to use in the simplified string model. The string model provides additional insight into the physics of the coherent LA phonons.

## II Microscopic theory

In this section, we derive the microscopic theory for coherent acoustic phonon generation in superlattices and multiple quantum wells, including the effects of (i) bandstructure, (ii) strain, (iii) piezoelectric fields, (iv) Coulomb interactions, and (v) laser optical excitation. In section III, we will show how this reduces to a simplified driven uniform string model with a nonuniform forcing term and a microscopic expression for the forcing term.

We model photogeneration of electrons and holes and the subsequant excitation of coherent acoustic phonons in a multi-quantum well (MQW) $`pin`$ diode shown schematically in Fig. 1. The intrinsic active region consists of a left $`GaN`$ buffer region, several pseudomorphically strained (0001) $`In_xGa_{1x}N`$ quantum wells sandwiched between $`GaN`$ barriers, and a right $`GaN`$ buffer region as indicated in the figure. The $`P`$ and $`N`$ regions are assumed to be abruptly terminated $`p`$\- and $`n`$-doped $`GaN`$ bulk layers separated by a distance, $`L`$, across which a voltage drop, $`\mathrm{\Delta }V=V_A`$, is maintained. Photoexcitation of carriers is achieved by means of an ultrafast laser pulse incident normally along the (0001) growth direction, taken to coincide with the $`z`$-axis.

### A Bulk Bandstructure

In bulk systems, the conduction and valence bands in wurzite crystals including the effects of strain are treated using effective mass theory. Near the band edge, the effective mass Hamiltonian for electrons is described by a $`2\times 2`$ matrix which depends explicitly on electron wavevector, k, and the strain tensor, $`ϵ`$. The electron Bloch basis states are taken to be

$$|c,1=|S$$
(2)
$$|c,2=|S.$$
(3)
The conduction band Hamiltonian is diagonal and we have (relative to the bottom of the conduction band) (Refs.  and )

$`H_{2\times 2}^c(𝐤,ϵ)=\{{\displaystyle \frac{\mathrm{}^2k_z^2}{2m_z^{}}}+{\displaystyle \frac{\mathrm{}^2k_t^2}{2m_{xy}^{}}}`$ (4)
$`+a_{c,z}ϵ_{zz}+a_{c,xy}(ϵ_{xx}+ϵ_{yy})\}𝐈_{2\times 2}.`$ (5)
where $`𝐈_{2\times 2}`$ is the identity matrix. The electron effective masses along z (taken parallel to the c-axis) and in the $`xy`$ plane are $`m_z^{}`$ and $`m_{xy}^{}`$, respectively, $`k_t^2=k_x^2+k_y^2`$, and $`ϵ_{xx}`$, $`ϵ_{yy}`$ and $`ϵ_{zz}`$ are strain tensor components, and $`a_{c,z}`$ and $`a_{c,xy}`$ are the deformation potentials.

The Hamiltonian for the valence bands is a $`6\times 6`$ matrix. Following Ref. , the Hamiltonian (relative to the top of the valence band) can be block diagonalized into two degnerate $`3\times 3`$ submatrices if we adopt the Bloch basis states

$$|v,1=\frac{\alpha ^{}}{\sqrt{2}}|(X+iY)+\frac{\alpha }{\sqrt{2}}|(XiY)$$
(7)
$$|v,2=\frac{\beta }{\sqrt{2}}|(XiY)\frac{\beta ^{}}{\sqrt{2}}|(X+iY)$$
(8)
$$|v,3=\beta ^{}|Z+\beta |Z$$
(9)
$$|v,4=\frac{\alpha ^{}}{\sqrt{2}}|(X+iY)\frac{\alpha }{\sqrt{2}}|(XiY)$$
(10)
$$|v,5=\frac{\beta }{\sqrt{2}}|(XiY)+\frac{\beta ^{}}{\sqrt{2}}|(X+iY)$$
(11)
$$|v,6=\beta ^{}|Z+\beta |Z$$
(12)
The phase factors, $`\alpha `$ and $`\beta `$, are functions of the angle $`\varphi =\mathrm{tan}^1(k_y/k_x)`$ and are given by

$$\alpha (\varphi )=\frac{1}{\sqrt{2}}e^{i(3\pi /4+3\varphi /2)}$$
(14)
$$\beta (\varphi )=\frac{1}{\sqrt{2}}e^{i(\pi /4+\varphi /2)}.$$
(15)
The block diagonalized Hamiltonian can be written as

$$H_{6\times 6}^v(𝐤,ϵ)=\left(\begin{array}{cc}H_{3\times 3}^U(𝐤,ϵ)& 0\\ 0& H_{3\times 3}^L(𝐤,ϵ)\end{array}\right),$$
(16)
where the upper and lower blocks of the Hamiltonian are

$$H_{3\times 3}^U(𝐤,ϵ)=\left(\begin{array}{ccc}F& K_t& iH_t\\ K_t& G& \mathrm{\Delta }iH_t\\ iH_t& \mathrm{\Delta }+iH_t& \lambda \end{array}\right)$$
(18)
and

$$H_{3\times 3}^L(𝐤,ϵ)=\left(\begin{array}{ccc}F& K_t& iH_t\\ K_t& G& \mathrm{\Delta }+iH_t\\ iH_t& \mathrm{\Delta }iH_t& \lambda \end{array}\right).$$
(19)
The elements appearing in the $`3\times 3`$ Hamiltonian matrices are

$$F=\mathrm{\Delta }_1+\mathrm{\Delta }_2+\lambda +\theta $$
(21)
$$G=\mathrm{\Delta }_1\mathrm{\Delta }_2+\lambda +\theta $$
(22)
$$K_t=\frac{\mathrm{}^2}{2m_0}A_5k_t^2$$
(23)
$$H_t=\frac{\mathrm{}^2}{2m_0}A_6k_tk_z$$
(24)
$$\mathrm{\Delta }=\sqrt{2}\mathrm{\Delta }_3$$
(25)
$$\lambda =\frac{\mathrm{}^2}{2m_0}\left(A_1k_z^2+A_2k_t^2\right)+D_1ϵ_{zz}+D_2\left(ϵ_{xx}+ϵ_{yy}\right)$$
(26)
$$\theta =\frac{\mathrm{}^2}{2m_0}\left(A_3k_z^2+A_4k_t^2\right)+D_3ϵ_{zz}+D_4\left(ϵ_{xx}+ϵ_{yy}\right)$$
(27)
In Eq. (II A), the $`A_i^{}s`$ are effective mass parameters, the $`D_i^{}s`$ are the Bir-Pikus deformation potentials, and the $`\mathrm{\Delta }^{}s`$ are related to the crystal field splitting, $`\mathrm{\Delta }_{cr}`$, and spin-orbit splitting, $`\mathrm{\Delta }_{so}`$, by $`\mathrm{\Delta }_1=\mathrm{\Delta }_{cr}`$ and $`\mathrm{\Delta }_2=\mathrm{\Delta }_3=\mathrm{\Delta }_{so}/3`$. $`m_0`$ is the free electron mass.

### B Quantized carrier states in MQW diodes

In quantum confined systems such as the $`pin`$ diode shown in Fig. 1, we must modify the bulk Hamiltonian. The finite MQW structure breaks translational symmetry along the $`z`$ direction but not in the $`xy`$ plane. Thus, quantum confinement of carriers in the MQW active region gives rise to a set of two-dimensional subbands. The wavefunctions in the envelope function approximation are

$$\psi _{n,𝐤}^\alpha (𝐫)=\underset{j}{}\frac{e^{i𝐤\rho }}{\sqrt{A}}F_{n,k,j}^\alpha (z)|\alpha ,j,$$
(28)
where $`\alpha =\{c,v\}`$ refers to conduction or valence subbands, $`n`$ is the subband index, $`𝐤=(k_x,k_y,0)=(k,\varphi )`$ is the two-dimensional wavevector, and $`j`$ labels the spinor component. For conduction subbands, ($`\alpha `$=c) $`j=1,2`$ while for valence subbands ($`\alpha `$=v) $`j=1\mathrm{}6`$. The slowly varying envelope functions $`F_{n,k,j}^\alpha (z)`$ are real and depend only on $`k=|𝐤|`$, while the rapidly-varying Bloch basis states $`|\alpha ,j`$ defined in Eqs. (II A) and (II A) depend on $`\varphi `$ in the case of valence subbands as given in Eq. (12). The area of the MQW sample in the $`xy`$ plane is $`A`$, and $`\rho =(x,y,0)`$ is the projection of r in the plane.

The envelope functions satisfy a set of effective-mass Schrödinger equations

$$\underset{j,j^{}}{}\left\{H_{j,j^{}}^\alpha (k)+\delta _{j,j^{}}\left[V_\alpha (z)E_n^\alpha (k)\right]\right\}F_{n,k,j^{}}^\alpha (z)=0,$$
(29)
subject to the boundary conditions

$$F_{n,k,j}^\alpha (z=0)=F_{n,k,j}^\alpha (z=L)=0.$$
(30)
where $`L`$ is again the length of the MQW diode structure (c.f. Fig. 1), $`V_\alpha (z)`$ are the quantum confinement potentials for conduction and valence electrons. $`E_n^\alpha (k)`$ are the energy eigenvalues for the $`n`$th conduction or valence subband. Note that in the envelope function approximation, the subband energy depends only on the magnitude $`k`$ of the transverse wavevector and not on the angle $`\varphi `$. For the quantum confined case, the matrix operators $`H_{j,j^{}}^\alpha (k)`$ depend explicitly on $`z`$ and are obtained by making the replacement $`k_zi\frac{}{z}`$ and letting all material parameters be $`z`$-dependent operators in the matrices $`H^\alpha (k,ϵ)`$ given in Eqs. (5) and (16). To ensure the Hermitian property of the Hamiltonian, we make the operator replacements

$$B(z)\frac{^2}{z^2}\frac{}{z}B(z)\frac{}{z},$$
(32)
and

$$B(z)\frac{}{z}\frac{1}{2}\left[B(z)\frac{}{z}+\frac{}{z}B(z)\right].$$
(33)
The quantum confinement potentials $`V_c(z)`$ and $`V_v(z)`$ arise from (i) bandgap discontinuities between well and barrier regions, (ii) the strain-induced piezoelectric field, and (iii) the time-dependent electric field due to photoexcited electrons and holes. Thus,

$$V_\alpha (z,t)=V_{\alpha ,\text{gap}}(z)+V_{\text{piezo}}(z)+V_{\text{photo}}(z,t),$$
(34)
and the bandstructure is explicitly time-dependent.

Material parameters for $`InN`$ and $`GaN`$ used in this work can be found in Table I. For $`In_xGa_{1x}N`$ alloys, we interpolate between the $`GaN`$ and $`InN`$ values listed in the table.

We obtain electron effective masses by linearly interpolating the reciprocals of the masses as a function of the indium concentration $`x`$, i.e., the concentration dependent effective masses are taken to be

$$\frac{1}{m_{xy}^{}(x)}=x\left(\frac{1}{m_{xy}^{}}\right)_{\text{InN}}+(1x)\left(\frac{1}{m_{xy}^{}}\right)_{\text{GaN}},$$
(36)
and

$$\frac{1}{m_z^{}(x)}=x\left(\frac{1}{m_z^{}}\right)_{\text{InN}}+(1x)\left(\frac{1}{m_z^{}}\right)_{\text{GaN}}.$$
(37)
For the alloy band gap, $`E_g(x)`$, we use an expression incorporating a bowing paremeter:

$$E_g(x)=xE_{g,InN}+(1x)E_{g,GaN}bx(1x),$$
(38)
where the bowing parameter, $`b=1.0eV`$.

For all other material parameters, we use linear interpolation in $`x`$ to obtain values for the alloy. Since we can’t find deformation potentials for $`InN`$, we use $`GaN`$ values by default. In the absence of values of $`ϵ_{\mathrm{}}`$ for either $`GaN`$ or $`InN`$, we use linear interpolation in $`x`$ to obtain $`ϵ_0`$ and simply take $`ϵ_{\mathrm{}}ϵ_0`$.

If the z-dependent band gap in the MQW is $`E_g(z)`$ as determined from Eq. (38) and the Indium concentration profile, and $`E_{g,\text{min}}=\text{min}_z[E_g(z)]`$ is the minimum band gap in the structure, then the confinement potentials for conduction and valence electrons are defined as

$$V_{c,\text{gap}}(z)=E_{g,\text{min}}+Q_c(E_g(z)E_{g,\text{min}})$$
(40)
$$V_{v,\text{gap}}(z)=(1Q_c)(E_g(z)E_{g,\text{min}})$$
(41)
where the conduction band offset is taken as $`Q_c=0.6`$.  With these definitions for the $`V_{\alpha ,\text{gap}}(z)`$, the zero of the gap confinement potential is placed at the top of the valence band profile.

The confinement potentials due to the strain-induced piezoelectric field are given by

$$V_{\text{piezo}}(z)=|e|E_z^0(z),$$
(42)
where $`|e|`$ is the electric charge and $`E_z^0(z)`$ is the strain-induced piezoelectric field. In a pseudomorphically strained MQW diode, the bulk source and drain (assumed to have identical composition) are unstrained while the in-plane MQW lattice constants adjust to the source and drain values. For a MQW grown along  (the z-direction) the z-dependent strain is

$$ϵ_{xx}(z)=ϵ_{yy}(z)=\frac{a_0a(z)}{a(z)}.$$
(43)
Here $`a_0`$ is the lattice constant in the source and drain and $`a(z)`$ is the z-dependent lattice constant in the MQW. Minimizing the overall strain energy, we find

$$ϵ_{zz}(z)=\frac{2C_{13}(z)}{C_{33}(z)}ϵ_{xx}(z),$$
(44)
where $`C_{13}(z)`$ and $`C_{33}(z)`$ are z-dependent elastic constants.

There are several issues concerning strain that one can worry about. One is the critical well thickness beyond which the strain relaxes. Studies have shown  that with 10 % Indium, and 6 nm well width, a 350 kV/cm field is measured, implying that the well is not fully relaxed and justifies using the pseudomorphic strain approximation as we do in this paper for wells having 6 % Indium and thickness near 4 nm. For thick wells, this pseudormorphic strain approximation clearly will begin to break down and a more detailed model will be needed. Interface roughness can also play a role. The roughness will not significantly affect the acoustic phonon modes (see introduction) but may affect the photogeneration of carriers . For simplicity, we do not consider interface roughness.

The strain-induced polarization directed along the z-direction is given by

$$P_z^0(z)=e_{31}(z)\left(ϵ_{xx}(z)+ϵ_{yy}(z)\right)+e_{33}(z)ϵ_{zz}(z),$$
(45)
where $`e_{31}(z)`$ and $`e_{33}(z)`$ are z-dependent piezoelectric constants. The unscreened piezoelectric field in the diode is obtained from the requirement that the electric displacement vanishes . Thus,

$$E_z^0(z)=\frac{4\pi }{\epsilon _0(z)}\left(P_z^0(z)+P_0\right),$$
(46)
where $`P_0`$ is a constant polarization induced by externally applied voltages and $`\epsilon _0(z)`$ is the position-dependent static dielectric constant. The value of $`P_0`$ is obtained from the voltage drop across the diode (of length $`L`$) in the unscreened limit, i.e. with no photoexcited carriers. In this limit, the voltage drop between source and drain due to the induced piezoelectric field is just

$$V_A=_0^L𝑑zE_z^0(z),$$
(47)
from which $`P_0`$ can be determined.

When photoexcited electrons and holes are generated by the laser, then there is an additional time-dependent confinement potential,

$$V_{\text{photo}}(z,t)=|e|E_z^{\text{photo}}(z,t).$$
(48)
This potential is obtained by solving the Poisson equation in the diode for $`E_z^{\text{photo}}(z,t)`$ subject to the boundary condition

$$V_A=_0^L𝑑zE_z^{\text{total}}(z,t).$$
(49)
Here, $`E_z^{\text{total}}(z,t)`$ is the total electric field and is just the sum of the strain induced electric field and the field due to photogenerated electrons and holes, i.e.,

$$E_z^{\text{total}}(z,t)=E_z^0(z)+E_z^{\text{photo}}(z,t).$$
(50)
Finally, we can write an effective-mass Schrödinger equation for the conduction electron envelope functions in terms of an effedctive electron potential, $`V_c^{eff}(z)`$:

$`{\displaystyle \frac{\mathrm{}^2}{2}}\left\{{\displaystyle \frac{}{z}}{\displaystyle \frac{1}{m_z^{}(z)}}{\displaystyle \frac{}{z}}\right\}F_{n,k,j}^c(z)`$ (51)
$`+\left\{V_c^{eff}(z)E_n^c(k)\right\}F_{n,k,j}^c(z)=0,`$ (52)
where the effective electron potential is

$`V_c^{eff}(z)=V_c(z)+{\displaystyle \frac{\mathrm{}^2k^2}{2m_{xy}^{}(z)}}`$ (53)
$`+a_{c,z}(z)ϵ_{zz}(z)+a_{c,xy}(z)\left(ϵ_{xx}(z)+ϵ_{yy}(z)\right).`$ (54)
Similar expressions can be derived for valence electrons. We arrive at a set of coupled ordinary differential equations (ODE’s) subject to the two-point boundary value condition of Eq. (30). These are solved for the envelope functions and subband energies.

In practice, we introduce a uniform grid, $`\{z_i\}`$, along the z-direction and finite-difference the effective mass Schrödinger equations to obtain a matrix eigenvalue problem which can be solved using standard matrix eigenvalue routines. The resulting eigenvalues are the subband energies, $`E_n^\alpha (k)`$, and the corresponding eigenvectors are the envelope functions, $`F_{n,k,j}^\alpha (z_i)`$, defined on the finite difference mesh.

### C Second quantized electron hamiltonians

We next describe the second quantized Hamiltonians for electrons moving freely in the MQW interacting via a screened Coulomb potential. We denote creation and destruction operators for electrons in conduction and valence subbands by $`c_{\alpha ,n,𝐤}^{}`$ and $`c_{\alpha ,n,𝐤}`$, respectively. The second quantized Hamiltonian for free electrons and holes is simply

$$_{e0}=\underset{\alpha ,n,𝐤}{}E_n^\alpha (k)c_{\alpha ,n,𝐤}^{}c_{\alpha ,n,𝐤}.$$
(55)
The Coulomb interaction Hamiltonian is given by

$`_{ee}`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle \underset{\alpha ,n,𝐤}{}}{\displaystyle \underset{\alpha ^{},n^{},𝐤^{}}{}}{\displaystyle \underset{\kappa 0}{}}𝒱_{\alpha ^{},n^{},𝐤^{}}^{\alpha ,n,𝐤}(\kappa )`$ (57)
$`\times c_{\alpha ,n,𝐤\kappa }^{}c_{\alpha ^{},n^{},𝐤^{}+\kappa }^{}c_{\alpha ^{},n^{},𝐤^{}}c_{\alpha ,n,𝐤}.`$
Eq. (57) describes two-body interactions where electrons in states $`|\alpha ,n,𝐤`$ and $`|\alpha ^{},n^{},𝐤^{}`$ scatter to subband states $`|\alpha ,n,𝐤\kappa `$ and $`|\alpha ^{},n^{},𝐤^{}+\kappa `$, respectively. Note that to simplify things, we have neglected terms corresponding to Coulomb induced interband transitions (the ”diagonal approximation”) since these are very unfavorable energetically.  The electrons thus stay in their original subbands (though they may scatter off of other electrons in different subbands) and exchange crystal momentum $`\kappa `$. The matrix elements describing the strength of these transitions are given by

$`𝒱_{\alpha ^{},n^{},𝐤^{}}^{\alpha ,n,𝐤}(\kappa )`$ $`=`$ $`{\displaystyle 𝑑z𝑑z^{}V_{|\kappa |}(zz^{})}`$ (60)
$`\times {\displaystyle \underset{j}{}}F_{n,|𝐤\kappa |,j}^\alpha (z)F_{n,k,j}^\alpha (z)`$
$`\times {\displaystyle \underset{j^{}}{}}F_{n^{},|𝐤^{}+\kappa |,j}^\alpha ^{}(z^{})F_{n^{},k^{},j^{}}^\alpha ^{}(z^{}).`$
From Eq. (60), it is apparent that the symmetry relation

$$𝒱_{\alpha ^{},n^{},𝐤^{}}^{\alpha ,n,𝐤}(\kappa )=𝒱_{\alpha ,n,𝐤}^{\alpha ^{},n^{},𝐤^{}}(\kappa )$$
(61)
must hold.

The Fourier transform in the $`xy`$ plane of the screened Coulomb potential depends only on $`\kappa |\kappa |`$ and $`|z|`$ and is given by

$$V_\kappa (z)=\frac{2\pi e^2}{ϵ_0A}\frac{e^{\kappa |z|}}{\kappa ϵ_s(\kappa )}.$$
(62)
To describe screening, we adopt an effective pseudodynamic dielectric function of the form

$$\frac{1}{ϵ_s(\kappa )}=\frac{\kappa }{\kappa +\kappa _s}.$$
(63)
In the pseudodynamic screening model, we completely neglect screening by the massive holes and treat screening by the lighter conduction electrons in the static screening limit. The screening wavevector, $`\kappa _s`$, is computed in the $`2D`$ limit. Thus

$$\kappa _s=\frac{2\pi e^2}{ϵ_0}\frac{N_{2D}}{\mu },$$
(64)
where $`N_{2D}`$, the two-dimensional conduction electron density, is related to an effective chemical potential, $`\mu `$, by

$$N_{2D}=\frac{m_{xy}^{}k_BT}{2\pi \mathrm{}^2}\underset{n}{}\mathrm{ln}\left[1+\mathrm{exp}\left(\frac{E_n^c(0)\mu }{k_BT}\right)\right].$$
(65)
In Eq. (65), $`E_n^c(0)`$ is the conduction subband energy evaluated at $`k=0`$. In our simulation, the value of $`\mu `$ is obtained by requiring that $`N_{2D}`$, evaluated using Eq. (65), be equal to the value

$$N_{2D}(t)=\underset{n,𝐤}{}f_n^c(k,t)$$
(66)
obtained from the time dependent conduction electron distribution functions, $`f_n^c(k,t)`$.

### D Photogeneration of carriers

Electron-hole pairs are created by the pump laser and we treat the electric field of the laser in the semiclassical dipole approximation. In this approximation, the electron-laser interaction Hamiltonian is

$$_{eL}=|e|𝐄(𝐭)\underset{n,n^{},𝐤}{}\left[𝐝_{n,n^{}}^{c,v}(𝐤)c_{c,n,𝐤}^{}c_{v,n^{},𝐤}+\text{h.c.}\right],$$
(67)
where h.c. denotes the Hermitian conjugate of the first term. The laser field is $`𝐄(𝐭)`$ and the dipole matrix elements are

$$𝐝_{n,n^{}}^{c,v}(𝐤)=\underset{j,j^{}}{}𝐃_{j,j^{}}^{c,v}(\varphi )𝑑zF_{n,k,j}^c(z)F_{n^{},k,j^{}}^v(z).$$
(68)
The vector operator, $`𝐃_{j,j^{}}^{c,v}(\varphi )`$, is a $`2\times 6`$ matrix with $`x`$, $`y`$, and $`z`$ components. Thus,

$$𝐃_{j,j^{}}^{c,v}(\varphi )D_X^{c,v}(\varphi )\widehat{𝐱}+D_Y^{c,v}(\varphi )\widehat{𝐲}+D_Z^{c,v}(\varphi )\widehat{𝐳},$$
(69)
where $`\widehat{𝐱}`$, $`\widehat{𝐲}`$, and $`\widehat{𝐳}`$ are unit vectors and

$$D_X^{c,v}(\varphi )=\frac{P_2}{\sqrt{2}E_g}\left[\begin{array}{cccccc}\alpha & \beta ^{}& 0& \alpha & \beta ^{}& 0\\ \alpha ^{}& \beta & 0& \alpha ^{}& \beta & 0\end{array}\right].$$
(71)
$$D_Y^{c,v}(\varphi )=\frac{iP_2}{\sqrt{2}E_g}\left[\begin{array}{cccccc}\alpha & \beta ^{}& 0& \alpha & \beta ^{}& 0\\ \alpha ^{}& \beta & 0& \alpha ^{}& \beta & 0\end{array}\right],$$
(72)
$$D_Z^{c,v}(\varphi )=\frac{P_1}{E_g}\left[\begin{array}{cccccc}0& 0& \beta & 0& 0& \beta \\ 0& 0& \beta ^{}& 0& 0& \beta ^{}\end{array}\right].$$
(73)
The $`6\times 2`$ vector operator, $`𝐃_{j^{},j}^{v,c}(\varphi )`$, is related to the $`2\times 6`$ operator, $`𝐃_{j,j^{}}^{c,v}(\varphi )`$, by

$$𝐃_{j^{},j}^{v,c}(\varphi )=\left(𝐃_{j,j^{}}^{c,v}(\varphi )\right)^{}.$$
(74)
In Eq. (69), $`\alpha `$ and $`\beta `$ are the $`\varphi `$-dependent phase factors defined in Eq. (12) and the Kane parameters, $`P_1`$ and $`P_2`$, for wurtzite materials are related to the effective masses and energy gaps by

$`P_1^2`$ $`=`$ $`{\displaystyle \frac{\mathrm{}^2}{2m_0}}\left({\displaystyle \frac{m_0}{m_z^{}}}1\right)`$ (76)
$`\times `$ $`{\displaystyle \frac{(E_g+\mathrm{\Delta }_1+\mathrm{\Delta }_2)(E_g+2\mathrm{\Delta }_2)2\mathrm{\Delta }_3^2}{E_g+2\mathrm{\Delta }_2}},`$ (77)
$`P_2^2`$ $`=`$ $`{\displaystyle \frac{\mathrm{}^2}{2m_0}}\left({\displaystyle \frac{m_0}{m_{xy}^{}}}1\right)`$ (78)
$`\times `$ $`{\displaystyle \frac{E_g\{(E_g+\mathrm{\Delta }_1+\mathrm{\Delta }_2)(E_g+2\mathrm{\Delta }_2)2\mathrm{\Delta }_3^2\}}{(E_g+\mathrm{\Delta }_1+\mathrm{\Delta }_2)(E_g+\mathrm{\Delta }_2)\mathrm{\Delta }_3^2}}.`$ (79)
For the semiclassical laser field, we write the real electric field as

$$𝐄(t)=\frac{1}{2}\left[\widehat{ϵ}(t)e^{i\omega t}+\widehat{ϵ}^{}(t)e^{i\omega t}\right]$$
(80)
where $`\omega `$ is the photon frequency, $`\widehat{ϵ}`$ is a complex unit polarization vector, and $`(t)`$ is the pulse shape envelope function. We assume a Gaussian pulse shape

$$(t)=_0\mathrm{exp}\left[\left(\frac{tt_0}{\tau \sqrt{\frac{1}{2\mathrm{ln}2}}}\right)^2\right]$$
(81)
centered at $`t=t_0`$ with an intensity full width at half maximum (FWHM) of $`\tau `$. The maximum electric field strength, $`_0`$, is related to the pump fluence, $``$, by

$$_0=\sqrt{\frac{16\pi }{cn_\omega \tau }\sqrt{\frac{\mathrm{ln}2}{\pi }}},$$
(82)
where $`n_\omega `$ is the index of refraction at the photon frequency.

For linearly polarized light incident normally on the MQW, the polarization vectors are real and given by either $`\widehat{𝐱}`$ or $`\widehat{𝐲}`$. For circularly polarized light, the polarization vectors are complex and given by

$$\widehat{ϵ}_\pm =\frac{\widehat{𝐱}\pm i\widehat{𝐲}}{\sqrt{2}}.$$
(83)
In Eq. (83), the upper sign refers to left circularly polarized light (positive helicity) and the lower sign refers to right circularly polarized light (negative helicity).

### E Coupling to LA phonons

We treat the acoustic phonons in the MQW as bulk-like plane-wave states with wavevector $`𝐪`$. Since the system exhibits cylindrical symmetry, only $`𝐪=q\widehat{𝐳}`$ longitudinal acoustic phonons are coupled by the electron-phonon interaction. The free LA phonon Hamiltonian can be written as

$$_{A0}=\underset{q}{}\mathrm{}\omega _qb_q^{}b_q.$$
(84)
where $`b_q^{}`$ and $`b_q`$ are creation and destruction operators for LA phonons with wavevector $`𝐪=q\widehat{𝐳}`$. The wave vector component, $`q`$, of LA phonons in the MQW is thus defined in an extended zone scheme where $`\mathrm{}<q<\mathrm{}`$. The phonon dispersion relation is given by a linear relation

$$\omega _q=C_s|q|=\sqrt{\frac{C_{33}}{\rho _0}}|q|,$$
(85)
where $`\rho _0`$ is the mass density and $`C_s`$ is just the LA phonon sound speed for propagation parallel to $`\widehat{𝐳}`$.  In computing the LA sound speed in the linear phonon dispersion relation of Eq. (85), we neglect the $`z`$-dependence of the material parmeters and use bulk $`GaN`$ values for $`C_{33}`$ and $`\rho _0`$.

The LA phonons in wurtzite MQW’s interact with the electrons through deformation potential and screened piezoelectric scattering. The electron-LA phonon interaction in an MQW is governed by the Hamiltonian

$$_{eA}=\underset{\alpha ,n,n^{},𝐤,q}{}_{n,n^{}}^\alpha (k,q)(b_q+b_q^{})c_{\alpha ,n,𝐤}^{}c_{\alpha ,n^{},𝐤}.$$
(86)
This Hamiltonian describes the scattering of an electron from subband state $`|\alpha ,n^{},𝐤`$ to subband state $`|\alpha ,n,𝐤`$ with either the emission or absorption of an LA phonon. We note that the electron wave vector, $`𝐤`$, in the $`xy`$ plane is conserved in this process since, as noted earlier, the phonon wave vector in the $`xy`$ plane is zero.

The interaction matrix elements describing deformation and screened piezoelectric scattering are

$`_{n,n^{}}^\alpha `$ $`(k,q)=\sqrt{{\displaystyle \frac{\mathrm{}^2}{2\rho _0(\mathrm{}\omega _q)V}}}`$ (88)
$`\times \left[iq𝒟_{n,n^{}}^\alpha (k,q){\displaystyle \frac{4\pi |e|e_{33}}{ϵ_{\mathrm{}}ϵ_s(q)}}𝒫_{n,n^{}}^\alpha (k,q)\right],`$
where $`V`$ is the crystal volume. The first term in Eq. (88) desribes deformation potential scattering while the second term describes screened piezoelectric scattering.

The relative strengths of the various transitions are determined by form factors for deformation potential and piezoelectric scattering. The form factor for screened piezoelectric scattering is given by

$$𝒫_{n,n^{}}^\alpha (k,q)=\underset{j}{}𝑑zF_{n,k,j}^\alpha (z)e^{iqz}F_{n^{},k,j}^\alpha (z),$$
(89)
while the form factor for deformation potential scattering is defined to be

$$𝒟_{n,n^{}}^\alpha (k,q)=\underset{j}{}\mathrm{\Theta }_j^\alpha 𝑑zF_{n,k,j}^\alpha (z)e^{iqz}F_{n^{},k,j}^\alpha (z).$$
(90)
The form factor for deformation potential scattering is similar to the form factor for piezoelectric scattering except that in summing over spinor components, $`j`$, the terms are weighted by $`j`$-dependent deformation potentials $`\mathrm{\Theta }_j^\alpha `$ which can be representd by the row vectors

$$\mathrm{\Theta }_j^c=\{a_{c,z},a_{c,z}\},$$
(92)
$$\mathrm{\Theta }_j^v=\{D_1+D_3,D_1+D_3,D_1,D_1+D_3,D_1+D_3,D_1\}.$$
(93)
for conduction and valence electrons, respectively.

### F Electron density matrices

We define statistical operators in terms of the electron and phonon eigenstates. The electron density matrix is

$$N_{n,n^{}}^{\alpha ,\alpha ^{}}(𝐤,t)c_{\alpha ,n,𝐤}^{}(t)c_{\alpha ^{},n^{},𝐤}(t),$$
(94)
where $``$ denotes the statistial average of the non-equlilbrium state of the system.

The interband components of the density matrix, $`N_{n,n^{}}^{c,v}(𝐤,t)`$ and $`N_{n^{},n}^{v,c}(𝐤,t)`$, describe the coherence between conduction and valence electrons in subbands $`n`$ and $`n^{}`$ and are related to the optical polarization. The intraband components of the density matrix, $`N_{n,n^{}}^{\alpha ,\alpha }(𝐤,t)`$ describe correlations between different subbands of the same carrier type if $`nn^{}`$. If $`n=n^{}`$, $`N_{n,n}^{\alpha ,\alpha }(𝐤,t)f_n^\alpha (𝐤,t)`$ is just the carrier distribution function for electrons in the subband state, $`\psi _{n,𝐤}^\alpha (𝐫)`$, defined in Eq. (28).

### G Coherent phonon amplitude

The coherent phonon amplitude of the $`q`$-th phonon mode, $`|q`$, is defined to be

$$D_q(t)b_q(t)+b_q^{}(t).$$
(95)
The coherent phonon amplitude is related to the macroscopic lattice displacement, $`U(z,t)`$, and velocity, $`V(z,t)`$, through the relations

$$U(z,t)=\underset{q}{}\sqrt{\frac{\mathrm{}^2}{2\rho _0(\mathrm{}\omega _q)V}}e^{iqz}D_q(t)$$
(96)
$$V(z,t)=\underset{q}{}\sqrt{\frac{\mathrm{}^2}{2\rho _0(\mathrm{}\omega _q)V}}e^{iqz}\frac{D_q(t)}{t}$$
(97)
The coherent phonon amplitude, $`D_q(t)`$, will vanish if there are a definite number of phonons in the mode, i.e., if the phonon oscillator is in one of its energy eigenstates, $`|q`$. In this case, there is no macroscopic displacement of the lattice.

The coherent phonon distribution is

$$𝒩_q^{coh}(t)b_q^{}(t)b_q(t)$$
(98)
and the total phonon distribution, $`𝒩_q(t)`$, can be separated into coherent and incoherent contributions as follows;

$`𝒩_q(t)`$ $`=b_q(t)b_q^{}(t)`$ (100)
$`𝒩_q^{coh}(t)+𝒩_q^{incoh}(t).`$
In general, a mode can have a number of both coherent and incoherent phonons, but only the coherent phonons contribute to the macroscopic lattice displacement.

We note that at the beginning of the experiment, there are no coherent phonons present, i.e. $`𝒩_q^{coh}(t)=0`$, and the incoherent phonon population is described by a thermal distribution, $`𝒩_q^{incoh}(t)e^{\mathrm{}\omega _q/k_BT}`$.

### H Equations of motion

In this section we develop equations of motion for the electron density matrices and coherent phonon amplitudes. The electron density matrices obey the general equations of motion

$$\frac{N_{n,n^{}}^{\alpha ,\alpha ^{}}(𝐤,t)}{t}=\frac{i}{\mathrm{}}[,c_{\alpha ,n,𝐤}^{}c_{\alpha ^{},n^{},𝐤}],$$
(101)
where $`[]`$ denotes the commutator and $``$ denotes the average over an initial ensemble. The density matrices are defined in the electron picture and initially the valence bands are filled while the conduction bands are empty. We have $`f_n^c(𝐤,t=\mathrm{})=0`$ and $`f_n^v(𝐤,t=\mathrm{})=1`$ which implies

$$N_{n,n^{}}^{\alpha ,\alpha ^{}}(𝐤,t=\mathrm{})=\delta _{n,n^{}}\delta _{\alpha ,v}\delta _{\alpha ^{},v}.$$
(102)
The total Hamiltonian, $``$, is the sum of the Hamiltonians described in the previous sections, i.e.

$$=_{e0}+_{ee}+_{eL}+_{A0}+_{eA}.$$
(103)
In deriving equations of motion for the density matrices, we make the ansatz that the density matrices depend only on $`k=|𝐤|`$. We use the rotating wave approximation (RWA) to factor out the rapid $`e^{i\omega t}`$ behavior of the interband density matrix elements, $`N_{n,n^{}}^{c,v}(k,t)`$. In the RWA, we have

$$N_{n,n^{}}^{c,v}(k,t)\stackrel{~}{N}_{n,n^{}}^{c,v}(k,t)e^{i\omega t},$$
(104)
where $`\stackrel{~}{N}_{n,n^{}}^{c,v}(k,t)`$ is a slowly varying envelope function. In addition, we treat the Coulomb interaction in the time-dependent Hartree-Fock approximation by factoring four-operator averages arising from $`_{eL}`$ into appropriate products of two-operator averages as described in Ref. .

The resulting equations of motion for the density matrices are

$`{\displaystyle \frac{N_{n,n^{}}^{c,c}(k,t)}{t}}={\displaystyle \frac{i}{\mathrm{}}}\left\{_n^c(k)_n^{}^c(k)\right\}N_{n,n^{}}^{c,c}(k)`$ (106)
$`i{\displaystyle \underset{m}{}}\left\{\mathrm{\Omega }_{n,m}^{c,v}(k)\stackrel{~}{N}_{m,n^{}}^{v,c}(k)\stackrel{~}{N}_{n,m}^{c,v}(k)\mathrm{\Omega }_{n,n^{}}^{v,c}(k)\right\}`$ (107)
$`+{\displaystyle \frac{i}{\mathrm{}}}{\displaystyle \underset{m}{}^{}}\left\{\mathrm{\Lambda }_{n,m}^c(k)N_{m,n^{}}^{c,c}(k)N_{n,m}^{c,c}(k)\mathrm{\Lambda }_{m,n^{}}^c(k)\right\}`$ (108)
$`{\displaystyle \frac{N_{n,n^{}}^{v,v}(k,t)}{t}}={\displaystyle \frac{i}{\mathrm{}}}\left\{_n^v(k)_n^{}^v(k)\right\}N_{n,n^{}}^{v,v}(k)`$ (109)
$`i{\displaystyle \underset{m}{}}\left\{\mathrm{\Omega }_{n,m}^{v,c}(k)\stackrel{~}{N}_{m,n^{}}^{c,v}(k)\stackrel{~}{N}_{n,m}^{v,c}(k)\mathrm{\Omega }_{n,n^{}}^{c,v}(k)\right\}`$ (110)
$`+{\displaystyle \frac{i}{\mathrm{}}}{\displaystyle \underset{m}{}^{}}\left\{\mathrm{\Lambda }_{n,m}^v(k)N_{m,n^{}}^{v,v}(k)N_{n,m}^{v,v}(k)\mathrm{\Lambda }_{m,n^{}}^v(k)\right\}`$ (111)
$`{\displaystyle \frac{\stackrel{~}{N}_{n,n^{}}^{c,v}(k,t)}{t}}={\displaystyle \frac{i}{\mathrm{}}}\left\{_n^c(k)_n^{}^v(k)\mathrm{}\omega \right\}\stackrel{~}{N}_{n,n^{}}^{c,v}(k)`$ (112)
$`i{\displaystyle \underset{m}{}}\left\{\mathrm{\Omega }_{n,m}^{c,v}(k)N_{m,n^{}}^{v,v}(k)N_{n,m}^{c,c}(k)\mathrm{\Omega }_{m,n^{}}^{c,v}(k)\right\}`$ (113)
$`+{\displaystyle \frac{i}{\mathrm{}}}{\displaystyle \underset{m}{}^{}}\left\{\mathrm{\Lambda }_{n,m}^c(k)\stackrel{~}{N}_{m,n^{}}^{c,v}(k)\stackrel{~}{N}_{n,m}^{c,v}(k)\mathrm{\Lambda }_{m,n^{}}^v(k)\right\}.`$ (114)
The equations of motion for $`\stackrel{~}{N}_{n,n^{}}^{v,c}(k,t)`$ are redundant since $`\stackrel{~}{N}_{n,n^{}}^{v,c}(k,t)=(\stackrel{~}{N}_{n^{},n}^{c,v}(k,t))^{}`$.

The first terms on the right-hand side of Eq. (II H) describe the free oscillation of the density matrices in the renormalized singel particle energy bands. The time-dependent single particle energies are

$$_n^\alpha (k,t)=E_n^\alpha (k)+\mathrm{\Lambda }_{n,n}^\alpha (k,t),$$
(115)
where $`E_n^\alpha (k)`$ are the single particle subband energies in the absence of conduction electrons and holes and $`\mathrm{\Lambda }_{n,n}^\alpha (k,t)`$ describes the time-dependent renormalization of the single particle subbands.

The renormalization energies, $`\mathrm{\Lambda }_{n,n}^\alpha (k,t)`$, are the diagonal elements of a generalized renormalization energy matrix (in the subband indices)

$$\mathrm{\Lambda }_{n,n^{}}^\alpha (k,t)=\mathrm{\Sigma }_{n,n^{}}^\alpha (k,t)+Q_{n,n^{}}^\alpha (k,t).$$
(116)
The first term in the renormalization energy matrix (116) is the generalized exchange self-energy matrix arising from the Coulomb interaction and is given by

$`\mathrm{\Sigma }_{n,n^{}}^\alpha (k,t){\displaystyle \underset{𝐤^{}𝐤}{}}𝒱_{\alpha ,n^{},k^{}}^{\alpha ,n,k}\left(|𝐤𝐤^{}|\right)`$ (117)
$`\times \left(N_{n,n^{}}^{\alpha ,\alpha }(k^{},t)\delta _{\alpha ,v}\delta _{n,n^{}}\right),`$ (118)
where $`𝒱_{\alpha ,n^{},k^{}}^{\alpha ,n,k}(|𝐤𝐤^{}|)`$ are angular averaged Coulomb interaction matrix elements. The second term in Eq. (116) accounts for renormalization due to coupling of carriers to coherent acoustic phonons. We have

$$Q_{n,n^{}}^\alpha (k,t)\underset{q}{}D_q(t)_{n,n^{}}^\alpha (k,q),$$
(119)
where $`D_q(t)`$ is the coherent phonon amplitude and the electron-phonon matrix elements, $`_{n,n^{}}^\alpha (k,q)`$, are defined in Eq. (88). The self energy corrections in Eq. (117) are small, though the can be important in some circumstances.

In computing the angular averaged Coulomb matrix elements in Eq. (117), we assume small momentum transfer, $`\kappa `$, and use the fact that the envelope functions, $`F_{n,k,j}^\alpha (z)`$, depend weakly on $`k`$ to obtain an effective interaction,

$`𝒱_{\alpha ^{},n^{},k^{}}^{\alpha ,n,k}(\kappa ){\displaystyle 𝑑z𝑑z^{}V_\kappa (zz^{})}`$ (120)
$`\times {\displaystyle \underset{j,j^{}}{}}|F_{n,(k,k^{}),j}^\alpha (z)|^2|F_{n^{},(k,k^{}),j^{}}^\alpha ^{}(z^{})|^2`$ (121)
where, by definition,

$$|F_{n,(k,k^{}),j}^\alpha (z)|^2\frac{|F_{n,k,j}^\alpha (z)|^2+|F_{n,k^{},j}^\alpha (z)|^2}{2}.$$
(122)
The effective Coulomb interaction, $`𝒱`$, defined in Eq. (121) is an even function of $`\kappa `$ and is symmetric in $`\alpha `$ and $`k`$, thus preserving the symmetry relation (61). Preserving this symmetry is essential in order to maintain conservation of carriers in the scattering process.

The second terms in Eq. (II H) describe photoexcitation of electron-hole pairs by the pump laser. The system reacts to an effective field which is the sum of the applied field and the dipole field of the electron-hole excitations. This gives rise to a matrix of generalized Rabi frequencies in the subband indices

$`\mathrm{}\mathrm{\Omega }_{n,n^{}}^{c,v}(k)={\displaystyle \frac{(t)}{2}}d_{n,n^{}}^{c,v}(k)+`$ (123)
$`{\displaystyle \underset{𝐤^{}𝐤}{}}𝒱_{v,n^{},k^{}}^{c,n,k}\left(|𝐤𝐤^{}|\right)\stackrel{~}{N}_{n,n^{}}^{c,v}(k^{},t),`$ (124)
which can be shown to satisfy the symmetry relations

$$\mathrm{}\mathrm{\Omega }_{n,n^{}}^{c,v}(k)=\left(\mathrm{}\mathrm{\Omega }_{n^{},n}^{v,c}(k)\right)^{}$$
(125)
The Gaussian pump envelope function, $`(t)`$, is defined in Eq. (81) and the optical dipole matrix elements

$$d_{n,n^{}}^{c,v}(k)=\left(d_{n^{},n}^{v,c}(k)\right)^{}_\pi ^\pi 𝑑\varphi \widehat{ϵ}𝐝_{n,n^{}}^{c,v}(𝐤)$$
(126)
are angular averages in the $`xy`$ plane of the vector dipole matrices dotted into the polarization vector. From Eq. (68), the $`\varphi `$-dependence of $`𝐝_{n,n^{}}^{c,v}(𝐤)`$ only appears in $`𝐃_{j,j^{}}^{c,v}(\varphi )`$ and we can get the angular averages by setting $`\alpha (\varphi )=\alpha _{avg}=(1i)/3\pi `$ and $`\beta (\varphi )=\beta _{avg}=(1+i)/\pi `$ in the $`2\times 6`$ matrices $`D_X^{c,v}(\varphi )`$, $`D_Y^{c,v}(\varphi )`$ and $`D_Z^{c,v}(\varphi )`$ defined in Eq. (69).

The last terms in Eq. (II H) are similar in structure to the renormalization corrections in the Hartree-Fock energies but are more complicated due to mixing among subbands and involve the off-diagonal components of $`\mathrm{\Lambda }_{n,n^{}}^\alpha `$. The prime on the summation sign indicates that terms containing factors of $`N_{n,n^{}}^{\alpha ,\alpha ^{}}(k)`$ are excluded from the sum since these terms have already been incorporated in the renormalized Hartree-Fock energies in Eq. (115).

The coherent phonon amplitudes, $`D_q(t)`$, satisfy the driven harmonic oscillator equations

$`{\displaystyle \frac{^2D_q(t)}{t^2}}+`$ $`\omega _q^2D_q(t)={\displaystyle \frac{2\omega _q}{\mathrm{}}}{\displaystyle \underset{\alpha ,n,n^{},𝐤}{}}_{n,n^{}}^\alpha (k,q)^{}`$ (128)
$`\times \left\{N_{n,n^{}}^{\alpha ,\alpha }(k,t)\delta _{\alpha ,v}\delta _{n,n^{}}\right\},`$
subject to the initial conditions

$$D_q(t=\mathrm{})=\frac{D_q(t=\mathrm{})}{t}=0.$$
(129)
The closed set of coupled partial differential equations, (II H) and (128), for the carrier density matrices and coherent phonon amplitudes are converted into a set of coupled ordinary differential equations (ODE’s) by discretizing $`k`$ and $`q`$ and solving for $`N_{n,n^{}}^{\alpha ,\alpha ^{}}(k_i)`$ and $`D(q_i)`$ for each of the mesh points $`k_i`$ and $`q_i`$ . The resulting initial value ODE problem is then solved using a standard adaptive step-size Runge Kutta routine .

The phonon distributions do not appear in the coupled set of equations (II H) and (128). If necessary, they can be determined from $`N_{n,n^{}}^{\alpha ,\alpha }(k)`$ and the pair of equations

$$\frac{𝒩_q^{coh}}{t}=\frac{2}{\mathrm{}}\text{Im}\underset{\alpha ,n,n^{},𝐤}{}_{n,n^{}}^\alpha (k,q)B_qN_{n,n^{}}^{\alpha ,\alpha }(k)$$
(131)
and

$$\frac{B_q}{t}+i\omega _qB_q=\frac{i}{\mathrm{}}\underset{\alpha ,n,n^{},𝐤}{}_{n,n^{}}^\alpha (k,q)^{}N_{n,n^{}}^{\alpha ,\alpha }(k).$$
(132)
In Eq. (132), $`B_q(t)b_q(t)`$ satisfies the initial condition $`B_q(t=\mathrm{})=0`$. For the incoherent phonon distribution,

$$\frac{𝒩_q^{incoh}}{t}=0$$
(133)
so no incoherent phonons are generated and the incoherent phonon population maintains its initial thermal equilibrium distribution.

## III Loaded string model

The microscopic equations are rather daunting and detailed. In this section, we show how they can be simplified (under certain conditions) to a more tractable model, manely that of a driven uniform string, provided one uses the appropriate driving function, $`S(z,t)`$, which is nonuniform. The microscopics, including details of the superlattice band structure and photogeneration process are included within the driving function.

In our detailed numerical simulations, we use the full microscopic formalism discussed in the previous sections. However, we gain a lot of insight if we can deal with the lattice displacement, $`U(z,t)`$, directly. If we assume that the acoustic phonon dispersion relation is linear as in Eq. (85), then we find that $`U(z,t)`$ satisfies the loaded string equation

$$\frac{^2U(z,t)}{t^2}C_s^2\frac{^2U(z,t)}{z^2}=S(z,t)$$
(134)
subject to the initial conditions

$$U(z,t=\mathrm{})=\frac{U(z,t=\mathrm{})}{t}=0.$$
(135)
The LA sound speed, $`C_s`$, is defined in Eq. (85), and the driving function, $`S(z,t)`$, is given by

$`S(z,t)=`$ $`{\displaystyle \frac{1}{\mathrm{}}}{\displaystyle \underset{\alpha ,n,n^{}}{}}{\displaystyle \underset{𝐤,q}{}}\sqrt{{\displaystyle \frac{2\mathrm{}C_s|q|}{\rho _0V}}}_{n^{},n}^\alpha (k,q)^{}`$ (137)
$`\times \left\{N_{n,n^{}}^{\alpha ,\alpha }(k,t)\delta _{\alpha ,v}\delta _{n,n^{}}\right\}e^{iqz}.`$
One may question whether a linear phonon dispersion relation is valid in a superlattice. For small wavevector, $`q`$, for which elasticity theory holds, the dispersion relation for LA phonons in a superlattice is linear with a dispersion, $`\omega =\overline{C}_sq`$, where $`\overline{C}_s`$ is the ”average” sound speed of LA phonons in the well and barriers . This, in fact, has been experimentally verified in InGaN/GaN superlattice samples studied by C.-K. Sun et. al. .

Note that coherent acoustic phonon generation in a superlattice is qualitatively different than coherent optical phonon generation in a bulk system where only the $`q0`$ optic mode can be excited. As a result, both the amplitude, $`U(z,t)`$, and the Fourier transform of the amplitude, $`D_q(t)`$, for an optic mode in bulk satisfy a forced oscillator equation. For the nonuniform, multiple quantum well case, one can excite acoustic modes with $`q0`$. The Fourier transform of the amplitude $`D_q(t)`$ of a coherent acoustic phonon obeys a forced oscillator equation, but owing to the linear dependence of $`\omega (q)`$ on $`q`$ ,the amplitude itself, $`U(z,t)`$, obeys a 1-D wave equation with a forcing term, $`S(z,t)`$.

Another important point is that Eq. (134) can be taken to be a uniform string with a non-uniform forcing function. This is because the speed of sound is approximately the same in both the $`GaN`$ and $`InGaN`$ layers (a more detailed theory would take into account differences in the sound velocities in each layer). For propagation of acoustic modes one can neglect, to lowest order, the differences between the different layers (this is not true for the optic modes). The non-uniformity of the forcing function $`S(z,t)`$ results from differences in the absorption (not sound velocity) in the well and barrier layers and is therefore $`z`$ dependent. We thus see from Eq. (134) that understanding coherent acoustic phonons in multiple quantum wells is equivalent to understanding a uniform string with an inhomogeneous forcing term, $`S(z,t)`$, containing the microscopics.

To simplify Eq. (137), we neglect valence band mixing and assume the effective masses, sound speeds, and coupling constants are uniform over regions where $`S(z,t)0`$, i.e., in regions where carriers are being photogenerated. We also assume the pump pulses are weak enough so that screening of the piezoelectric interaction can be neglected. Finally, if the pump duration is long enough so that transient effects associated with photogeneration of virtual carriers can be ignored, then the off-diagonal elements of the carrier density matrices in Eq. (137) can be dropped. In this case, the driving function takes on the simple form

$$S(z,t)=\underset{\nu }{}S_\nu (z,t),$$
(138)
where the summation index, $`\nu `$, runs over carrier species, i.e., conduction electrons, heavy holes, light holes, and crystal field split holes.

Eq. (138) says that each carrier species makes a separate contribution to the driving function. The partial driving functions, $`S_\nu (z,t)`$, are

$$S_\nu (z,t)=\pm \frac{1}{\rho _0}\left\{a_\nu \frac{}{z}+\frac{4\pi |e|e_{33}}{ϵ_{\mathrm{}}}\right\}\rho _\nu (z,t),$$
(139)
where the plus sign is used for conduction electrons and the minus sign is used for holes. Here $`\rho _\nu (z,t)`$ is the photognerated electron or hole number density, which is real and positive, and $`\rho _0`$ is the mass density. We note that the loaded string equation for the propagation of coherent phonons together with the simplified driving function in Eqs. (138) and (139) have also been independently derived by other authors in the limit $`e_{33}=0`$.

In Eq. (139), the partial driving function for a given species is obtained by applying a simple operator to the photogenerated carrier density. This operator is a sum of two terms, the first due to deformation potential scattering and the second to piezoelectric scattering. The piezoelectric coupling constant, $`e_{33}`$, is the same for all carrier species, while the deformation potential, $`a_\nu `$, depends on the species. For conduction electrons, $`a_\nu =a_{c,z}`$, for heavy or light holes, $`a_\nu =D_1+D_3`$, and for crystal field split holes, $`a_\nu =D_1`$.

It is interesting to note that Planck’s constant does not appear in either the loaded string equation, (134), or in its associated driving function defined in Eqs. (138) and (139). Thus, we find that coherent LA phonon oscillations in MQW’s can be viewed as an essentially classical phenomenon, an observation that was made in the context of coherent LO phonon oscillations in bulk semiconductors by Kuznetsov and Stanton in Ref. .

The driving function, $`S(z,t)`$, satisfies the sum rule

$$_{\mathrm{}}^{\mathrm{}}𝑑zS(z,t)=0.$$
(140)
This is most easily seen from Eqs. (138) and $`(\text{139})`$, but it also holds for the general expression in Eq. (137). The significance of the sum rule is readily appreciated. After the pump dies away, the carrier density in Eq. (139), neglecting tunneling between wells, is essentially constant and thus $`S(z,t)`$ is time-independent. In the loaded string analogy, the integral of the driving function over position is proportional to the average force per unit length on the string. If this integral were non-zero, then the center of mass of the string would undergo a constant acceleration resulting in the buildup of an infinite amount of kinetic energy. Such an alarming result in the context of coherent LA phonons is precluded by the sum rule in Eq. (140).

For a given driving function, the wave equation (134), together with the initial conditions (135), can be solved for the coherent phonon lattice displacement using the Green’s function method. Thus,

$$U(z,t)=_{\mathrm{}}^{\mathrm{}}𝑑t^{}_{\mathrm{}}^{\mathrm{}}𝑑z^{}G(zz^{},tt^{})S(z^{},t^{}).$$
(141)
In our MQW diode model, the substrate is assumed to be infinite and the Green’s function in this case is just

$$G(z,t)=\frac{\mathrm{\Theta }(t)}{2C_s}\left\{\mathrm{\Theta }(z+C_st)\mathrm{\Theta }(zC_st)\right\},$$
(142)
where $`\mathrm{\Theta }(x)`$ is the Heaviside step function.

We note that the loaded string model described above is not restricted to the special case of an infinite substrate and can be extended to study the generation and propagation of coherent LA phonons in more complicated heterostructures. If the driving function, $`S(z,t)`$, due to photoexcited carriers is localized, then the assumptions leading to Eqs. (138) and (139) need only hold in those regions where $`S(z,t)`$ in nonvanishing. The wave equation applies in regions where the LA sound speed, $`C_s`$, is constant. Heterostructure in which the LA sound speed is piecewise constant have abrupt acoustic impedance mismatches which can be handled by introducing more complicated Green’s functions or by using other standard techniques. An example of such a problem would be a MQW structure embedded in a free standing substrate in which coherent LA phonons generated in the MQW could bounce back and forth between two parallel substrate-air interfaces.

## IV Results

In this section, we discuss simulations based on our microscopic theory of coherent LA phonon generation in a $`pin`$ diode structure with four periods of $`InGaN/GaN`$ MQW’s photoexcited by a Gaussian pump normally incident along the (0001) $`z`$-direction. The parameters for our numerical example are listed in Table II. The MQW dimensions and Gaussian pump parameters were chosen to match those typically encountered in room temperature pump-probe differential transmission measurements of coherent LA phonon oscillations carried out by C.-K. Sun, et. al.  on $`In_{0.06}Ga_{0.94}N/GaN`$ MQW structures having fourteen periods.

### A Bulk wurtzite bandstructure

Bulk wurtzite $`GaN`$ and $`InN`$ are direct gap materials with band gaps of $`3.4`$ and 1.95 eV, respectively. The bulk bandstructure of unstrained wurtzite $`GaN`$ is shown in Fig. 2. As can be seen from equations (5) and (16), the bandstructure is anisotropic and depends on $`k_z`$, the wavevector along the (0001) $`z`$-axis, and $`k_t`$, the wavevector within the $`xy`$ plane perpendicular to the $`z`$-axis. The effective mass conduction band is two-fold degenerate and has a parabolic dispersion with anisotropic effective masses $`m_z^{}=0.19`$ along the $`z`$-direction and $`m_{xy}^{}=0.18`$ in the $`xy`$ plane.

The two-fold degenerate valence bands are mixtures of heavy hole (HH), light hole (LH), and crystal-field splitoff hole (CH) character. At the zone center, the off-diagonal components of the $`3\times 3`$ upper and lower Hamiltonians in Eq. (16) vanish and the valence bands can be labeled according to their pure state wavefunctions at $`k=0`$. For the zone-center HH state, the degenerate wavefunctions are the basis states $`|v,1`$ and $`|v,4`$ defined in Eq. (II A). For the zone-center LH state, the wavefunctions are $`|v,2`$ and $`|v,5`$, and for the CH band the zone-center wavefunctions are $`|v,3`$ and $`|v,6`$. The heavy-hole effective masses along $`z`$ and $`xy`$ are $`m_z^{HH}=|A_1+A_3|^1=1.96`$ and $`m_{xy}^{HH}=|A_2+A_4A_5|^1=1.92`$ for heavy holes, $`m_z^{LH}=m_z^{HH}=1.96`$ and $`m_{xy}^{LH}=|A_2+A_4+A_5|^1=0.14`$ for light holes, and $`m_z^{CH}=|A_1|^1=0.14`$ and $`m_{xy}^{CH}=|A_2|^1=1.96`$ for crystal-field splitoff holes.

### B Pseudomorphic strain

Bulk $`GaN`$ and $`InN`$ have different lattice constants so when an (0001) $`InGaN`$ MQW structure is grown, a significant lattice mismatch occurs between the $`In_xGa_{1x}N`$ wells and $`GaN`$ barriers. For the $`InGaN`$ MQW diode specified in Table II, we assume pseudomorphic strain conditions. In a pseudomorphically strained device, the lattice constant throughout the MQW adjusts to the value of the lattice constant in the bulk $`N`$ and $`P`$ substrates in order to minimize the overall strain energy. In our simulated diode, the substrates are $`n`$\- and $`p`$\- doped $`GaN`$, so the lattice constant throughout the device takes on the $`GaN`$ value, i.e. $`a_0=3.189\AA `$. The non-vanishing position-dependent strain tensor components, $`ϵ_{xx}`$, $`ϵ_{yy}`$, and $`ϵ_{zz}`$, for the MQW diode, as computed from Eqs. (43) and (44), are shown in Fig. 3 as a function of $`z`$. Clearly, the $`GaN`$ barriers are unstrained since the $`N`$ and $`P`$ substrates are composed of $`GaN`$ and all the strain from the lattice mismatch is accomodated in the $`In_{0.06}Ga_{0.94}N`$ wells.

### C Built-in piezoelectric field

The presence of strain in the MQW’s results in the creation of a strain-induced polarization, $`P_z^0(z)`$, directed along $`z`$ as described by Eq. (45). The strain-induced polarization, in turn, results in a strong bult-in piezoelectric field which can be computed from Eqs. (46) and (47), given the strain field and the dc bias, $`V_A`$, applied across the diode. The computed strain-induced piezoelectric field, $`E_z^0(z)`$, and the piezoelectric confinement potential, $`V_{\text{piezo}}(z)`$, which result from the strain field in Fig. 3 are shown in Fig. 4. Prior to the application of the pump pulse, we assume the applied dc bias, $`V_A`$, has been adjusted so flat-band biasing in the diode is achieved, i.e. $`V_A`$ is such that the band edges seen in Fig. 4 are periodic functions of position.

Given the piezoelectric field and confinement potentials, position dependent band edges for the MQW can be computed. The conduction and valence band edges for our pseudomorphically strained MQW diode are shown as functions of position in Fig. 5. These are just the confinement potentials, $`V_\alpha (z)=V_{\alpha ,\text{gap}}(z)+V_{\text{piezo}}(z)`$, in the diode prior to photoexcitation. It is clear from Fig. 5 that the confinement of elecrons and holes in the MQW is mostly due to strong built-in piezoelectric fields which result in the triangular confinement potentials seen in each well.

### D Photogeneration of carriers

In our numerical example, we simulate photoexcitation of electrons and holes and the generation and subsequent propagation of coherent LA phonons in the hypothetical MQW diode when a Gaussian pump laser pulse is normally incicent along the $`z`$-axis. As seen in Table II, the Gaussian pump pulse is assumed to be left circularly polarized with a photon energy of $`3.21\text{eV}`$. The pump fluence is taken to be $`100.0\mu \text{J}/\text{cm}^2`$ and the Gaussian FWHM is taken to be $`180.0\text{fs}`$. The expermient is assumed to take place at room temperature.

In Fig. 6, the computed conduction and valence subband energies are shown as functions of $`k`$ for the $`In_xGa_{1x}N`$ diode. At the chosen pump energy of $`3.21eV`$, electrons from the first two valence subbands are excited into the lowest lying conduction subband.

The computed densities of photoexcited electrons and holes, neglecting and including Coulomb interaction effects, are shown as functions of position and time in Figs. 7 and 8, respectively and the total photoexcited electron density per unit area as a function of time is shown in Fig. 9. In Fig. 9 the pulse shape is shown for comparison. We find that including Coulomb effects decreases the total photogenerated carrier density. The electrons and holes screen the built-in piezoelectric field widening the effective band gap. This quantum confined Stark effect acts to suppress the photogeneration of carriers.

### E Generation of coherent phonons

The driving function, $`S(z,t)`$, for the driven string equation (134) is shown in Fig. 10 as a function of position and time. The driving function has units of acceleration and in Fig. 10, we compute $`S(z,t)`$ using the full microscopic formalism of Eq. (137).

For comparison, we also computed the driving function in the simplified loaded string model of Eqs. (138) and (139) using the carrier densities shown in Fig. 8 to facilitate the comparison. Since the photoexcited holes are predominantly a mixture of heavy- and light- holes, we use $`a_\nu =D_1+D_3`$ in computing hole deformation potential contributions in Eq. (139). The sum over species, $`\nu `$, then yields the total driving function

$`S(z,t)`$ $`={\displaystyle \frac{1}{\rho _0}}\left\{a_{c,z}{\displaystyle \frac{}{z}}+{\displaystyle \frac{4\pi |e|e_{33}}{ϵ_{\mathrm{}}}}\right\}\rho _{\text{elec}}(z,t)`$ (144)
$`{\displaystyle \frac{1}{\rho _0}}\left\{(D_1+D_3){\displaystyle \frac{}{z}}+{\displaystyle \frac{4\pi |e|e_{33}}{ϵ_{\mathrm{}}}}\right\}\rho _{\text{hole}}(z,t),`$
where $`\rho _{\text{elec}}(z,t)`$ and $`\rho _{\text{hole}}(z,t)`$ are the total conduction electron and valence hole densities plotted in Fig. 8. The resulting $`S(z,t)`$ is shown in Fig. 11.

By comparing Figs. 10 and 11, we see that for the diode structure and Gaussian pump used in our simulation the simplified loaded string model produces essentially the same results as those obtained using the full microscopic formalism.

Acoustic LA phonon generation due to the piezoelectric effect depends on the piezoelectric constant, $`e_{33}`$, the number of photogenerated electrons and holes, as well as the spatial separation of electron and hole densities brought about by the strong built-in piezoelectric field in the MQW’s. From Eq. (144), the piezoelectric contribution to the driving function is given by

$$S_{\text{piezo}}(z,t)=\frac{1}{\rho _0}\frac{4\pi |e|e_{33}}{ϵ_{\mathrm{}}}\left\{\rho _{\text{elec}}(z,t)\rho _{\text{hole}}(z,t)\right\}$$
(145)
In the absence of a built-in piezoelectric field (such as found in a square well with infinite barriers) we would have $`\rho _{\text{elec}}(z,t)\rho _{\text{hole}}(z,t)`$ and hence $`S_{\text{piezo}}(z,t)0`$ even for relatively large values of $`e_{33}`$. The built-in piezoelectric field serves to spatially separate the electrons and holes so that $`\rho _{\text{elec}}(z,t)\rho _{\text{hole}}(z,t)`$ and hence $`S_{\text{piezo}}(z,t)0`$. However, if the built-in piezoelectric field is too strong and the spatial separation of electrons and holes too large, then $`\rho _{\text{elec}}(z,t)\rho _{\text{hole}}(z,t)0`$. This is because the overlap between the conduction and valence envelope functions enters into the optical dipole matrix elements in Eq. (68). If there’s negligable overlap between electron and hole envelope functions due to strong piezoelectric fields then $`𝐝_{n,n^{}}^{cv}(𝐤)0`$, no electron-hole pairs are photogenerated and once again $`S_{\text{piezo}}(z,t)0`$.

The deformation potential contribution to the driving function is given by

$$S_{\text{def}}(z,t)=\frac{a_{c,z}}{\rho _0}\frac{\rho _{\text{elec}}(z,t)}{z}\frac{(D_1+D_3)}{\rho _0}\frac{\rho _{\text{hole}}(z,t)}{z}.$$
(146)
From Table I, the conduction electron deformation potential, $`a_{c,z}`$, is roughly twice the valence hole deformation potential, $`D_1+D_3`$. Thus, the two terms in Eq. (146) are of comparable magnitude. The first term, due to conduction electrons, gives rise to a contribution to $`S_{\text{def}}(z,t)`$ which is localized on the right side of each MQW while the second, due to valence holes, gives rise to contribution which is localized on the left hand side of each MQW.

In our simulation, we find that piezoelectric and deformation potential contriburtions to the driving function are comparable. This is seen in Fig. 12 where $`S_{\text{piezo}}(z)`$ and $`S_{\text{def}}(z)`$, along with their sum, are plotted at $`t=2\text{ps}`$. In this example, we find that $`S_{\text{piezo}}(z)`$ makes the dominant contribution to $`S(z,t)`$ as can be seen in Fig. 12.

The macroscopic lattice displacement, $`U(z,t)`$, and velocity field, $`V(z,t)`$, can be obtained from the coherent phonon amplitudes, $`D_q(t)`$, using Eqs. (96) and (97). Alternatively, they can also be obtained from the driving function, $`S(z,t)`$, through the Green’s function solution of the driven string equation. In Fig. 13, we plot the macroscopic lattice displacement, $`U(z,t)`$, and velocity field, $`V(z,t)`$, for coherent LA phonon modes generated by the driving function shown in Fig. 10.

### F Coherent phonon energy

From the lattice displacement, $`U(z,t)`$, we can obtain the total energy density per unit volume associated with coherent LA phonons,

$$_{LA}(z,t)=𝒯_{LA}(z,t)+𝒱_{LA}(z,t),$$
(148)
as the sum of a kinetic energy density term,

$$𝒯_{LA}(z,t)=\frac{\rho _0(z)}{2}\left(\frac{U(z,t)}{t}\right)^2,$$
(149)
and a potential energy density term,

$$𝒱_{LA}(z,t)=\frac{C_{33}(z)}{2}\left(\frac{U(z,t)}{z}\right)^2.$$
(150)
The coherent LA phonon energy per unit area, $`E_{LA}(t)`$, is obtained by integrating $`_{LA}(z,t)`$ over position, $`z`$:

$$E_{LA}(t)=_{\mathrm{}}^{\mathrm{}}𝑑z_{LA}(z,t).$$
(151)
The total energy density, $`_{LA}(z,t)`$, for coherent LA phonons as a function of position and time is shown in Fig. 14 and in Fig. 15 the total energy density is plotted as a function of position for equally spaced values of the time ranging from $`t=0`$ to $`t=8\text{ps}`$ in increments of $`2\text{ps}`$. The curves in Fig. 15 for different times have been offset to avoid overlapping. At short times, the evolution of the total phonon energy density is complicated, but the long times behavior, $`t6\text{ps}`$, can be easily understood.

As $`t\mathrm{}`$, a localized energy density appears in the MQW region due almost entirely to the potential energy term in Eq. (150). This is due to near steady-state loading by the driving function at long times. Assuming the driving function, $`S(z,t)`$, is approximately constant at long times, the loaded string equation, (134), can be integrated once in the steady state limit. We find the steady state solution

$$\frac{U(z)}{z}_{\mathrm{}}^z𝑑z^{}\frac{S(z^{})}{C_s^2}.$$
(152)
from which the long time behavior of the LA phonon energy density per unit volume in the MQW’s,

$$_{LA}^{\mathrm{}}(z)\frac{C_{33}(z)}{2A}\left(_{\mathrm{}}^z𝑑z^{}\frac{S(z^{})}{C_s^2}\right)^2,$$
(153)
can be obtained. The fact that the energy density in Eq. (153) is localized in the MQW’s follows directly from the sum rule (140) and is clearly seen in Figs. 14 and 15.

In addition to the localized energy density, which remains behind in the MQW’s, two propagating wave trains consisting of four pulses each are seen to exit the MQW region and travel off to infinity at the acoustic phonon sound speed, $`C_s`$. The distance between the pulses is just the inter-well separation distance. In these radiating wave trains, the kinetic and potential energy densities, $`𝒯_{LA}(z,t)`$ and $`𝒱_{LA}(z,t)`$, are found to be equal as one would expect.

The power spectrum of the coherent LA phonon energy density in $`q`$-space can be written in terms of the coherent phonon amplitudes $`D_q(t)`$. The power spectrum for the total coherent LA phonon energy density,

$$_{LA}(q,t)=𝒯_{LA}(q,t)+𝒱_{LA}(q,t),$$
(155)
is again the sum of a kinetic energy term,

$$𝒯_{LA}(q,t)=\frac{1}{2A}\frac{\mathrm{}}{\omega _q}\left|\frac{D_q(t)}{t}\right|^2,$$
(156)
and a potential energy term,

$$𝒱_{LA}(q,t)=\frac{1}{2A}\mathrm{}\omega _q\left|D_q(t)\right|^2.$$
(157)
The phonon energy density per unit area is obtained by summing the power spectrum over positive phonon wavevectors, $`q`$. Thus,

$$E_{LA}(t)=\underset{q>0}{}_{LA}(q,t).$$
(158)
The total energy density power spectrum for coherent LA phonons as a function of phonon wavevector, $`q`$, and time is shown in Fig. 16. The peak near $`q=0`$ is associated with buildup of the steady state energy density localized in the MQW region. Secondary peaks are seen near $`q_0=0.59nm^1`$ and twice this wavevector, i.e. $`q_1=2q_0=1.18nm^1`$. The wavevector $`q_0`$ corresponds to the the wavevector of the MQW period,

$$q_0=\frac{2\pi }{L_w+L_b},$$
(159)
where $`L_w`$ and $`L_b`$ are the well and barrier widths.

The total coherent LA phonon energy per unit area can be obtained from either Eqs. (151) or (158). In Fig. 17, we show the total phonon energy per unit area, $`E_{\text{LA}}(t)`$, as a function of time for the coherent LA phonons generated by the driving function shown in Fig. 10. The total energy per unit area is the sum of kinetic- and potential- energy terms. For comparison, the pulse shape is shown as a dotted line. It is clear from the figure that the buildup of energy in coherent LA phonons takes place on a time scale that is much longer than the pump duration. In addition, we see that the total energy buildup in the phonons saturates at around $`5\text{ps}`$ and that some strong but rapidly decreasing oscillations are superimposed on top of an increasing trend.

The saturation phenomenon results from the fact that we have a finite number of quantum wells and not an infinite superlattice. The results can best be explained in terms of the driven string equation. In general, the rate at which energy is fed into the phonon field per unit area is described by the energy equation,

$$\frac{E_{\text{LA}}}{t}=\rho _0_0^L𝑑zS(z,t)V(z,t),$$
(160)
in which $`S(z,t)`$ and $`V(z,t)`$ are the phonon driving function and velocity fields defined in Eqs. (137) and (97), and $`\rho _0`$ is the $`GaN`$ mass density used for computing the sound speed, $`C_s`$, in Eq. (85). The energy equation simply says that the rate at which energy is added to a driven string is proportional to the local force times velocity integrated over the length of the string. The integral in Eq. (160) vanishes when the transient velocity field, $`V(z,t)`$, exits the MQW region in which the driving function, $`S(z,t)`$, is localized. Thus, the time, $`t_{\text{sat}}`$, required for $`_{LA}`$ to saturate is just the time it takes for an LA sound wave to cross the MQW, i.e. $`t_{\text{sat}}W/C_s`$, where $`W`$ is the width of the MQW region over which the driving function is localized. In our example the LA sound speed is $`C_s=80\AA /\text{ps}`$ in $`GaN`$ and the MQW width (four well and three barrier layers) is $`W=381\AA `$, from which we obtain $`t_{\text{sat}}=4.8\text{ps}`$.

The oscillations of $`E_{LA}`$ observed in Fig. 17 reflect the number and periodicity of the diode MQW’s. The pump laser generates spatially periodic electron and hole distributions, as seen in Fig. 7, due to the fact that the pump photoexcites carriers in the wells but not the barriers. From each of the wells in the MQW, two sound pulses emerge traveling in opposite directions thus giving rise to an outwardly propagating velocity field pattern, $`V(z,t)`$, with four peaks traveling outward in each direction as seen in Fig. 13. The driving function, $`S(z,t)`$, on the other hand, is localized in the MQW’s and is relatively constant in time after the pump pulse dies out. The driving functions localized in each well do work on four traveling velocity disturbances, the one generated in the well itself as well as the ones generated in the three neighboring wells that subsequently pass by. This gives rise to the three peaks plus saturation plateau seen in Fig. 17. The time interval between coresponding peaks in adjacent wells is just the time it takes LA sound waves to travel between wells, i.e. $`t_{\text{period}}=(L_w+L_b)/C_s`$, where $`L_w`$ and $`L_b`$ are the well and barrier thicknesses. For the MQW structure, $`L_w=63.0\AA `$ and $`L_b=43.0\AA `$ and we have $`t_{\text{period}}=1.325\text{ps}`$ which agrees with the peak-to-peak time between the first and second peaks seen in Fig. 17.

## V Summary and Conclusions

We have developed a microscopic theory for the generation and propagation of coherent LA phonons in pseudomorphically strained wurzite (0001) $`InGaN/GaN`$ multi-quantum well (MQW) $`pin`$ diodes. Both $`GaN`$ and $`InN`$ have different lattice constants so that a significant mismatch occurs between the wells and barriers. The presence of strain in the MQW’s results in the creation of strain-induced built-in piezoelectric fields on the order of sevaral MV/cm which significantly alter the electronic and optical properties of the diode structure. In particular, the effective band gap can be lower than the band gap in unstrained $`InGaN`$ wells due to the presence of triangular piezoelectric potentials.

To a first approximation, the generation of coherent LA phonons is driven by optical photoexcitation of electron-hole pairs by an ultrafast Gaussian pump laser. Under typical experimental conditions, the propagation of coherent LA phonons is described by a uniform loaded string equation for the lattice displacement where the time- and position-dependent driving force on the string is a function of the photoexcited carrier density. This differs from coherent LO phonon oscillations in bulk systems where the coherent LO phonons obey a forced oscillator equation. Both deformation potential and piezoelectric coupling mechanisms contribute to the driving force. We find that deformation potential coupling contributes a driving force proportional to the derivative of the carrier density while piezoelectric coupling contributes a driving force proportional to the photoexcited carrier density.

We found that the driving term in the loaded string equation is suddenly turned on by rapid generation of electron-hole pairs by the pump and remains approximately constant theafter. This sudden displaceive loading of the string results in a new static equilibrium lattice displacement. This new static equilibrium displacement corresponds to a population of coherent LA phonons with $`q0`$. As the lattice adjusts to the new equilibrium, coherent LA phonons are transmitted in the positive and negative $`z`$-directions at the LA sound speed. These traveling coherent LA phonons are characterized by $`q2\pi /L`$ where $`L`$ is the superlattice period.

The formalism described here can be applied to the analysis of more complicated device geometries as well as more complicated laser pulse sequences. This gives a simpler method for calculating the coherent LA phonon generation in more complicated geometries and gives additional insight into the acoustic coherent response.

###### Acknowledgements.

The work of GDS and CJS was supported by the National Science Foundation through Grant No. DMR-9817828 and INT-9414072. The work of CSK was supported by Chonnam National University through a grant in the year 1999.